A TAILORED APPROACH TO ENGINEERING SOLID STATE SINGLE PHOTON 
SOURCES 
 
 
 
 
 
 
 
 
 
 
 
 
by 
 
RACHAEL KLAISS 
 
 
 
 
 
 
 
 
 
 
 
 
A DISSERTATION 
 
Presented to the Department of Physics 
and the Division of Graduate Studies of the University of Oregon 
in partial fulfillment of the requirements 
for the degree of 
Doctor of Philosophy  
 
June 2023 
 
  
DISSERTATION APPROVAL PAGE 
 
Student: Rachael Klaiss 
 
Title: A Tailored Approach to Engineering Solid State Single Photon Sources 
 
This dissertation has been accepted and approved in partial fulfillment of the requirements for 
the Doctor of Philosophy degree in the Department of Physics by: 
 
Hailin Wang Chairperson 
Benjamín Alemán Advisor 
Benjamin McMorran Core Member 
Christopher Hendon Institutional Representative 
 
and 
 
Krista Chronister Vice Provost for Graduate Studies  
 
Original approval signatures are on file with the University of Oregon Division of Graduate 
Studies.  
 
Degree awarded June 2023 
  
2 
  
 
 
 
 
 
 
 
 
 
 
 
 
 
© 2023 Rachael Klaiss  
  
3 
  
DISSERTATION ABSTRACT 
Rachael Klaiss 
Doctor of Philosophy 
Department of Physics 
June 2023 
Title: A Tailored Approach to Engineering Solid State Single Photon Sources 
Integrated quantum information technologies such as photonic circuits, quantum 
transducers, and magnetic sensors require robust single-photon sources in precise locations. 
Solid-state single photon emitters (SPEs) hosted by mid-bandgap defects in 2D material 
hexagonal boron nitride (hBN) are bright and stable at room temperature and demonstrate strong 
coupling to external fields, making them desirable candidates for quantum device applications. 
However, the specific atomic structure of hBN SPEs remains unidentified, making deterministic 
engineering a challenge. While recent studies have narrowed the range of possible defect 
candidates by demonstrating the role of carbon in hBN SPEs, the methods to engineer carbon-
based defects in hBN either produce randomly located emitters or require bottom-up crystal 
growth on structured substrates.  
We achieved patterned arrays of SPEs via focused ion beam (FIB) milling followed by 
chemical vapor deposition (CVD) of nanocrystalline graphite source for carbon diffusion, and 
found that both techniques are necessary for significant and repeatable creation of SPEs. This 
technique creates localized emitters with ten times the yield of carbon annealing alone. 
Furthermore, by adjusting the parameters of FIB exposure time and carbon annealing time, we 
found multiple different parameter combinations that successfully created SPEs, demonstrating 
the adjustability of this technique based on device application requirements. 
4 
  
Additionally, we performed atomic force microscopy to characterize the surface 
morphology of hBN regions patterned by Ga+ FIB to create SPEs at a range of ion doses and 
found that material swelling is strongly correlated to successful SPE creation. Furthermore, we 
simulated vacancy and impurity profiles to elucidate how Ga+ FIB patterning induces lattice 
damage in the form of vacancies, structural voids, and amorphous layers, creating a diffusion 
barrier to control the introduction of carbon impurities to engineer isolated SPEs with high 
resolution of process control. Our results provide novel insight into the formation of hBN SPEs 
created by high-energy, heavy-ion FIB that can be leveraged for monolithic hBN photonic 
devices and a wide range of low-dimensional solid-state SPE hosts. 
This dissertation includes previously published and unpublished coauthored material. 
5 
  
CURRICULUM VITAE 
 
NAME OF AUTHOR:  Rachael Klaiss 
 
 
GRADUATE AND UNDERGRADUATE SCHOOLS ATTENDED: 
 
 University of Oregon, Eugene, OR 
 University of California Berkeley, Berkeley, CA 
 
DEGREES AWARDED: 
 
 Doctor of Philosophy, Physics, 2023, University of Oregon 
 Bachelor of Arts, Physics, 2017, University of California Berkeley 
 
AREAS OF SPECIAL INTEREST: 
 
 Color centers in wide bandgap semiconductors 
 Defect engineering 
 2D materials 
 Materials characterization 
 
 
PROFESSIONAL EXPERIENCE: 
 
 Graduate Research Assistant, Alemán Lab, University of Oregon, Eugene, 
 Oregon, 2017-2023 
 
 Graduate Teaching Assistant, Department of Physics, University of Oregon, 
 Eugene, Oregon, 2017-2022 
 
 Graduate Teaching Assistant, Graduate Internship Program Photovoltaic and 
Semiconductor Device Processing Track, Knight Campus for Accelerating Scientific 
Impact, Eugene, Oregon 2020-2022 
 
 Diversity Equity and Inclusion Graduate Employee, Department of Physics, University of 
 Oregon, Eugene, Oregon, 2021 
 
 
GRANTS, AWARDS, AND HONORS: 
 
Future Leaders Society, Finalist for K. Patricia Cross Future Leaders Award, American 
Association of Colleges & Universities, 2023 
 
 Weiser Leadership Award, Department of Physics, University of Oregon, 2022 
 
6 
  
 Karfilis Leadership Award, University of Oregon Women in Graduate Sciences, 2021 
 
 Best Student Poster, “Hexagonal Boron Nitride Defect Investigation via Super Resolution 
      Imaging of Quantum Emission”, Oregon Center for Optical, Molecular, and Quantum 
      Science, 2019 
 
 Summer Research Institute Fellowship, CalTeach, University of California Berkeley, 
2015 
 
 
PUBLICATIONS: 
 
 R. Klaiss, J. Ziegler, D. Miller, K. Zappitelli, K. Watanabe, T. Taniguchi, and B.J. 
Alemán. Uncovering the morphological effects of high energy Ga+ focused ion beam     
milling on hBN single-photon emitter fabrication. J. Chem. Phys. 157, 074703 (2022). 
 
 J. Ziegler, R. Klaiss, A. Blaikie, D. Miller, V.R. Horowitz, and B.J. Alemán. 
     Deterministic quantum emitter formation in hexagonal boron nitride via controlled 
     edge creation. Nano Lett. 19, 2121-2127 (2019). 
 
 
7 
  
ACKNOWLEDGMENTS 
 
-- 
I am a blank canvas 
and the Mona Lisa 
if we were all born with a paintbrush in hand 
-- 
I owe the deepest gratitude to my advisor, Benjamín Alemán, for a phone call back in the 
Fall of 2016 that kickstarted this whole thing, and for all of the moments along the way when you 
were my biggest cheerleader and advocate. Thank you for making this work possible, for trusting 
me as I took this project down the path that spoke to me, and for validating my skills as a scientist, 
educator, and communicator. I move on to the next stage of my career with full confidence in my 
abilities because of your mentorship. 
To my lab mates and collaborators who taught me that science is a group effort both inside 
and outside the lab: Brittany Carter, for chem room cleaning sessions, office chats, and soup 
nights; Uriel Hernandez, for shared experiences and putting on so many hats; Rudy Resch, David 
Miller, Andrew Blaike, Kara Zappitelli, and Josh Ziegler, for laying foundations and teaching me 
how to conquer mountains and graduate school; Davis Austin, Skippy Clairmont, Paulina Nelson, 
Josh D’Mello, Halle Foster, Silkie Melloul, and Alena McVicker, for learning along with me; Sai 
Krishna Narayanan and Pratibha Dev, for calculations and reaching out; Viva Horowitz, for 
showing me what’s possible; Hailin Wang, Ben McMorran, and Chris Hendon, for serving on my 
committee and providing valuable insight beyond what the role requires. 
To the educators and scientists who shaped my passions for physics, unfurling its wonders 
for others, and ending the active exclusion of and harm done to many in the field: Mrs. Revelle, for 
“Best Girl Scientist”; Mr. Emmanuel and Professor Haefner, for seeing my interest and intuition 
for physics before I could see it myself; Frank Ferrone, for answering a cold email from a college 
8 
  
freshman and making my research career possible; Geoff Wehmeyer and the Dames Lab, for 
including me in part of a broader research community; Matt Gilbert and the Zettl Lab, for teaching 
me that I can always figure it out and try again; Robbie McGehee and Eli Mizrahi, for still being 
my models for a great mentorship group; Nima Dinyari, Fuding Lin, Stacey York, and Lynde 
Ritzow, for showing me that I can bring my whole self to work every day; Joel Corbo, Gina Quan, 
and the Access Network, for helping me find my niche and doing what you do; Tien-Tien Yu and 
Steven van Enk, for giving me agency, encouragement, and structure in the classroom; Stephanie 
Majewski, Laura Jeanty, and Tim Cohen, for listening; My mentees Francis Duey, Jasper Marcum, 
and Cristian Garcia, for sticking with me; The student leaders and organizers involved with 
UOWiP/PAGES, DuckREFS, North Star, and GTFF, for showing up and getting it done. 
To my friends that span states and stages of my life, I’m so glad I found you all when I 
did: Erin Coffman, Megan Bottino, and Sofia McCline, for being there through all evolutions; 
Lauren Capelluto, for being my Berkeley first and longest; Mike Lu and Jack Twilling, for Math 
1B, Napa wineries, and everything in between; Siena Ramirez and Andy Meislin, for helping us 
get it together; Lexi Goldwyn, Maddy Spiering, Megan Perley, Rayna Kim, and so many other 
Zeta sisters, for showing me the true meaning of #ZLAM and #4Years4Life; Nicole O’Shea, for 
keeping it classy and catty; Christine Abelyan, Melissa Hazelwood, Ramsey Karim, Joe Costello, 
Danny Tuthill, Rob Bentley, Mac Robertson, and Chris Agostino, for the reading room, study 
snacks, and support; The many iterations of Big Pink: Justin Kittell, Ethan Turner, Marija Glisic, 
Claire Pokorny, Andrew Lesak, and Cassidy Wagner, for always tuning in for Rachael’s Opinion 
Time and turning our little old house into a little old home; Michael Hurst, for sharing the 
longest and hardest leg of this journey with me; Grace Waddell and Steph VanBeuge, for 
showing me that physicists, chemists, and biologists aren’t so different after all; Daniel Moore 
9 
  
and Dylan Heim, for your campaign shenanigans and unique perspectives; Amy Soudachanh and 
Sean Brudney, for the chats; UO Phys Grads, especially the Fall 2017 cohort, for reading all my 
emails and showing me that grad school can be fun; Grace Lindquist, Kendra Lady, Chelsea 
Godon, and the many dancers I’ve had the privilege of working with over the last 18 months, for 
introducing me to your community and reminding me about balance. 
To my whole extended family that I know I’m making proud, and in particular: Aunt GG, 
for remembering more than I ever expected; Aunt Carole, for being a safe and comfortable place; 
Shannon, for popping in with a cat video and smile every so often; Aunt Kathy, for making a 
place at the table for me. 
To my parents, AnnaMarie and Harry: I have always believed I could do anything I set 
my mind to because of you, and I knew that you would support me every step of the way. This 
wouldn’t have been possible without either of you. I know you are proud of me, and I’m proud 
of me—so please give yourselves a space to be proud of yourselves. This is your 
accomplishment too. 
To my late baby brother, Harrison: I never doubted for a second that you’d want me to 
finish this thing, and honestly, I never even considered not finishing it. I couldn’t suddenly 
change what you thought was going to happen, so I kept going. I wish we could’ve had those 
thesis writing parties like we talked about, but I’m grateful for all you taught me about writing. 
Thank you for taking notes every time I tried to teach you about my work, because you said you 
wanted to understand my defense when it came around. Thank you for your unending curiosity, 
admiration, and love. I will carry them with me the rest of my life, in everything I do. 
 
If you are feeling alone and having thoughts of suicide—whether or not you are in crisis—or know 
someone who is, don’t remain silent. Call or text the 988 Suicide & Crisis Lifeline at 988. The Lifeline 
provides free and confidential emotional support to people in suicidal crisis or emotional distress 24 hours a 
day, 7 days a week, across the United States.   
10 
  
 
 
 
 
 
 
 
 
 
 
 
 
 
For Harrison. Thanks for stopping by. 
 
“Go confidently in the direction of your dreams. Live the life you have imagined” 
-Henry David Thoreau 
 
 
 
11 
  
TABLE OF CONTENTS 
 
Chapter       Page 
 
I.    INTRODUCTION .................................................................................................................. 20 
1.1. Solid State Single Photon Sources for Integrated Photonics ..................................... 20 
1.2. Two-Dimensional Wide-Bandgap Insulator Hexagonal Boron Nitride SPEs ........... 21 
1.3. Current Challenges in Engineering hBN SPEs .......................................................... 23 
 
II.   ENGINEERING CARBON-BASED QUANTUM EMITTERS IN HEXAGONAL 
      BORON NITRIDE.................................................................................................................. 25 
2.1. Introduction ................................................................................................................ 25 
2.2. Results and discussion ............................................................................................... 26 
2.3. Conclusions ................................................................................................................ 35 
2.4. Methods...................................................................................................................... 35 
2.5. Supplementary information ....................................................................................... 38 
2.6. Bridge ......................................................................................................................... 56 
 
III. UNCOVERING THE MORPHOLOGICAL EFFECTS OF HIGH ENERGY GA+ 
      FOCUSED ION BEAM MILLING ON HBN SINGLE-PHOTON EMITTER 
      FABRICATION ...................................................................................................................... 58 
3.1. Introduction ................................................................................................................ 58 
3.2. Results and discussion ............................................................................................... 60 
12 
  
Chapter       Page 
 
3.3. Conclusions ................................................................................................................ 76 
3.4. Methods...................................................................................................................... 77 
3.5. Supplementary information ....................................................................................... 88 
3.6. Bridge ......................................................................................................................... 90 
 
IV. MODELING EXTRINSIC IMPURITY DIFFUSION FOR COLOR CENTER 
      FORMATION ......................................................................................................................... 91 
4.1. Introduction ................................................................................................................ 91 
4.2. Results and discussion ............................................................................................... 92 
4.3. Conclusions .............................................................................................................. 100 
4.4. Methods.................................................................................................................... 101 
4.5. Supplementary information ..................................................................................... 103 
 
V.  CONCLUDING REMARKS ................................................................................................ 114 
 
APPENDICIES ........................................................................................................................... 116 
A. CORRELATIVE SCANNING NEAR-FIELD OPTICAL FLUORESCENCE 
     AND TOPOGRAPHY OF HBN COLOR CENTERS .............................................. 116 
B. SAMPLE IMAGES .................................................................................................... 121 
C. HOME-BUILT CONFOCAL MICROSCOPE USE INSTRUCTIONS ................... 125 
13 
  
Chapter       Page 
 
REFERENCES CITED ............................................................................................................... 130 
 
 
14 
  
LIST OF FIGURES 
 
Figure      Page 
 
Figure 2.1. SPE fabrication process flow diagram ....................................................................... 27 
Figure 2.2. Characterizing SPEs patterned in pristine hBN.......................................................... 28 
Figure 2.3. Four factor two level screening SPE yield response .................................................. 30 
Figure 2.4. Parameters and interactions that influence SPE yield ................................................ 31 
Figure 2.5. Optimizing SPE yield via FIB exposure and carbon annealing interaction ............... 33 
Figure 2.6. Tailoring antibunching at 𝑔𝑔2(0) and emission energy ............................................... 34 
Figure 2.7. Process flow diagrams for control experiments ......................................................... 39 
Figure 2.8. Raman characterization of diffusion source deposition ............................................. 41 
Figure 2.9. SEM and PL centroid analysis to localize patterned SPEs......................................... 41 
Figure 2.10. Fabrication parameter factor screening details ......................................................... 42 
Figure 2.11. Plots for single factor effects on SPE yield .............................................................. 47 
Figure 2.12. Plots for two-factor interaction effects on SPE yield ............................................... 48 
Figure 2.13. Plots for three-factor interaction effects on SPE yield ............................................. 49 
Figure 2.14. FIB exposure and carbon annealing interaction optimization parameters ............... 50 
Figure 2.15. Second order photon correlation fit parameter distributions .................................... 55 
Figure 2.16. Example zero phonon line spectra shapes for 𝑔𝑔2(0) > 0.5 ..................................... 56 
 
Figure 3.1. Testing FIB exposure on a single hBN flake .............................................................. 60 
Figure 3.2. AFM characterization of FIB exposure levels............................................................ 61 
Figure 3.3. Morphological characteristics calculated from AFM height data .............................. 62 
Figure 3.4. Average morphological characteristics for each FIB exposure level ......................... 63 
Figure 3.5. Summary of photophysical characterization by FIB exposure ................................... 66 
15 
  
Figure      Page 
 
Figure 3.6. Correlating material damage to SPE yield and photoluminescence ........................... 67 
Figure 3.7. Simulated FIB exposure dose profile ......................................................................... 68 
Figure 3.8. SRIM vacancy concentration profile .......................................................................... 69 
Figure 3.9. Estimating damage thresholds from simulated vacancies .......................................... 70 
Figure 3.10. Simulated swelling and milling via vacancy concentration thresholds .................... 72 
Figure 3.11. Illustration of vacancy distribution, void nucleation and growth, and milling ........ 73 
Figure 3.12. Analyzing AFM height data to calculate morphological characteristics .................. 80 
Figure 3.13. Stopping Range of Ions in Matter vacancy concentration ........................................ 83 
Figure 3.14. Simulating FIB patterning beam path....................................................................... 84 
Figure 3.15. AFM height scan data of pre-exposure hBN flake. .................................................. 88 
Figure 3.16. Mean calculated values for average heights and exposed edge surface area ........... 89 
Figure 3.17. Measured values for ZPL and 𝑔𝑔2(0) ........................................................................ 89 
 
Figure 4.1. Extended atomic force microscopy and photoluminescence characterization. .......... 93 
Figure 4.2. Illustrations demonstrating multilayer diffusion model formulation ......................... 95 
Figure 4.3. Evaluated model at tested anneal times ...................................................................... 96 
Figure 4.4. Model parameter dependence of SPE supersaturation characteristics ....................... 97 
Figure 4.5. Full AFM characterization at 4 different anneal times ............................................. 104 
Figure 4.6. Multi-layer diffusion carbon concentration profiles ................................................. 105 
Figure 4.7. Illustration of Boltzmann distribution model for defect formation .......................... 106 
Figure 4.8. Density functional theory (DFT) calculations .......................................................... 113 
 
Figure A.1. Schematic of WITec SNOM adapted for fluorescence ............................................ 118 
16 
  
Figure      Page 
 
Figure A.2. Improving confocal resolution with SNOM ............................................................ 118 
Figure A.3. Simultaneous correlation of hBN color center SNOM to AFM .............................. 119 
Figure A.4. Correlative imaging of an hBN quantum emitter .................................................... 120 
 
Figure B.1. Optical and SEM of factor screening samples ......................................................... 121 
Figure B.2. Confocal PL scans of factor screening samples ....................................................... 122 
Figure B.3. Optical and SEM of optimization samples .............................................................. 123 
Figure B.4. Confocal PL of optimization samples ...................................................................... 124 
 
 
17 
  
LIST OF TABLES 
 
Table      Page 
 
Table 2.1. Control Experiments Results ....................................................................................... 39 
Table 2.2. Control Experiments – SPE Yield Mean Comparison Tests ....................................... 40 
Table 2.3. Factor Screening Experiments Results ........................................................................ 45 
Table 2.4. Thin, 500 nm, 3.8 ms, 5 hr – Mean SPE Yield Comparison Tests .............................. 45 
Table 2.5. Factor Screening Full Model Effects Test for Number of SPEs.................................. 46 
Table 2.6. Ion Fluences for Optimization FIB Exposure Levels .................................................. 50 
Table 2.7. Optimization Experiments Results .............................................................................. 51 
Table 2.8. Optimization Experiments – Mean SPE Yield Comparison Tests .............................. 52 
Table 2.9. Optimization Experiments – Mean Region PL Comparison Tests .............................. 53 
Table 2.10. Top 5 Optimization Experiments – Mean 𝑔𝑔2(0) Comparison Tests ......................... 54 
Table 2.11. ZPL by FIB Exposure – Mean Comparison Tests ..................................................... 55 
Table 2.12. ZPL by Carbon anneal time – Mean Comparison Tests ............................................ 55 
 
Table 3.1. Morphological Characteristics Calculated from AFM ................................................ 63 
 
Table 4.1. Average Mill Depth by FIB Exposure time– Mean Comparison Tests ..................... 103 
Table 4.2. Average Swell Height by FIB Exposure time– Mean Comparison Tests ................. 103 
Table 4.3. Model Parameters Used for Figure 4.3 ...................................................................... 105 
Table 4.4. Model Parameter 24 Results –Time of Max SPE Probability .................................... 109 
Table 4.5. Model Parameter 24 Results –SPE Onset Time ......................................................... 109 
Table 4.6. Model Parameter 24 Results –Color Center Saturation Time .................................... 110 
Table 4.7. Factor Screening Full Model Effects Test for Time of Max SPE Probability ........... 110 
18 
  
Table      Page 
 
Table 4.8. Factor Screening Full Model Effects Test for SPE Onset Time ................................ 111 
Table 4.9. Factor Screening Full Model Effects Test for Color Center Saturation Time ........... 111 
19 
 
 
CHAPTER I 
 
INTRODUCTION 
 
1.1. Solid State Single Photon Sources for Integrated Photonics 
 Advances in quantum technology that promise more complex computation, secure 
communication, and accessible sensing rely on scaling the unique capabilities of quantum 
emitters (QEs)—atom-like systems such as quantum dots and solid-state defects— to practical 
devices, while connecting sources that operate at multiple wavelength ranges in a quantum 
network.1 Promising qubit candidates for communicating between quantum technologies over 
long distances are single photon emitter (SPE) polarization qubits with tunable wavelength and 
stable fluorescence2, while spin qubits with SPE readout via optically detectable magnetic 
resonance (ODMR) are desirable for quantum sensing and transduction applications. Currently, 
classical light mediates communication between various systems, for example in the 
transmission of information in fiber optic cables for television or internet. SPEs can be coupled 
to optical fibers for transmission, demonstrating an avenue for sending information in a quantum 
network.3 
 However, a requirement for telecommunication is microwave photon conversion, which is 
yet to be demonstrated in a tunable and stable quantum light system.4 For example, silicon-based 
QEs (e.g. G, W, and T-centers) are desirable for photonics applications due to the mature 
understanding of silicon-based electronics fabrication and emission in the telecom band, but 
these systems require cryogenic temperatures for operation and photon conversion platforms 
require increasingly complex fabrication, making practical applications elusive.5–7  
2 0 
  
 A way to convert QEs outside the telecom band to telecommunication frequencies is through 
coupling to mechanical modes, which operate at high frequencies. Solid-state QEs that operate at 
room temperature in the visible band are an attractive way to couple optical readout signals to 
spin, electronic, or phononic degrees of freedom when applied to hybrid and heterogenous 
integrated photonics, allowing for communication and translation between different components 
in a range of applications.8,9 These quantum applications require stable, indistinguishable SPEs 
and scalable placement, the limitations of which depend highly on the solid-state source and 
fabrication technique.  
 Optomechanical devices have demonstrated coupling between mechanical modes and light in 
quantum dots and NV centers in diamond via strain mediation and electron spin centers in 
diamond via surface acoustic waves.10–12 Two-dimensional (2D) materials such as graphene, 
hexagonal boron nitride (hBN), and transition metal dichalcogenides (TMDCs) (e.g. MoS2, 
WSe2, etc.) have a wide parameter space in the mechanical regime due to their low mass and 
high capacity for strain, with mechanical resonances in the telecommunication frequency 
regime.13 SPEs from 2D material hosts therefore promise strong single photon optomechanical 
coupling, with an ideal source operating at room temperature to decrease the reliance on 
cryogenics.4 
1.2. Two-Dimensional Wide-Bandgap Insulator Hexagonal Boron Nitride SPEs 
 hBN SPEs are a desirable platform for integrated photonics and quantum technologies 
applications due to the utility of hBN in many existing technologies as an insulating layer, as 
well as unique SPE properties when compared to both bulk and other 2D hosts. hBN is a wide 
bandgap semiconductor and emitters are hosted in mid-bandgap states caused by atomic vacancy 
defects.14 States that are within the band gap can be excited without also mediating exciton 
21 
  
recombination optical emission. Therefore, hBN single photon sources do not require cooling to 
cryogenic temperatures to exhibit bright, stable, and narrow spectral lines of single photon 
emission in contrast to SPEs from TMDCs.15,16  
 The characteristic spectral line shape of a SPE is called the zero-phonon line (ZPL) due to 
emission not requiring phonon assistance to reach the energy between excited and ground states. 
The ZPL of different hBN SPEs can range over hundreds of nanometers.17,18 The range in ZPL 
has been attributed to various vacancy-antisite and vacancy-interstitial defect complexes such as 
NBVN and CNVB as well as to a range of vacancy defect complexes with anywhere from one to 
sixteen vacancies per defect.19,20 Substrate interactions are also expected to induce ZPL variation 
even for the same defect due to local strain and charge gradients.21 
 hBN SPEs can be heavily influenced by the host environment. Blinking emitters have been 
observed and characterized in numerous systems.22–24 Consistent with other quantum systems 
such as quantum dots and single molecules, the dynamics of blinking follow a power law fit, 
which attributes the blinking to a population trap or accessing a triplet dark state by spontaneous 
fluctuations.25–27 Spectral diffusion has also been attributed to random changes to the charge 
surrounding the single photon source, inducing a strong local Stark shift.28 These random 
fluctuations are hypothesized as due to charge diffusing from the sample substrate or the 
environment.17-23   
 Similar to the NV-center in diamond, V -B  in hBN exhibits optically detectable magnetic 
resonance (ODMR) that has been demonstrated as a spin sensor for ferromagnetic materials.29,30 
Additionally, hBN defects not specifically identified as V -B  are sensitive to strain31 and external 
electric fields, exhibiting strong stark shifts32, while carbon-based defects of various complex 
types and sizes are expected to be spin sensitive33,34. Furthermore, recent work has shown that 
22 
  
twisted bilayers of hBN can host emitters, which are also magnetically sensitive and tunable 
when exposed to electron beams35,36. 
1.3. Current Challenges in Engineering hBN SPEs 
 To be useful for applications, e.g. into arrays for photonic circuits, or into structures like 
nanoelectromechanical systems or optical cavities, it is essential that the SPEs be patterned at 
well-defined positions8,9. There is some success in large-scale activation of hBN SPEs with bulk 
focused ion beam (FIB)37,38 and plasma irradiation39, but these methods generate randomly 
positioned emitters in a large area and tend to produce ensembles that cannot be isolated into 
single photons. Moreover, carbon implantation (and deposition) of hBN revealed the likely role 
of carbon in emitting defects—carbon-vacancy centers—but was also only successful at 
producing randomly located emitters40. Similar results have been attained by bottom-up growth40 
or annealing of hBN crystals in the presence of atomic carbon41.  
 Patterning of SPEs in hBN has been realized via low electron beam (e-beam) irradiation42 
and nanoindentation with an atomic force microscopy (AFM) tip43, and was achieved by the 
authors in earlier work23 by milling nanoscale holes and edges into hBN followed by annealing, 
but the replication success and efficiency highly depended on material properties and even local 
variations of a given sample; for example, SPE nanofabrication was not possible in high-quality 
exfoliated hBN. Ultimately, challenges in engineering and localizing hBN SPEs persist because 
the specific atomic structure of defects varies greatly and the physical mechanism for defect 
formation enabled by patterning and doping has not been thoroughly investigated.   
 In order to address these challenges, we developed an adjustable technique using FIB 
patterning to create edges and vacancies that allow for diffusive impurity doping to create 
visible-range, room temperature SPEs in hBN, likely in the form of carbon-vacancy defect 
23 
  
complexes. By characterizing both photophysical properties of the SPEs and morphological 
features of the patterned hBN regions for multiple FIB and impurity diffusion parameter 
combinations, we propose a model for the physical mechanism of defect formation with this 
technique. Understanding the physical mechanism of defect formation provides a framework to 
optimize fabrication of not just hBN SPEs, but other impurity defect-based solid state SPEs that 
can be used in a wider range of photonics applications. Furthermore, this adjustable engineering 
technique can be leveraged to minimize negative effects of substrate interactions and local 
charge fluctuations on spectral diffusion and photon indistinguishability, thereby improving the 
viability of solid-state SPEs for scalable integrated photonics. 
 The work described in this dissertation is based largely on co-authored material. Chapters II 
and IV contain unpublished co-authored material with contributions from Josh Ziegler, David 
Miller, Kara Zappitelli, Kenji Watanabe, Takashi Tanaguchi, Sai Krishna Narayanan, Pratibha 
Dev, and Benjamín Alemán. Chapter III contains published co-authored material with 
contributions from Josh Ziegler, David Miller, Kara Zappitelli, Kenji Watanabe, Takashi 
Tanaguchi, and Benjamín Alemán. 
 
 
24 
 
 
CHAPTER II 
 
ENGINEERING CARBON-BASED QUANTUM EMITTERS IN HEXAGONAL  
BORON NITRIDE 
 
This chapter contains in-preparation co-authored material; it has been reproduced with 
permission from: R. Klaiss, J. Ziegler, D. Miller, K. Zappitelli, K. Watanabe, T. Tanaguchi, S.K. 
Narayanan, P. Dev, B. J. Alemán, “An adjustable focused ion beam – impurity diffusion 
approach to pattern hexagonal boron nitride quantum emitters.” In this work, I contributed to 
designing the research, performing the research, analyzing the data, and writing the paper. 
 
2.1. Introduction 
 Solid-state single photon sources are essential platforms for all-optical quantum computing, 
quantum communication, and as sources for quantum optics experiments, and are desirable due 
to ease of on-chip integration and coupling to the material host as a means of control and 
tunability1,8,9,44,45. Single photon emitters (SPEs) have been discovered in various 2D materials, 
which offer particularly straightforward integration, strong surface-sensitivity, and a wider range 
of strain-tunability than their bulk counterparts3,44,46. Hexagonal boron nitride (hBN) SPEs are 
particularly effective due to 2D material host and defect states hosted within the wide bandgap, 
leading to robustness at room temperature47,48.   
 Practical application of SPEs for photonics applications requires straightforward, 
deterministic placement of defects, but many techniques currently require structured substrates49 
or niche ion implantation systems50. In our previous work utilizing focused ion beam (FIB) 
2 5 
  
patterning to address this challenge, milling into poorer quality chemical vapor deposition 
(CVD) grown hBN worked the best23. This CVD grown hBN likely had high levels of 
amorphous carbon from the multi-step polymer transfer process, further indicating that the 
presence of carbon was essential to creating the emitters. Guided by these observations, we 
conjecture that nanofabrication was enabled by diffusive driving of carbon into vacancies located 
near edges of holes.   
 To test this idea, we combine FIB nanopatterning with CVD of nanocrystalline-graphite 
(NC-Gr) as a carbon diffusion source to nanofabricate quantum emitters in hBN. We find SPEs 
in positions localized to our FIB patterning process and use factor screening to determine which 
process parameters are important for SPE creation. We postulate that this process works when 
elemental carbon diffuses into vacancies created by FIB, entering at the edges from the surface 
of the hBN and perform process optimization experiments to test this model. 
2.2. Results and discussion 
 The steps to fabricate samples used in this work are outlined in Figure 2.1, with more details 
in 2.4 Methods. To create edge defects at deterministic locations, we patterned arrays of 500 nm 
circles separated by a 1 μm pitch with beam parameters determined in our previous work as 
optimal for SPE creation.23,51 We then performed CVD of NC-Gr via non-catalytic ethylene 
decomposition (confirmed via Raman spectroscopy [Supplementary information Figure 2.8]) 
at 1000°C, the annealing temperature used by Lyu et al. to activate hBN SPEs in a carbon rich 
environment41. To remove residual surface NC-Gr prior to characterization, we annealed the 
samples under flowing oxygen.  
 We identified spots of photoluminescence (PL) and measured the degree of second-order 
coherence (𝑔𝑔2(𝑡𝑡)) of emitted photons via Hanbury-Brown and Twiss experiment in a home-built 
26 
  
Figure 2.1. SPE fabrication process flow diagram starting with hBN transfer, followed by 
cleaning anneal, Ga+ FIB patterning, carbon annealing shown via ethylene decomposition, 
cleaning of carbon diffusion source shown via oxidation, and ending with characterization of 
single photon emission via 532 nm excitation in home-built confocal microscope 
confocal microscope with 532 nm excitation [Figure 2.2(a), Supplementary information 
Figure 2.9, 2.4 Methods]. These measurements were used to determine if PL sources exhibited 
antibunching, and also yielded their average photon number, nonradiative and radiative lifetimes, 
and bunching amplitude. These parameters were extracted from antibunching measurements by 
fitting our data to the following second order correlation function model for a three-level system: 
𝑔𝑔2(𝑡𝑡) = 1 − 𝜌𝜌2
|𝑡𝑡| |𝑡𝑡|
+ 𝜌𝜌2 �1 − (1 − 𝑎𝑎)𝑒𝑒− �𝜏𝜏1 + 𝑎𝑎𝑒𝑒− �𝜏𝜏2� (2.1) 
where 𝑎𝑎 is the bunching amplitude, 𝜏𝜏1 is the nonradiative lifetime, 𝜏𝜏2 is the radiative lifetime, 
and 𝜌𝜌2 = 1 − 𝑔𝑔2(0) where 𝑔𝑔2(0) is the degree of antibunching. Photon number in a given milled 
hole was determined by binning 𝑔𝑔2(0) values according to: 
𝑔𝑔2
1
(0) = 1 − (2.2) 
𝑛𝑛
27 
  
Figure 2.2. Characterizing SPEs patterned in pristine hBN (a) Representative confocal PL 
map of a color center array on an hBN flake patterned with the described process. Dashed lines 
indicate hBN flake edges. Scale bar: 4 µm (b) Zero-phonon line (ZPL) at ~2.15 eV and 
corresponding phonon sideband (PSB) emission of 160 meV for a representative single photon 
emitter created with this process. Inset: 𝑔𝑔2(𝑡𝑡) for the emitter with 𝑔𝑔2(0) =  0.24 , well below 
the 𝑔𝑔2(0) < 0.5 threshold for single photon emission. (c) Yield of patterned sites that hosted 
single photon emitters on samples subjected to 4 different treatment combinations. Error bars 
show standard error on mean. 
We use the accepted antibunching threshold for SPE (𝑔𝑔2(0) < 0.5) to categorize emitting sites 
as SPE [Figure 2.2(b) Inset]. We denote any emitting sites that exhibit antibunching (with 
0.5 < 𝑔𝑔2(0) < 1) as quantum emitters (QE). We categorize any bright emitting sites that do not 
exhibit antibunching (𝑔𝑔2(0) = 1) and do not photobleach as color centers. We further 
characterized the quantum emitting nature of the spots by acquiring spectra, where we observe 
28 
  
the expected zero-phonon line (ZPL) and corresponding phonon sideband (PSB) for quantum 
emitters [Figure 2.2(b)]. 
 FIB patterning followed by carbon annealing produced SPEs with higher yield than only FIB 
patterning or only carbon annealing [Figure 2.2(c)]. The combined process (optimized 
parameters) yielded 0.1146±0.0094 SPEs/patterned site (0.1146±0.0094 SPEs/μm2), 
approximately 9.2 times more SPEs than carbon annealing alone (0.0125±0.0071 SPEs/ μm2) 
[Supplementary information Table 2.2 Student’s t test: p < 0.0001]. This observed difference 
in yield is in qualitative agreement with our density functional theory (DFT)-based calculations 
of the defect formation energies for a representative set of carbon-based defect complexes, which 
consist of substitutional carbons adjacent to different vacancies. The formation energy of these 
defect complexes in pristine hBN are much larger than the energies required to incorporate 
carbon once defects have already been created in hBN [see Chapter IV Supplementary 
information Figure 4.8(a)]. The latter scenario emulates experimental process of FIB followed 
by CVD. 
 To determine which fabrication parameters of the combined FIB-CVD technique contributed 
to the improved SPE yield, we performed a two-level factor screening on four targeted factors: 
hBN thickness (14.59 ± 4.81 nm, 28.47 ± 5.42  nm), FIB pattern size (250 nm, 500 nm), FIB 
exposure time (2.5 ms, 3.8 ms), and carbon annealing time (20 min, 300 min) [Supplementary 
information Figure 2.10]. The parameter combination from this factor screening that produced 
the highest SPE yield was: 500 nm pattern diameter holes exposed for 3.8 ms on a 16 𝑛𝑛m flake 
with 5 hours of carbon annealing time [Figure 2.3]. Samples fabricated with this combination of 
parameters yielded 0.0673±0.0058 SPEs/patterned site, which was 3.4 times higher than the 
29 
  
Figure 2.3. Four factor two level screening SPE yield response (a) Confocal PL scan for a 
sample processed under the best conditions for SPE creation during factor screening. SPEs are 
circled in solid white line, QEs are circled in dashed white line. Scale Bar: 4 µm (b) All factor 
combinations tested, showing that the low thickness, high FIB pattern size, high FIB exposure 
time, high carbon annealing time combination lead to the highest SPE yield (Tukey-Kramer HSD) 
second-best parameter combination [Supplementary information Table 2.4 Student’s t test: 
p = 0.0022]. 
 We found that hBN thickness had the strongest single-factor effect on SPE yield. Thinner 
flakes hosted 9.1 times higher SPE yields than thicker flakes, regardless of the chosen FIB 
exposure time, pattern size, and carbon annealing time [Supplementary information Table 2.5 
F-test: p = 0.004] [Figure 2.4(a)]. The FIB exposure time was the only other single factor that 
influenced the SPE yield, with “high” exposure producing 4 times higher SPE yields than “low” 
exposure [Supplementary information Table 2.5 F-Test: p = 0.028] [Figure 2.4(b)]. 
Somewhat surprisingly, annealing samples for 5 hours did not, on average, produce a higher or 
lower SPE yield than annealing samples for 20 minutes [Supplementary information Table 2.5 
F-test: p = 0.544].  
 When we accounted for parameter interactions, we found that increasing the anneal time 
increased the SPE yield only in the high exposure case, whereas it decreased the SPE yield in the 
low exposure case, indicating a strong FIB exposure time and carbon anneal time interaction 
effect [Supplementary information Table 2.5 F-test, p = 0.042] [Figure 2.4(c)]. Furthermore, 
30 
  
Figure 2.4. Parameters and interactions that influence SPE yield (c) Effect of hBN thickness 
on SPE yield. Closed yellow circles show least squares mean SPE yield for thin hBN and thick 
hBN. (d) Effect of FIB exposure time on SPE yield. Closed yellow circles show least squares 
mean SPE yield for low exposure and high exposure. Error bars show standard error on the 
mean. (e) Effect of FIB exposure-Carbon annealing interaction on SPE yield. Closed yellow 
circles show least squares mean SPE yield for low FIB exposure interactions and closed green 
circles show least squares mean SPE yield for high FIB exposure interactions. (f) Effect of hBN 
thickness-FIB exposure-Carbon annealing interaction on SPE yield.  Left panel shows thin hBN 
interactions and right panel shows thick hBN interactions. Closed yellow circles show least 
squares mean SPE yield for low FIB exposure interactions and closed green circles show least 
squares mean SPE yield for high FIB exposure interactions. Error bars show standard error on 
the mean. 
on thinner flakes, this effect was enhanced, while on thicker flakes, this effect was reversed, 
indicating a strong hBN thickness, FIB exposure time, and carbon annealing time interaction 
effect [Supplementary information Table 2.5 F-test, p = 0.005] [Figure 2.4(d)].  
31 
  
 Based on the results of this factor screening, we established that the optimal hBN thickness to 
fabricate SPEs with this process is 10-19 nm. Furthermore, since SPE yield depends strongly on 
the interaction of FIB exposure and carbon annealing times, we can more finely tune this 
technique by sweeping those parameters. We postulate that balancing the concentration of 
vacancy defects and carbon impurity atoms is essential to engineering carbon-based quantum 
emitters in pristine hBN and increasing SPE yield. 
 To optimize emitter yield with this technique, we explored the interaction of FIB exposure 
and carbon annealing for four FIB exposure times (1.3 ms, 2.5 ms, 3.8 ms, 5 ms) and five 
carbon annealing times (20 min, 75 min, 150 min, 300 min, 720 min) in order to cover the 
range within the two extreme values chosen for the factor screening. For more details about the 
experimental design, refer to Supplementary information Figure 2.14. 
 Our technique successfully created color centers for 13 of the 18 FIB exposure-carbon 
annealing time combinations. We observed that longer carbon annealing decreased SPE yield for 
shorter FIB exposure times and increased SPE yield for longer exposure [Figure 2.5(a)]. This 
behavior agrees with the factor screening results. We found the highest SPE yields on regions 
exposed for 5.0 ms and annealed for 5 hr, with average SPE yield of 0.1146±0.0094 
SPEs/patterned site—approximately twice the yield produced by the best factor screening 
parameters [Supplementary information Table 2.8 Student’s t test: p < 0.0001]. Four 
additional FIB exposure-carbon annealing time combinations yielded quantum emitters at 
patterned sites [Supplementary information Table 2.7 Dunett test with no FIB no CVD 
control: p < 0.0001], (B) 3.8 ms:2.5 hr, (C) 3.8 ms:5 hr, (D) 2.5 ms:20 min, (E) 5.0 ms:12 hr 
[Figure 2.5(b)]. 
32 
  
Figure 2.5. Optimizing SPE yield via FIB exposure and carbon annealing interaction (a) 
Interaction plots showing mean SPE and mean QE yield for samples processed under the same 
optimization conditions. The conditions that on average produce the highest number of SPEs and 
QEs are labeled with A being the conditions that produces the highest number of emitters on 
average, B the second highest, and so on. (b) Representative PL maps for the conditions that 
produce the 4 highest number of emitters on average. Scale Bar: 2 µm 
 While the current, most persistent need for integrating quantum emitters into photonic 
technologies is to optimize the yield of SPEs at deterministic locations, it is also essential to 
characterize the influence of the fabrication technique on optical characteristics such as single-
photon purity and ZPL energy that determine the utility of the emitters in these technologies. To 
determine if the particular combination of FIB exposure and carbon annealing times impact the 
single-photon purity, we characterized 157 emitters with 𝑔𝑔2(0) < 1 on regions subjected to 12 
of the 18 possible exposure and annealing combinations, with 65.6% of those emitters on regions 
subjected to the top five combinations listed above as producing the highest QE yields. We found 
that emitters fabricated with conditions (A) and (B) had average 𝑔𝑔2(0) < 0.5, while regions 
subjected to conditions (C), (D), and (E) had average 𝑔𝑔2(0) > 0.5 [Figure 2.6(a)]. Furthermore, 
we observed that for constant FIB exposure level, longer carbon annealing time increased the 
average 𝑔𝑔2(0) [Supplementary information Table 2.10 B-C Student’s t test 0.436 ± 0.038 – 
0.631 ± 0.033: p = 0.0002, A-E Student’s t test 0.477 ± 0.028 – 0.698 ± 0.028: p < 0.0001]. 
33 
  
Figure 2.6. Tailoring antibunching at g2(0) and emission energy (a) 𝑔𝑔2(0) value distribution 
for parameter combinations that yielded the most QEs, where A, B, C, D, E labels correspond to 
the same labels in Figure 2.3. Diamond with error bars for each experiment represent the mean 
𝑔𝑔2(0) and standard error on mean. Dashed lines mark single-photon emitter threshold  
𝑔𝑔2(0) =  0.5  (b) ZPL distribution by FIB exposure level. Purple diamond with error bars for 
each level represent mean ZPL and standard error on the mean. 
These results demonstrate that we can tailor the quantum purity of emitters fabricated with this 
technique by adjusting the exposure time and carbon annealing time. 
 To determine if we can similarly tailor the ZPL, we collected spectra from 143 color centers 
and observed ZPL within the expected visible range for hBN QEs. We found that ZPL and ZPL 
distributions varied due to FIB exposure time [Figure 2.6(b)]. Although the average ZPL of 
emitters patterned with the 2 highest FIB exposure times were similar [Supplementary 
information Table 2.11 5.0 ms-3.8 ms Student’s t test 2.22 ± 0.02 – 2.19 ± 0.01: p = 0.312], the  
average ZPL for emitters patterned with the second FIB level shifted to 2.16 ± 0.02 eV 
[Supplementary information Table 2.11 3.8 ms-2.5 ms Student’s t test: p = 0.0164, 5.0 ms-2.5 
ms Student’s t test: p = 0.0022]. This difference could indicate variations in the types of defects 
and defect complexes, local strain gradients, or substrate effects that are dependent on the 
specific exposure and annealing combinations.  
34 
  
2.3. Conclusions 
 In this work, we engineered single photon sources in pristine hBN by combining FIB 
patterning and carbon diffusion via CVD of nanocrystalline graphite as a carbon source. By 
exploring a portion of the rich parameter space of this technique, we found multiple 
combinations of FIB exposure dose and diffusion time that successfully created SPEs and 
determined that crystals < 20 nm, FIB exposure for 5.0 ms, and carbon annealing for 5 hours lead 
to the highest yield of SPEs, with other combinations providing a means to produce SPEs with 
different average purity and ZPL range. Further work to characterize the influence of these 
parameters on material properties such as sputtering damage and vacancy distributions could 
help identify the classes of defects created with this technique, thereby leading to greater control 
over SPE creation and photophysical properties. This adjustable fabrication technique provides 
an avenue to tunable SPE properties that can be applied to a wider range of integrated quantum 
photonics applications. 
2.4. Methods 
Sample Fabrication 
 All experiments: The samples studied were prepared on Si wafers with 300 nm of thermal 
wet oxide from University Wafer. Wafers were diced into 7 mm by 7 mm chips and patterned 
with a number array via laser writing and reactive ion etching for labeling purposes only. After 
the labeling process, chips were cleaned in acetone sonication bath for 5 min and IPA sonication 
bath for 5 min. Cleaned chips were exposed to 300 mTorr O2 plasma at 300 W RF power in a 
March plasma cleaner for 1 min to create a hydrophilic surface for material transfer adhesion. 
HPHT hBN crystals were prepared for exfoliation transfer on thermal release tape. The thermal 
release tape with hBN crystals was placed on the chip surface immediately after oxygen plasma 
35 
  
exposure. The chip-hBN-tape stack was placed on a 120°C hot plate for 5 minutes to release the 
tape from the chip-hBN stack. The sample was cleaned of any tape residue in a tube furnace 
open to air at 500°C for 3 hours with no gas flow. 
 FIB only, FIB+CVD experiments: FIB exposure to create vacancy and edge defects 
occurred in a FEI Helios Dual Beam with gallium ion source. Prior to sample loading, the 
chamber was plasma cleaned with air for 5 min. The FIB was set to acceleration voltage of 20 
keV and beam current of 59 pA. The circular pattern exposure mask was drawn directly in the 
FIB window of the execution software interface with desired circle diameter of either 250 nm or 
500 nm depending on the experiment. Once the pattern was drawn in the software, an array of 
the circle patterns with 1 μm pitch was made. The exposure time was set in the pattern settings of 
the execution software to 1.3 ms, 2.5 ms, 3.8 ms, or 5.0 ms depending on the experiment. When 
the desired hBN flake was located via SEM, the FIB pattern was executed, exposing only the 
regions of the material as determined in the circular array pattern exposure mask.  
 CVD only, FIB+CVD experiments: Samples immediately underwent the CVD process to 
introduce carbon into the crystal environment. Samples were loaded into the tube furnace 
connected to vacuum system and the system was evacuated to a base pressure of ~30 mTorr. 
When the base pressure was reached, the C2H4, Ar, and H2 gas lines were purged. Once the base 
pressure was reached again, the furnace temperature setpoint was ramped to 1000°C over 30 
minutes under 50 sccm Ar. After 40 minutes, the temperature was held at 1000°C for 20 minutes, 
75 minutes, 150 minutes, or 300 minutes depending on the experiment. During this time, the gas 
flow was 35 sccm ethylene, 2 sccm H2, and 5 sccm Ar at 1 Torr. At the end of the designated 
time, the furnace temperature setpoint was set to 0°C and the sample was allowed to gradually 
cool under 50 sccm Ar. When the temperature readout reached 400°C, the furnace lid was 
36 
  
opened to speed up the rest of the cooling process. The CVD process was ended when the 
temperature readout reaches 200°C. The vacuum system was vented, and samples were removed.  
 FIB only, CVD only, and FIB+CVD experiments: To remove residual carbon on the 
surface, the samples underwent a cleaning anneal at 850°C for 30 min under 5 sccm O2.  
Characterization 
 Atomic Force Microscopy: After material transfer and 500°C cleaning anneal, the hBN 
flake thickness was characterized by taking a line cut measurement of a few micron square area 
AFM image taken via Bruker Dimension FastScan AFM in Tapping Mode.  
 Optical Characterization: Confocal PL images were taken to identify spots of 
photoluminescence in a home-built confocal microscope with 532 nm solid state laser 
(OptoEngine) excitation, 100× 0.7 NA objective, dichroic beam splitter, and 550 nm long pass 
filter. Emission was coupled to a ThorLabs fiber optic 50:50 beam splitter and photon counts for 
each path were detected by avalanche photodiodes (Micro Photon Devices). PL images show the 
signal from one APD. Time correlated photon counting was performed in a Hanbury-Brown and 
Twiss experiment using the signal from both APDs sent to a PicoQuant TimeHarp 260.  
PL spectra of emitting sites and Raman of material were collected in a commercial WITec 
alpha300 spectrometer equipped with a Peltier-cooled Andor iDus CCD.  
 Statistical Analysis: Characterized regions were divided into 4 𝜇𝜇m2 segments, which 
corresponds to 4 FIB patterned locations for patterned samples because we patterned exactly one 
circle per 1 𝜇𝜇m2. To calculate mean SPE (QE) density or yield, we recorded the number of 
emission sites with 𝑔𝑔2(0) < 0.5 (𝑔𝑔2(0) < 1) in each 4 𝜇𝜇m2 region and averaged all regions 
processed with a given combination of experimental input factors. To calculate the mean PL map 
intensity, we converted all confocal image points from kcps to cts by the detection bin size of 0.1 
37 
  
s and plotted in grayscale where 0 corresponded to the minimum counts value and 255 
corresponded to the maximum counts value for the given PL map. We recorded the mean pixel 
height in ImageJ for each region, converted it to the mean PL intensity in counts, and averaged 
all regions processed with a given combination of experimental input factors. To elucidate which 
input factors influence successful SPE creation alone or combined with other input factors and 
obtain significance values, we performed a Full Model Effects Test on the number of SPEs 
counted, number of QEs counted, and average PL in each 4 𝜇𝜇m2 region. 
 All statistical tests were performed in SAS JMP 15.1 statistical analysis software. Factor 
Screening (4-factor 2-level) Full Model Effects Test and all ANOVA tests (test specified in main 
text i.e. Student’s t-test, Tukey-Kramer HSD, or Dunnett test) were performed with significance 
level 𝛼𝛼 = 0.05. Prior to all statistical tests for SPE Yield, QE Yield, and Region PL Counts, we 
performed a variance stabilizing transform for Poissonian data: 
𝑥𝑥� = √𝑥𝑥 (2.3) 
2.5. Supplementary information 
Control experiments 
 To test if both FIB milling and carbon annealing are necessary to create localized quantum 
emitters, we compared the combined process to the following 2 treatments: (1) FIB milled only 
without carbon anneal and (2) carbon anneal only with no FIB milling, as well as to a pristine 
hBN control [Figure 2.7]. We used a 5-hour carbon deposition and anneal time for both the 
CVD only process and the combined process. Our CVD process was performed at 1000°C due to 
the temperature dependent diffusion constant of carbon in semiconductor and insulating 
materials52. Recent works have also shown that this temperature is optimal for carbon-based 
CVD for emitter creation in hBN41. hBN flakes characterized for CVD only experiments were on 
38 
  
Figure 2.7. Process flow diagrams for control experiments (i) pristine hBN control, (ii) 
carbon annealing (CVD) only control, (iii) FIB only control. All samples underwent the same 
transfer and cleaning process. 
the same chips as those characterized for FIB+CVD experiments so may have been exposed to 
residual charging effects, but were not directly patterned via FIB exposure nor imaged via SEM. 
Samples processed with FIB milling that skip the carbon anneal still undergo the 30 min O2 
anneal at 850°C even though there is no need to remove carbon from the surface because we 
found this to be required to activate emitters in our previous work.  
  Table 2.1 shows the average SPE yield per μm2 for each control experiment, including a 
pristine hBN control. Carbon annealing alone produced yields comparable to those found by Lyu 
et. al. using methane feedstock gas in a similar carbon annealing process41. In Table 2.2, we 
report p-values for Student’s t and Tukey-Kramer HSD Means Comparison Tests for the 
optimized combined process parameters with each and all control experiments. 
Table 2.1. Control Experiments Results for average SPE yield. 1 SPE per μm2 on pristine and 
CVD only controls is analogous to 100% yield on FIB only control and combined process 
samples because we patterned 1 site per μm2. 
Experiment Studied region SPE Yield (Mean±Std Err) 
Pristine hBN Control 192 μm2 0 
CVD Only Control 240 μm2  0.0125±0.0071 /μm2 
FIB Only Control 1200 sites 0.0008±0.0024 /μm2 (0.08 %) 
 
 
39 
  
 
Table 2.2. Control Experiments – SPE Yield Mean Comparison Tests 𝛼𝛼 = 0.05 significance 
level, unbalanced, unequal variance 
Experiment 1 – Difference Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Experiment 2 Difference Each Pair) HSD All Pairs) 
Combined Process – 0.3973 0.0410 < 0.0001 < 0.0001 
Pristine hBN Control 
Combined Process – FIB 0.3940 0.0313 < 0.0001 < 0.0001 
Only Control 
Combined Process – 0.3473 0.0389 < 0.0001 < 0.0001 
CVD Only Control 
CVD Only Control – 0.0500 0.0389 0.1997 0.5734 
Pristine hBN Control 
CVD Only Control – FIB 0.0467 0.0284 0.1014 0.3565 
Only Control 
FIB Only Control – 0.0033 0.0313 0.9151 0.9996 
Pristine hBN Control 
 
Raman characterization of graphite diffusion source deposition 
 We confirmed the introduction of carbon via Raman spectroscopy with 532 nm excitation 
[Figure 2.8]. The observed nanocrystalline graphite D and G Raman peaks demonstrates that 
during the carbon annealing process, ethylene successfully decomposed, depositing carbon on 
the surface53. Therefore, we assert that elemental carbon was present during the anneal process. 
The visibility of the hBN Raman peak54 following the full fabrication process demonstrates that 
the oxygen anneal removed the deposited graphite and confirms that the process did not change 
the overall hBN crystal composition or structure.   
Pattern localization 
On samples patterned using the combined process, we observed spots separated by an 
average pitch of 0.934±0.013 μm, found via average centroid calculation of a background 
subtracted PL scan [Figure 2.9]. This value is slightly out of range of the 1 μm input FIB 
patterning pitch. This discrepancy could be due to FIB stage calibration, sample tilt during 
milling, confocal fast scanning mirror calibration, or sample tilt during confocal data collection. 
However, these scans show qualitative agreement to the periodicity of the pattern, showing  
40 
  
Figure 2.8. Raman characterization of diffusion source deposition Peaks showing hBN peak 
(1370 rel. cm-1) for the cleaned pristine flake (black line), nanocrystalline-graphite D (1360 rel. 
cm-1) and G (1610 rel. cm-1) peaks after carbon annealing (red line), and hBN peak again after 
removing NC-Gr in the oxygen anneal 
Figure 2.9. SEM and PL centroid analysis to localize patterned SPEs (a) SEM image of FIB 
patterned 500 nm diameter holes in hBN flakes, separated by 1 µm pitch. (b) Grayscale confocal 
PL scan image showing fluorescent emission at FIB patterned locations after the combined FIB 
and carbon annealing process. Inset: PL image was background subtracted and ImageJ software 
found all “particles”, outlined the closed circles, and recorded the centroid position of each spot. 
The “x” of each spot marks the centroid position and the dashed line shows an example 
displacement measurement to give emission location pitch. Scale bar: 2 µm 
 
41 
  
 
𝑛𝑛 × 𝑛𝑛 emission sites when we mill an 𝑛𝑛 × 𝑛𝑛 pattern. Therefore, we assert this technique as a 
method to produce localized, patterned emission sites. 
Factor screening experiments 
 Parameters: We chose four experimental factors that are straightforward to control and 
chose qualitative “high” and “low” regimes justified by past work or observations prior to SPE 
characterization. Demonstrations of these parameter values are shown in Figure 2.10. 
Figure 2.10. Fabrication parameter factor screening details (a) Optical microscope image 
showing hBN flakes with different optical contrasts due to thickness dependence of thin film 
interference (High: labeled “thick” during optical inspection with mean ± s.d. = 28.47 ± 5.42  
nm, Low: labeled “thin” with mean ± s.d = 14.59 ± 4.81 nm). Scale bar: 30 µm), (b) SEM 
images showing differences in FIB pattern size (high: 500 nm, low: 250 nm), and FIB exposure 
time (high: 3.8 ms, low: 2.5 ms). Scale bars: 2µm, (c) CVD process parameters as a function of 
time for carbon annealing step to facilitate the deposition and diffusion of carbon impurities 
(high: 300 min, low: 20 min). 
42 
  
The hBN thickness is a natural parameter to vary because our straightforward transfer 
method produces a wide range of hBN flake thicknesses. Furthermore, previous studies cite 
material thickness as an influential factor in QE stability55. We placed hBN flakes into “high” or 
“low” regime for our factor screening experiment based on qualitative optical contrast .Then, we 
measured the thickness of individual hBN flakes via Atomic Force Microscopy (AFM), finding 
that flakes placed in the “high” regime were > 20 nm and those placed in the “low” regime were 
< 20 nm. 
 We chose to vary the FIB pattern size to test the lower limit of our QE localizing abilities 
with this method. We set this parameter in the FIB pattern software by varying the pattern 
diameter. A 500 nm pattern size was chosen as our “high” based on [NL paper] and 250 nm 
pattern size was chosen as our “low” based on the expected resolution limits of the FIB 
patterning and confocal microscope characterization systems.  
 To control the amount of damage done to the material in the form of vacancies, Ga+ 
implantation, and sputtering to create milled edges, we chose to vary the FIB exposure time, 
thereby controlling the ion fluence for a constant energy and current. We defaulted to the FIB 
energy and current determined as optimal from our previous work23, detailed in Methods. 
Furthermore, recent work utilizing ion implantation to create SPEs in hBN showed that ion 
fluence has a larger influence on the SPE density than ion energy56. We chose “high” and “low” 
values of 3.8 ms and 2.5 ms respectively, because these exposure times produced distinct levels 
of material damage during our preliminary experiments.  
 CVD involves many parameters that could influence the creation of SPEs, we chose to vary 
the carbon annealing time in order to determine if a steady supply of carbon introduced for a 
long time is necessary for the creation of QEs in hBN or if a steady supply of carbon introduced 
43 
  
for a short time is sufficient. Furthermore, the carbon annealing time controls both the 
concentration of elemental carbon and the diffusion time for all defects and therefore contains 
valuable information within one parameter. We expected the carbon annealing time to 
independently and strongly influence SPE yield because the probability of forming carbon-based 
defects depends explicitly on the concentration of impurity atoms and therefore implicitly on the 
annealing time 𝑡𝑡 by: 
𝑥𝑥
𝐶𝐶(𝑥𝑥, 𝑡𝑡) =  𝐶𝐶𝑠𝑠𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑒 � � (2.3) 2√𝐷𝐷𝑡𝑡
where 𝑥𝑥 is the depth,  𝐶𝐶𝑠𝑠 is the surface concentration of carbon impurity atoms, which is 
dependent on the deposition of NC graphite we observed via Raman, and 𝐷𝐷 is the diffusion 
constant. We kept all other CVD parameters constant for experimental simplicity, but we expect 
that further study will elucidate the role of these parameters such as chosen carbon feedstock, gas 
ratios, and process pressure. We chose a “high” value of 5 hr and a “low” value of 20 min. 
 Data: The mean SPE Yield for all 16 factor screening combinations are shown in Table 2.3. 
The table also lists the number of sites across multiple samples patterned with each combination. 
Additionally, we performed Dunnett Means Comparison Test with pristine hBN control and 
report the p-value for each combination to determine if that particular combination successfully 
creates SPEs with higher yield than what would be found in unprocessed material. Furthermore, 
we report the results of Student’s t and Tukey-Kramer HSD Mean SPE Yield Comparison Tests 
for pairs that include the thin, 500 nm, 3.8 ms, 5 hr parameter combination in Table 2.4. 
 
 
 
 
44 
  
Table 2.3. Factor Screening Experiments Results for average SPE yield. 
Experiment (thickness, pattern Studied region SPE Yield (Mean±Std p-value (Dunnett with 
size, exposure, anneal time) Err) Pristine hBN Control) 
Thin, 500 nm, 3.8 ms, 5 hr 208 sites 0.0673±0.0058 < 0.0001 
Thin, 250 nm, 3.8 ms, 5 hr 48 sites 0.0208±0.0120 0.8784 
Thin, 500 nm, 2.5 ms, 20 min 100 sites 0.0200±0.0083 0.6300 
Thin, 250 nm, 2.5 ms, 20 min 204 sites 0.0098±0.0058 0.9856 
Thin, 500 nm, 3.8 ms, 20 min 104 sites 0.0096±0.0082 0.9983 
Thick, 250 nm, 3.8 ms, 20 min 120 sites 0.0083±0.0076 0.9994 
Thin, 250 nm, 3.8 ms, 20 min 132 sites 0.0076±0.0073 0.9997 
Thick, 500 nm, 3.8 ms, 20 min 296 sites 0.0068±0.0068 0.9993 
Thin, 500 nm, 2.5 ms, 5 hr 196 sites 0 1 
Thin, 250 nm, 2.5 ms, 5 hr 32 sites 0 1 
Thick, 500 nm, 2.5 ms, 5 hr 84 sites 0 1 
Thick, 250 nm, 2.5 ms, 5 hr 92 sites 0 1 
Thick, 500 nm, 2.5 ms, 20 min 124 sites 0 1 
Thick, 250 nm, 2.5 ms, 20 min 60 sites 0 1 
Thick, 500 nm, 3.8 ms, 5 hr 144 sites 0 1 
Thick, 250 nm, 3.8 ms, 5 hr 148 sites 0 1 
 
Table 2.4. Thin, 500 nm, 3.8 ms, 5 hr – Mean SPE Yield Comparison Tests 𝛼𝛼 = 0.05 
significance level, unbalanced, unequal variance. Experiment labels given according to key in 
Figure 2.3(b) 
Experiment 1 – Experiment 2 Difference Standard Error p-value (Student’s t p-value (Tukey-Kramer 
(thickness, pattern size, Difference Each Pair) HSD All Pairs) 
exposure, anneal time) 
-+++ – ---+ 0.2692 0.0717 0.0002 0.0197 
-+++ – -+-+ 0.2692 0.0375 < 0.0001 < 0.0001 
-+++ – +--- 0.2692 0.0553 < 0.0001 0.0002 
-+++ – +--+ 0.2692 0.0472 < 0.0001 < 0.0001 
-+++ – +-++ 0.2692 0.0406 < 0.0001 < 0.0001 
-+++ – ++-- 0.2692 0.0428 < 0.0001 < 0.0001 
-+++ – ++-+ 0.2692 0.0488 < 0.0001 < 0.0001 
-+++ – ++++ 0.2692 0.0409 < 0.0001 < 0.0001 
-+++ – +++- 0.2422 0.0341 < 0.0001 < 0.0001 
-+++ – --+- 0.2389 0.0420 < 0.0001 < 0.0001 
-+++ – +-+- 0.2358 0.0432 < 0.0001 < 0.0001 
-+++ – -++- 0.2307 0.0453 < 0.0001 < 0.0001 
-+++ – ---- 0.2300 0.0372 < 0.0001 < 0.0001 
-+++ – -+-- 0.1892 0.0459 < 0.0001 0.0050 
-+++ – --++ 0.1858 0.0604 0.0022 0.1565 
 
 
 
 
45 
  
 Full Model Effects Test: We report results from the SPE Yield Full Model Effects test in 
Table 2.5, including both the effect size and p-value, listed in order by p-value. Interaction plots 
for all single factor, 2-factor, and 3-factor effects for SPE Yield are displayed in Figures 2.11-
2.13. 
Table 2.5. Factor Screening Full Model Effects Test for Variance Transformed Number of 
SPEs per 4 patterned sites 𝛼𝛼 = 0.05 significance level, unbalanced, unequal variance. 
Input Factor or Interaction Effect Size (F) p-value 
hBN thickness 8.24 𝟎𝟎.𝟎𝟎𝟎𝟎𝟎𝟎 
hBN thickness × FIB exposure time × Carbon annealing time 7.89 𝟎𝟎.𝟎𝟎𝟎𝟎𝟎𝟎 
FIB exposure 4.85 𝟎𝟎.𝟎𝟎𝟎𝟎𝟎𝟎 
FIB exposure × Carbon annealing time 4.17 𝟎𝟎.𝟎𝟎𝟎𝟎𝟎𝟎 
hBN thickness × FIB exposure time 2.06 0.152 
hBN thickness × FIB pattern size 2.05 0.153 
FIB pattern size 1.95 0.164 
hBN thickness × Carbon annealing time 1.89 0.170 
FIB exposure time × FIB pattern size × Carbon annealing time 1.78 0.183 
hBN thickness × FIB exposure time × FIB pattern size × Carbon annealing 1.69 0.195 
time 
hBN thickness × FIB exposure time × FIB pattern size 0.89 0.347 
FIB exposure × FIB pattern size 0.82 0.366 
FIB pattern size × Carbon annealing time 0.71 0.400 
hBN thickness × FIB pattern size × Carbon annealing time 0.65 0.421 
Carbon annealing time 0.37 𝟎𝟎.𝟎𝟎𝟎𝟎𝟎𝟎 
 
  
46 
  
Figure 2.11. Plots for single factor effects on SPE yield for 2 level factor screening showing least 
squares means by level and parameter combination averages for each factor (a) hBN thickness (b) 
FIB exposure time, (c) FIB pattern size, (d) carbon anneal time. Error bars are standard error on 
mean. 
 
  
47 
  
 
Figure 2.12. Plots for two-factor interaction effects on SPE yield for 2 level factor screening 
showing least squares means by level combination for each interaction (a) FIB exposure and 
carbon anneal time (b) hBN thickness and FIB exposure (c) hBN thickness and pattern size (d) 
hBN thickness and carbon anneal (e) pattern size and FIB exposure (f) pattern size and carbon 
anneal time. Points are least squares mean of parameter combinations that include the labeled 
interaction. Error bars are standard error on mean. 
 
  
48 
  
 
Figure 2.13. Plots for three-factor interaction effects on SPE yield for 2 level factor screening 
showing least squares means by level combination for each interaction (a) hBN thickness, FIB 
exposure, and carbon anneal time (b) pattern size, FIB exposure, and carbon anneal time (c) 
pattern size, hBN thickness, and FIB exposure (d) pattern size, hBN thickness, and carbon anneal 
time.  Points are least squares mean of parameter combinations that include the labeled interaction. 
Error bars are standard error on mean. 
Optimization 
 Parameters: Carbon annealing and FIB exposure levels and processes are shown in Figure 
2.14. Based on the factor screening results, we only exposed flakes with thickness in the range of 
11-18 nm and only patterned circles of 500 nm diameter. For each substrate, we found four 
flakes of suitable area in the prescribed thickness range determined as optimal during Factor 
Screening. On each flake, we patterned one 7×7 array of circles at each dose. We randomized 
the order of the four doses patterned on each flake to minimize any pattern order-dependent 
49 
  
variance in measured response variables. The lowest value of 1.3 ms was the lowest time 
possible given the fixed pattern size, current, and dwell time specified for the Helios Dual Beam 
FIB system used for this work. Exposure times longer than 5 ms in past experiments milled 
through hBN within this thickness range completely and began milling SiO2. Given the constant 
pattern area for this set of fabrication, we can also refer to the FIB exposure time parameter by 
the ion fluence [Table 2.6] calculated by: 
𝐼𝐼 × 𝑡𝑡
𝐹𝐹𝐹𝐹𝐹𝐹𝑒𝑒𝑛𝑛𝑒𝑒𝑒𝑒 =  𝑏𝑏𝑏𝑏𝑏𝑏𝑏𝑏 𝑏𝑏𝑒𝑒𝑒𝑒𝑒𝑒𝑠𝑠𝑏𝑏𝑒𝑒 (2.4)
𝑒𝑒+ × 𝐴𝐴  𝑒𝑒𝑏𝑏𝑡𝑡𝑡𝑡𝑏𝑏𝑝𝑝𝑝𝑝
2
Where 𝐼𝐼𝑏𝑏𝑏𝑏𝑏𝑏𝑏𝑏 = 59 𝑝𝑝𝐴𝐴, 𝑒𝑒+ = 1.6𝐸𝐸 − 7 𝑝𝑝𝐶𝐶, and 𝐴𝐴
500 𝑝𝑝𝑏𝑏
𝑒𝑒𝑏𝑏𝑡𝑡𝑡𝑡𝑏𝑏𝑝𝑝𝑝𝑝 = 𝜋𝜋 � �  2
Figure 2.14. FIB exposure and carbon annealing interaction optimization parameters, (a) 
Optical microscope image of hBN flake marked with array locations for each FIB level. Pattern 
order was random to randomize the changes to the FIB exposure time (b) SEM image and 
corresponding confocal PL scan of patterned flake shown in (a). Scale bars for (a) and (b): 10 
µm. (c) CVD process showing the 4 carbon anneal time levels 
 
Table 2.6. Ion Fluences for Optimization FIB Exposure Levels calculated according to 
Equation 2.4 
FIB Exposure Level Exposure Time (ms) Ion Fluence (cm-2) 
1 1.25 2.35×1014 
2 2.5 4.70×1014 
3 3.75 7.04×1014 
4 5.0 9.39×1014 
 
50 
  
 Data: The mean SPE Yield, QE Yield, and Region PL for all 18 optimization combinations 
are shown in Table 2.7. The table also lists the number of sites across multiple samples patterned 
with each combination. Additionally, we performed Dunnett Means Comparison Test with 
pristine hBN control and report the p-value for each combination to determine if that particular 
combination successfully creates SPEs and QEs (PL) with higher yield (intensity) than what 
would be found in unprocessed material. Furthermore, we report the results of Student’s t and 
Tukey-Kramer HSD Mean SPE Yield Comparison Tests for pairs that include parameter 
combinations where the SPE yield was greater than zero in Table 2.8. We also report the results 
of Student’s t and Tukey-Kramer HSD Mean Region PL Grouped by Anneal Time Comparison 
Tests for pairs that include parameter combinations in Table 2.9. 
 
Table 2.7. Optimization Experiments Results for average SPE yield, average QE yield, and 
Region PL intensity. 𝛼𝛼 = 0.05 significance level 
Experiment Number SPE Yield p-value QE Yield p-value Region PL p-value 
(FIB exposure, of (Mean±St (Dunnett (Mean±St (Dunnett [cts] (Dunnett 
carbon anneal patterned d Err) with d Err) with (Mean±St with 
time) sites Pristine Pristine d Err) Pristine 
hBN hBN hBN 
Control) Control) Control) 
5.0 ms, 5 hr 192 0.1146±0. < 0.0001 0.1875±0. < 0.0001 46.03±2.5 0.9716 
0094 0144 6 
3.8 ms, 2.5 hr 184 0.0598±0. < 0.0001 0.0978±0. < 0.0001 86.15±10. 0.2305 
0096 0147 30 
3.8 ms, 5 hr 192 0.0469±0. 0.0014 0.1094±0. < 0.0001 186.37±25 0.0004 
0094 0144 .08 
2.5 ms, 20 196 0.0305±0. 0.1318 0.1453±0. < 0.0001 239.80±48 0.0001 
min 0093 0142 .36 
5.0 ms, 12 hr 388 0.0155±0. 0.7758 0.0799±0. < 0.0001 108.63±20 0.0525 
0054 0083 .79 
2.5 ms, 1.25 hr 184 0.0109±0. 0.9877 0.0380±0. 0.2691 393.22±46 < 0.0001 
0096 0147 .33 
5.0 ms, 20 min 176 0.0057±0. 1.0000 0.0057±0. 1.0000 28.62±0.7 1.0000 
0098 0150 7 
3.8 ms, 1.25 hr 184 0.0053±0. 1.0000 0.0053±0. 1.0000 36.16±0.6 1.0000 
0095 0145 7 
3.8 ms, 20 min 132 0.0050±0. 1.0000 0.0147±0. 0.911 39.06±1.6 0.9998 
0091 0140 7 
2.5 ms, 2.5 hr 192 0 1.0000 0.0052±0. 1.0000 769.81±57 < 0.0001 
0144 .56 
1.3 ms, 20 min 192 0 1.0000 0 1.0000 468.32±21 < 0.0001 
.62 
51 
  
 
Table 2.7. (continued) 
Experiment Number SPE Yield p-value QE Yield p-value Region PL p-value 
(FIB exposure, of (Mean±St (Dunnett (Mean±St (Dunnett [cts] (Dunnett 
carbon anneal patterned d Err) with d Err) with (Mean±St with 
time) sites Pristine Pristine d Err) Pristine 
hBN hBN hBN 
Control) Control) Control) 
1.3 ms, 1.25 hr 192 0 1.0000 0 1.0000 298.54±26 < 0.0001 
.19 
1.3 ms, 2.5 hr 184 0 1.0000 0 1.0000 385.09±42 < 0.0001 
.60 
1.3 ms, 5 hr 192 0 1.0000 0 1.0000 573.44±26 < 0.0001 
.55 
2.5 ms, 5 hr 192 0 1.0000 0 1.0000 1097.49±5 < 0.0001 
5.80 
3.8 ms, 12 hr 376 0 1.0000 0 1.0000 233.19±21 < 0.0001 
.12 
5.0 ms, 1.25 hr 204 0 1.0000 0 1.0000 29.35±0.4 1.0000 
7 
5.0 ms, 2.5 hr 192 0 1.0000 0 1.0000 31.46±0.5 1.0000 
7 
 
Table 2.8. Optimization Experiments – Mean SPE Yield Comparison Tests for pairs where 
SPE yield > 0. 𝛼𝛼 = 0.05 significance level, unbalanced, unequal variance 
Experiment 1 – SPE Difference SPE Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Experiment 2 Difference Each Pair) HSD All Pairs) 
4:4 – 3:1 0.3777 0.0438 < 0.0001 < 0.0001 
4:4 – 3:2 0.3760 0.0447 < 0.0001 < 0.0001 
4:4 – 4:1 0.3745 0.0455 < 0.0001 < 0.0001 
4:4 – 2:2 0.3538 0.0449 < 0.0001 < 0.0001 
4:4 – 4:5 0.3414 0.0384 < 0.0001 < 0.0001 
4:4 – 2:1 0.2868 0.0442 < 0.0001 < 0.0001 
4:4 – 3:4 0.2220 0.0445 < 0.0001 0.0001 
3:3 – 3:1 0.1940 0.0443 < 0.0001 0.002 
3:3 – 3:2 0.1923 0.0452 < 0.0001 0.0034 
3:3 – 4:1 0.1909 0.0459 < 0.0001 0.0051 
4:4 – 3:3 0.1836 0.0449 < 0.0001 0.0068 
3:3 – 2:2 0.1701 0.0454 0.0002 0.0242 
3:3 – 4:5 0.1578 0.0390 < 0.0001 0.0079 
3:4 – 3:1 0.1556 0.0438 0.0004 0.0464 
3:4 – 3:2 0.1540 0.0447 0.0006 0.0653 
3:4 – 4:1 0.1525 0.0455 0.0008 0.0856 
3:4 – 2:2 0.1318 0.0449 0.0035 0.2559 
3:4 – 4:5 0.1194 0.0384 0.002 0.1687 
3:3 – 2:1 0.1031 0.0447 0.0214 0.7067 
2:1 – 3:1 0.0908 0.0436 0.0375 0.8464 
2:1 – 3:2 0.0892 0.0445 0.0454 0.8851 
2:1 – 4:1 0.0877 0.0453 0.0529 0.9121 
2:1 – 2:2 0.0670 0.0447 0.1347 0.9933 
3:4 – 2:1 0.0648 0.0442 0.1437 0.9948 
52 
  
 
Table 2.8 (continued) 
Experiment 1 – SPE Difference SPE Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Experiment 2 Difference Each Pair) HSD All Pairs) 
2:1 – 4:5 0.0546 0.0382 0.1529 0.9961 
3:3 – 3:4 0.0383 0.0449 0.3941 1.0000 
4:5 – 3:1 0.0362 0.0377 0.3373 1.0000 
4:5 – 3:2 0.0345 0.0387 0.3731 1.0000 
4:5 – 4:1 0.0330 0.0396 0.404 1.0000 
2:2 – 3:1 0.0238 0.0443 0.5905 1.0000 
2:2 – 3:2 0.0222 0.0452 0.6236 1.0000 
2:2 – 4:1 0.0207 0.0459 0.6519 1.0000 
4:5 – 2:2 0.0123 0.0390 0.752 1.0000 
4:1 – 3:1 0.0031 0.0448 0.9446 1.0000 
3:2 – 3:1 0.0016 0.0440 0.9698 1.0000 
4:1 – 3:2 0.0014 0.0457 0.9747 1.0000 
 
Table 2.9. Optimization Experiments – Mean Region PL Grouped by Anneal Time 
Comparison Tests 𝛼𝛼 = 0.05 significance level, unbalanced, unequal variance 
Experiment 1 – PL Difference PL Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Experiment 2 Difference Each Pair) HSD All Pairs) 
2:1 – 1:1 8.66 1.04 < 0.0001 < 0.0001 
2:1 – 3:1 6.37 1.02 < 0.0001 < 0.0001 
2:1 – 4:1 7.24 1.06 < 0.0001 < 0.0001 
1:1 – 3:1 15.03 1.03 < 0.0001 < 0.0001 
1:1 – 4:1 15.90 1.06 < 0.0001 < 0.0001 
3:1 – 4:1 0.87 1.05 0.4085 0.8409 
2:2 – 1:2 1.32 0.999 0.1877 0.5500 
2:2 – 3:2 11.93 1.00 < 0.0001 < 0.0001 
2:2 – 4:2 12.52 0.98 < 0.0001 < 0.0001 
1:2 – 3:2 10.61 0.99 < 0.0001 < 0.0001 
1:2 – 4:2 11.20 0.97 < 0.0001 < 0.0001 
3:2 – 4:2 0.59 0.98 0.5455 0.9302 
2:3 – 1:3 8.40 1.08 < 0.0001 < 0.0001 
2:3 – 3:3 18.21 1.08 < 0.0001 < 0.0001 
2:3 – 4:3 21.23 1.07 < 0.0001 < 0.0001 
1:3 – 3:3 9.81 1.09 < 0.0001 < 0.0001 
1:3 – 4:3 12.82 1.08 < 0.0001 < 0.0001 
3:3 – 4:3 3.02 1.08 0.0059 0.0297 
2:4 – 1:4 8.93 0.98 < 0.0001 < 0.0001 
2:4 – 3:4 20.34 0.98 < 0.0001 < 0.0001 
2:4 – 4:4 25.90 0.98 < 0.0001 < 0.0001 
1:4 – 3:4 11.41 0.98 < 0.0001 < 0.0001 
1:4 – 4:4 16.97 0.98 < 0.0001 < 0.0001 
3:4 – 4:4 5.56 0.98 < 0.0001 < 0.0001 
 
 
53 
  
𝒈𝒈𝟎𝟎(𝟎𝟎) Analysis 
 We report the results of Student’s t and Tukey-Kramer HSD Mean 𝑔𝑔2(0) Comparison Tests 
for the optimization parameter combinations with the 5 highest SPE/QE yields in Table 2.10. 
 
Table 2.10. Top 5 Optimization Experiments – Mean 𝒈𝒈𝟎𝟎(𝟎𝟎) Comparison Tests 𝛼𝛼 = 0.05 
significance level, unbalanced, unequal variance 
Experiment 1 – 𝑔𝑔2(0) Difference 𝑔𝑔2(0) Standard p-value (Student’s t p-value (Tukey-Kramer 
Experiment 2 Error Difference Each Pair) HSD All Pairs) 
4:5 – 3:3 0.262 0.047 < 0.0001 < 0.0001 
4:5 – 4:4 0.221 0.040 < 0.0001 < 0.0001 
3:4 – 3:3 0.194 0.050 0.0002 0.0015 
2:1 – 3:3 0.170 0.050 0.0009 0.0076 
3:4 – 4:4 0.154 0.043 0.0005 0.0048 
2:1 – 4:4 0.129 0.043 0.0032 0.0261 
4:5 – 2:1 0.092 0.043 0.0325 0.2014 
4:5 – 3:4 0.067 0.043 0.1215 0.5270 
4:4 – 3:3 0.041 0.047 0.3894 0.9097 
3:4 – 2:1 0.025 0.046 0.5913 0.9832 
 
 Figure 2.15 shows the distribution of fit parameters for the antibunching data collected from 
the QEs characterized for this work. While 𝑔𝑔2(0) values were process dependent, we did not 
find any significant variations in average bunching amplitude, lifetime 𝜏𝜏1, and decoherence time 
𝜏𝜏2 across the parameter combinations that produced the highest emitter yields. 
ZPL characterization 
 We report the results of Student’s t and Tukey-Kramer HSD Mean average ZPL Comparison 
Tests for each FIB exposure time averaged over all anneal times in Table 2.11 and for each 
carbon annealing time averaged over all FIB exposure times in Table 2.12. 
54 
  
Figure 2.15. Second order photon correlation fit parameter distributions (a) Histograms for 
the fit parameter values extracted from fitting Equation 2.1 to photon correlation data for each 
color center characterized for this work where 𝑔𝑔2(0) < 1.0 (b) Example antibunching curves 
from 𝑔𝑔2(𝑡𝑡) with green points showing rebinned data and black lines showing fit (c) Violin plots 
showing distributions of bunching amplitude, 𝜏𝜏1, and 𝜏𝜏2 for the five best QE fabrication 
parameter combinations. 
 
Table 2.11. ZPL by FIB Exposure – Mean Comparison Tests 𝛼𝛼 = 0.05 significance level, 
unbalanced, unequal variance 
Exposure 1 – Exposure 2 Difference Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Difference Each Pair) HSD All Pairs) 
5.0 ms – 2.5 ms 0.064 0.021 0.0022 0.0063 
3.8 ms – 2.5 ms 0.043 0.018 0.0164 0.0430 
5.0 ms  – 3.8 ms 0.021 0.021 0.3121 0.5690 
 
Table 2.12. ZPL by Carbon anneal time – Mean Comparison Tests 𝛼𝛼 = 0.05 significance 
level, unbalanced, unequal variance 
Anneal Time 1 – Anneal Difference Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Time 2 Difference Each Pair) HSD All Pairs) 
2.5 hr – 1.25 hr 0.068 0.029 0.0213 0.0964 
5 hr  – 1.25 hr 0.063 0.022 0.0049 0.0252 
2.5 hr  – 20 min 0.036 0.028 0.1885 0.5509 
20 min – 1.25 hr 0.032 0.024 0.1871 0.5481 
5 hr – 20 min 0.031 0.020 0.1172 0.3954 
2.5 hr – 5 hr 0.006 0.026 0.8247 0.9961 
 
55 
  
 Figure 2.16 demonstrates examples of spectra that were similar to the ZPL shape but had 
𝑔𝑔2(0) > 0.5. From this data, we conclude that at these pattern sites where we did not find SPE, 
we fabricated clusters of defects, with multiple excited at once by a single laser excitation spot. 
Figure 2.16. Example zero phonon line spectra shapes for 𝒈𝒈𝟎𝟎(𝟎𝟎) > 𝟎𝟎.𝟎𝟎 Many of the color 
centers characterized for this work exhibited spectral lineshapes similar to ZPL, with a distinct 
peak and secondary phonon sideband ~160 meV, but with either (a) high background, (b) multiple 
ZPL, or (c)-(d) broad linewidths characteristic of ensemble emission. 
 
2.6. Bridge 
 In this chapter, I have reported a technique we developed to deterministically and repeatably 
engineer carbon-based defects in hexagonal boron nitride using a combination of focused ion 
beam patterning and carbon diffusion. By designing experiments to optimize the yield of single 
photon emitters and characterize photophysical properties including emission energy and 
quantum purity, we found that these outputs depended strongly on the combination of FIB 
56 
  
exposure dose and the length of the carbon diffusion step chosen for fabrication. However, we 
were unable to improve our SPE yields over ~11.5% and our results are inconclusive in 
identifying a particular defect candidate responsible for this emission. To improve our SPE 
yields and potentially reduce variations in defect type produced, we need to better control defect 
creation as well as the material conditions surrounding the desired defects.  
 FIB patterning is a ballistic process known to cause structural damage to materials at the 
beam energy we employ. Because the FIB exposure dose had strong independent effects on 
optical outputs, we therefore set out to investigate how the FIB dose, independent of carbon 
annealing time, impacts the host material and therefore color center creation. We achieve this via 
atomic force microscopy (AFM) height topography measurements and Stopping Range of Ions in 
Matter (SRIM) simulations to quantify damage and correlate FIB damage characteristics to 
optical outputs. The following chapter details our AFM characterization, observed damage 
characteristics, FIB exposure simulations, and correlation analysis, as well as proposes a 
preliminary physical mechanism for how FIB exposure dose influences color center creation that 
is consistent with the observed morphological damage features and optical characteristics. 
 
57 
 
 
CHAPTER III 
 
UNCOVERING THE MORPHOLOGICAL EFFECTS OF HIGH ENERGY GA+  
FOCUSED ION BEAM MILLING ON HBN SINGLE-PHOTON EMITTER  
FABRICATION 
 
This chapter contains previously published co-authored material; it has been reproduced with 
permission from: R. Klaiss, J. Ziegler, D. Miller, K. Zappitelli, K. Watanabe, T. Tanaguchi, B. J. 
Alemán, “Deterministic Quantum Emitter Formation in Hexagonal Boron Nitride via Controlled 
Edge Creation.” J Chem Phys. 2022, 157 (074703), 1-10. Copyright 2022 American Institute of 
Physics. In this work, I contributed to designing the research, performing the research, analyzing 
the data, and writing the paper. 
 
3.1. Introduction 
 The advancement of quantum information technologies (QITs) such as quantum transducers 
and sensors require solid-state qubit platforms that host optically addressable single-photon or 
spin states that couple strongly to external degrees of freedom and can be easily engineered for 
integration into hybrid systems1,4,57,58. Low dimensional materials are desirable solid-state single 
photon emitter (SPE) hosts because they exhibit unique electronic, magnetic, and mechanical 
properties impossible for their bulk counterparts32,35,59–62. Two-dimensional wide bandgap 
insulator hexagonal boron nitride (hBN) hosts bright and stable room temperature SPEs that 
couple strongly to applied electric32,59 and strain31,63 fields, as well as optically detectable 
magnetic resonance (ODMR) with spin state splitting sensitive to external magnetic 
5 8 
  
fields29,38,64,65. Furthermore, hBN functions well as both a dielectric substrate for other 2D 
materials such as graphene62,66,67 and transition metal dichalcogenides (TMDCs)46,68,69 and as a 
monolithic platform for photonic circuits70–72, making it a versatile material to integrate into 
QITs.  
 Fabrication challenges such as patterning hBN SPEs23,43,50,73–77 and direct-writing 
nanostructures in hBN78–80 have been resolved using focused ion beam (FIB) exposure. However, 
the high energies required to pattern emitters at precise locations and mill nanostructures cause 
changes to the morphology of thin hBN crystals that may influence further fabrication steps and 
impact the material’s ability to host SPEs46,80.  Moreover, while SPEs found in hBN after FIB 
exposure have been attributed to either vacancies created deep within the crystal or edges created 
due to milling, it is not understood how changes to surface morphology that occur due to FIB 
exposure influence the probability of SPE creation.   
 To uncover the morphological life-story of hBN crystals to support SPEs while undergoing 
high energy Ga+ FIB exposure, we fabricated samples with a range of exposure doses and 
characterized changes to the surface roughness, milling depth, swelling, and single photon 
emission photophysics. Then, we analyzed the correlation between morphological factors and 
photophysical responses. Furthermore, we performed Stopping Range of Ions in Matter (SRIM) 
calculations81 to investigate how vacancy concentration influences morphology and SPE 
creation. This work provides a necessary framework to elucidate how pristine hBN crystals 
become damaged in controlled ways at the atomic and nano- scales to perform as needed for 
quantum technologies. 
59 
  
3.2. Results and discussion 
To create samples that contain many 
visible hBN crystals, we performed an 
exfoliation transfer of high-pressure high 
temperature (HPHT) grown hBN crystals onto 
Si/SiO2 substrates using heat-release tape and 
annealed samples at 500°C to clean residual 
tape from the surface. Once we identified 
crystals of area ~200 μm2 [Figure 3.1(a)] and 
thickness in the range of 10-17 nm [see 
Figure 3.1. Testing FIB exposure on a single 
Supplementary information Figure 3.15], we hBN flake (a) Optical microscope image of an 
exfoliated HPHT hBN flake studied in this 
performed Ga+ FIB patterning tests. Based on work prior to FIB exposure. (b) Scanning 
electron microscope (SEM) image of an hBN 
our previous work23 that created SPEs at flake after patterning arrays of 500 nm circles 
with Ga+ FIB at the four studied exposure 
precise locations, we patterned arrays of 500 times, from left to right, labeled by FIB Level: 
(3) 3.8 ms FIB Level, (4) 5.0 ms, (1) 1.3 ms, 
nm diameter circles with beam acceleration and (2) 2.5 ms. Scale bar: 5 μm 
voltage 20 kV and beam current 59 pA. The Ga+ FIB used for this work had a minimum 
exposure time of 1.3 ms for the specified pattern size and current, corresponding to the shortest 
loop time the FIB can direct-write the circle pattern. We adjusted the ion fluence by changing the 
exposure time, thereby changing the number of passes performed by the software to pattern the 
circle and tested 1.3 ms (FIB Level 1: 1 pass), 2.5 ms (FIB Level 2: 2 passes), 3.8 ms (FIB Level 
3: 3 passes), and 5.0 ms (FIB Level 4: 4 passes) [Figure 3.1(b)]. We annealed the exposed 
samples in a home-built chemical vapor deposition (CVD) system used for vertically-aligned 
carbon nanotube (VACNT) growth with C2H4:H2:Ar flowing gas mixture at 1000°C for 5 hr in 
60 
  
order to create a carbon-rich environment that is expected to activate hBN quantum 
emitters37,40,41. For more information about the emitter fabrication process, refer to 3.4 Methods.  
 After the full emitter fabrication process, we performed tapping-mode atomic force 
microscopy (AFM) to measure the topography of nanoscale surface features and observed 
qualitatively different morphological features at each exposure level [Figure 3.2].  
Figure 3.2. AFM characterization of FIB exposure levels AFM height scan images for circles 
patterned at each of the studied exposure levels marked with black dashed lines where line cuts 
(below) were taken. AFM height color bar range: [-15 nm, 5 nm]. Line cut: 1 μm 
 From the AFM height data, we calculated the following morphological characteristics of 1 
μm × 1 μm regions containing one exposed circle each: average mill depth ∆��𝑧𝑧��−�, maximum mill 
depth ∆𝑧𝑧− −𝑏𝑏𝑏𝑏𝑒𝑒, mill volume ∆𝑉𝑉 , average swell height ∆��𝑧𝑧��+�, maximum swell height ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒, 
swell volume ∆𝑉𝑉+, and exposed edge surface area 𝑆𝑆𝐴𝐴 [Figure 3.3]. Local unexposed surface 
height 𝑧𝑧0 and RMS roughness 𝑅𝑅𝑅𝑅𝑆𝑆0 for each 1 μm × 1 μm region were calculated from points 
outside of the exposed pattern area and used as the zero height for calculations and as the 
threshold for identifying FIB induced damage, respectively. Milling values ∆��𝑧𝑧��−�,  ∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒, and 
∆𝑉𝑉− were calculated from data points with height below 𝑧𝑧0 − 𝑅𝑅𝑅𝑅𝑆𝑆0 and characterized the depth 
of milled locations and total amount of material removed from the exposed regions. Swelling 
values ∆��𝑧𝑧��+�,  ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒, and ∆𝑉𝑉+ were calculated from data points above 𝑧𝑧0 + 𝑅𝑅𝑅𝑅𝑆𝑆0 and 
61 
  
characterized the height and total amount of 
FIB induced swelling observed in exposed 
regions. Exposed edge surface area 
calculations were performed using all data 
points in the region, with baseline area of 1 
μm2 subtracted. For more details about the 
calculations, refer to 3.4 Methods.  
 We report the average values for all 
Figure 3.3. Morphological characteristics 
calculated morphological characteristics for calculated from AFM height data AFM height 
profile showing 𝑅𝑅𝑅𝑅𝑆𝑆0 threshold for swell and 
each exposure level in Table 3.1 and show mill ∆𝑧𝑧 and ∆𝑉𝑉 calculations from AFM height 
data with corresponding AFM height image 
the average values of ∆𝑧𝑧+  and ∆𝑧𝑧−  showing points above (top) and below (bottom) 𝑏𝑏𝑏𝑏𝑒𝑒 𝑏𝑏𝑏𝑏𝑒𝑒 the calculated roughness threshold included in 
for each exposure level in Figure 3.4(a), swelling (top) and milling (bottom) 
calculations. AFM height image scale bar: 200 
with average values for ∆𝑉𝑉+ and ∆𝑉𝑉− for nm 
each exposure in Figure 3.4(b). Plots showing the average values for all calculated 
morphological characteristics as a function of FIB exposure level can be found in 3.5 
Supplementary information. Based on these calculated values as well as observations made 
from the data visualized as AFM scan plots in Figure 3.2, we found that swelling characteristics 
and milling characteristics evolve differently with increasing levels of FIB exposure.  
  Starting at the lowest FIB exposure level (1.3 ms), exposed regions contained small and 
randomly located areas of milling and swelling within the patterned area with �∆�𝑧𝑧��−�1 = −3.26 ±
0.15 nm and ∆��𝑧𝑧��+�1 = 1.92 ± 0.09 nm. Regions exposed at the second FIB level (2.5 ms) showed 
qualitatively similar milling and swelling characteristics and calculated values were similar to the 
first exposure level, with �∆�𝑧𝑧��−�2 = −2.89 ± 0.11 nm (Table 3.1 Student’s t-test: p = 0.1566) and 
62 
  
 
Table 3.1. Morphological Characteristics Calculated from AFM mean ± standard error 
FIB − + +����−�
Level ∆𝑧𝑧  (nm) ∆𝑧𝑧
−  (nm) ∆𝑉𝑉   ����+� ∆𝑧𝑧 𝑏𝑏𝑏𝑏𝑒𝑒 ∆𝑉𝑉   𝑆𝑆𝐴𝐴  𝑏𝑏𝑏𝑏𝑒𝑒 (105 nm3) ∆𝑧𝑧  (nm) (nm) (105 nm3) (103 nm2) 
1 −3.26±0.15 −8.63±0.33 −3.87±0.26 1.92±0.09 2.92±0.18 1.37±0.15 2.03±0.11 
2 −2.89±0.11 −7.92±0.22 −3.33±0.22 1.81±0.06 2.57±0.11 1.33±0.14 1.77±0.06 
3 −3.25±0.14 −8.19±0.27 −2.53±0.15 2.20±0.08 3.70±0.23 1.89±0.22 1.67±0.07 
4 −6.23±0.29 −13.23±0.32 −9.79±0.67 2.12±0.08 3.96±0.25 2.44±0.23 4.12±0.29 
 
�∆�𝑧𝑧��+�2  = 1.81 ± 0.06 nm (Table 3.1 Student’s t-test: p = 0.3245). The largest change observed 
was an 8% decrease in the magnitude of ∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒, with ∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒1 = −8.63 ± 0.33 nm and 
∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒2 = −7.92 ± 0.22 nm, but the difference between these values was not strong enough to 
draw any conclusions about the evolution of morphological characteristics from the first 
exposure level to the second (Table 3.1 Student’s t-test: p = 0.0837). Therefore, we determine  
that there are no significant morphological differences between regions exposed for 1.3 ms and 
those exposed for 2.5 ms. 
Figure 3.4. Average morphological characteristics for each FIB exposure level (a) Mean 
calculated values for maximum swell height and maximum mill depth at each FIB exposure level 
of all samples included in this work. (b) Mean calculated values for swell volume and mill 
volume at each FIB exposure level of all samples included in this work. The lines connecting 
data points are a guide to the eye. Error bars on show standard error on the mean. 
63 
  
 At the third exposure level (3.8 ms), we began to observe significant changes to 
morphological characteristics. While regions exposed at the third exposure level were milled 
similarly to those at the second with ∆��𝑧𝑧��−�3  = −3.25 ± 0.14 nm (Table 3.1 Student’s t-test ∆��𝑧𝑧��−�: 
p = 0.1742), both the average and maximum swell heights increased from the second to the third  
FIB level (t-test ∆��𝑧𝑧��+�: p = 0.0006, t-test ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒: p < 0.0001), with the mean ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒 increasing 
44% from ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒2  = 2.57 ± 0.11 nm to ∆𝑧𝑧
+
𝑏𝑏𝑏𝑏𝑒𝑒3= 3.70 ± 0.23 nm.  
 Regions exposed at the fourth level (5.0 ms) showed swell heights similar to the values 
calculated at the third (t-test ∆��𝑧𝑧��+�: p = 0.4653, t-test ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒: p = 0.3754). However, ∆𝑉𝑉+ 
increased 23% from ∆𝑉𝑉+3= 1.89 ± 0.22 × 105 nm3 to ∆𝑉𝑉+4= 2.44 ± 0.23 × 105 nm3 (t-test: p 
= 0.0444). The areas of swelling for regions exposed at the fourth level appeared more localized 
to the pattern edges, away from the pattern center whereas swelling for regions exposed at the 
third level was observed within the pattern center. While the swelling was pushed out to the 
pattern edges, we observed deeper milling more localized to the circle pattern at this level, with 
the magnitude of ∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒 increasing 62% from ∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒3= −8.19 ± 0.27 nm to ∆𝑧𝑧
−
𝑏𝑏𝑏𝑏𝑒𝑒4  = 
−13.23 ± 0.32 nm (t-test: p < 0.0001). Furthermore, we observed similar behavior in the 
exposed surface area 𝑆𝑆𝐴𝐴, where 𝑆𝑆𝐴𝐴3= 1.67 ± 0.07 × 103 nm2 increased 147% to 𝑆𝑆𝐴𝐴4= 4.12 ±
0.29 × 103 nm2 (t-test: p < 0.0001) most likely due to the increase in directional milling 
forming connected sidewalls at the fourth exposure level. We conclude that swelling is the 
dominant morphological characteristic starting at the third exposure level and is a precursor to 
the deeper directed milling that began at the fourth level.  
 To determine how the targeted morphological characteristics influenced emitter creation, we 
identified fluorescent locations via confocal microscopy with 532 nm excitation [Figure 3.5(a)] 
and characterized quantum emission via Hanbury Brown and Twiss experiment, fitting 𝑔𝑔2(𝑡𝑡) 
64 
  
data to the three-level system model for antibunching and using the accepted threshold of 
𝑔𝑔2(0) < 0.5 to denote SPE.46 To further confirm the quantum nature of emitters, we measured 
zero-phonon line (ZPL) spectra [see Figures 2.2(b),  2.16, and Supplementary information 
Figure 3.17] within the expected range of 1.8 eV-2.2 eV for room-temperature visible range 
emission in hBN under 532 nm (2.33 eV) excitation15,23,40. We calculated the average PL of 
exposed regions at each level by image analysis on confocal PL scan plots segmented into 
regions of 4 μm2, with 1 patterned site per 1μm2. We also calculated the average SPE yield of 
each segmented 4 μm2 region by counting the number of emitters with 𝑔𝑔2(0) < 0.5 and dividing 
by 4, the number of patterned sites. For more details on the average PL and SPE yield 
calculations, refer to 3.4 Methods.  
 We observed PL emission from successfully patterned sites exposed at all levels with region 
PL and SPE yield dependent on exposure level. Patterned sites that bordered or overlapped with 
an hBN flake edge accounted for ~1% of all patterned sites and exhibited low PL and SPE yield 
regardless of exposure level. We expect that partially patterned circles created structural 
environments inconsistent with the fully contained patterned sites intended. Therefore, we 
considered circle patterns bordering or overlapping with hBN edges to be unsuccessfully 
fabricated and excluded those sites from our analysis. Of the regions that hosted non-zero SPE 
yields, we found that patterned sites exposed at the third level hosted SPEs with a yield of 
0.047 ±  0.016 (mean comparison test to zero: p = 0.0027) and the fourth hosted SPEs with a 
yield of 0.115 ±  0.025 (mean comparison test to zero: p < 0.0001) [Figure 3.5(c)]. The regions 
exposed at the two lowest levels, that did not host SPEs, exhibited brighter regions on confocal 
scan plots, with 𝑃𝑃𝑃𝑃1 = 573 ± 27 cts and 𝑃𝑃𝑃𝑃2 = 1097 ± 56 cts, while the third level exhibited 
65 
  
Figure 3.5. Summary of photophysical characterization by FIB exposure (a) Confocal PL map 
showing all 4 levels of FIB exposure. Dashed line marks the flake edges. Scale bar: 4 𝜇𝜇m. (b) 
Mean values Region PL intensity as a function of FIB exposure time, showing statistically 
significant differences between each level when compared via Tukey-Kramer HSD all-pairs 
comparison test. (c) Mean values SPE Yield (number of SPEs per 4 milled sites) as a function of 
FIB exposure time, showing statistically significant differences between the highest exposure 
times and all others when compared via Tukey-Kramer HSD all-pairs comparison test. The lines 
connecting data points in (c) and (d) are guides to the eye. Error bars in both plots show standard 
error on the mean. 
lower PL on average with 𝑃𝑃𝑃𝑃3 = 186 ± 25 cts (Tukey test: p < 0.0001) and the fourth even 
lower on average with 𝑃𝑃𝑃𝑃4 = 46 ± 3  cts (Tukey test: p < 0.0001) [Figure 3.5(b)]. 
 We conjecture that brighter PL indicates that many optically active vacancy defects remained 
within the crystal at exposed regions, so excited locations within the diffraction limited spot 
exhibited ensemble emission rather than single photon emission. We further expect that higher 
SPE yields were possible when many optically active defects were removed by milling, reducing 
the overall background PL. This conjecture is consistent with the highest SPE yield, lowest 
66 
  
region PL, and largest milling observed at the highest FIB exposure level. However, it is 
inconsistent with the non-zero SPE yield, lower region PL, and lower milling observed at the 
third exposure level. Therefore, the morphological characteristic of deep and directed milling 
cannot be fully responsible for the optical characteristics of FIB exposed hBN crystals. 
When comparing all calculated morphological characteristics to the measured optical 
characteristics [Figure 3.6(a)], we found that the PL intensity at each FIB exposure level was 
strongly correlated only to the maximum swell height at each FIB level, with c(𝑃𝑃𝑃𝑃, ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒) =
 −0.969 (p = 0.0309), while the SPE yield at each FIB level was strongly correlated only to the 
average volume of swell areas at each FIB level with c(𝑆𝑆𝑃𝑃𝐸𝐸, ∆𝑉𝑉+) = 0.9955 (p = 0.0045)  
[Figure 3.6(b)]. Furthermore, the observed increases in material swell height and volume were 
the only morphological features that correlated even weakly (0.05 < p < 0.20) to increased SPE 
Figure 3.6. Correlating material damage to SPE yield and photoluminescence (a) 
Correlation matrix for calculated morphological characteristics and measured optical 
characteristics, with p-value denoting the significance of the correlation strength (b) Plots 
showing the correlation by 1𝜎𝜎 density ellipse between the morphological characteristics and 
optical characteristics. Each data point represents the mean calculated value of the specified 
optical characteristic at a given exposure level plotting at the corresponding mean calculated 
value of the specified morphological characteristic at the same exposure level.  Error bars show 
standard error on the mean for the calculated values at each exposure level along the specified 
axis. 
67 
  
yield and diminished PL simultaneously. We expect that the morphological characteristic of 
swell is indicative of damage levels within the crystal that decrease PL and increase the 
likelihood of hosting SPEs prior to deep and directed milling.  
 To elucidate possible shared causes for both the observed swelling due to FIB damage and 
the observed optical characteristics on the atomic level, we performed SRIM calculations for 
vacancy concentrations. We performed the Monolayer Collision Steps/Surface Sputtering 
calculation, which was considered optimal for the thickness range of flakes studied in this 
work81. From the SRIM calculations, we acquired the vacancy number per ion as a function of 
depth and lateral position 𝑁𝑁𝑣𝑣/𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧)for a 20 keV Ga+ beam entering the material at a single 
point, where SRIM assumes a beam diameter of 0. 
To determine the ion exposure profile 𝐷𝐷𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑦𝑦) for 
the full 500 nm circular pattern, we modeled the 
beam as a 2D gaussian82 and summed over all points 
in the pattern following a circular beam path, with 
the beam pitch set constant as in the patterning 
software for all patterning experiments [Figure 3.7]. 
 We simulated the vacancy concentration 
𝑛𝑛𝑣𝑣(𝑥𝑥, 𝑧𝑧) in exposed regions by convolving 
𝑁𝑁𝑣𝑣/𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧) at a given point with the simulated 
Figure 3.7. Simulated FIB 
exposure profile 𝐷𝐷𝐺𝐺𝑏𝑏+(𝑥𝑥) at that point over a 1 μm exposure dose profile calculated 
using a 2D Gaussian beam 
position range. Lastly, to mimic the material [Equation 3.15] scanned in a 
circular trajectory for a 500 nm circle 
conditions upon which optical data was collected, pattern as set in the FIB patterning 
software with corresponding profile 
we further convolved 𝑛𝑛𝑣𝑣(𝑥𝑥, 𝑧𝑧) with the height data 2showing ions/nm  as a function of 
position 
68 
  
collected by AFM at each point, removing areas of the simulated material at locations where the 
experimental material was milled a depth −∆𝑧𝑧. For more information about simulation 
parameters and how we accounted for milling at each increased FIB exposure level, refer to 3.4 
Methods.  
Based on the results of these simulations, shown in Figure 3.8, we observed that the 𝑛𝑛𝑣𝑣(𝑥𝑥, 𝑧𝑧) 
increased steadily with each increase in FIB exposure level. Furthermore, while the significant 
increase in milling at the fourth exposure level removed a large volume of regions with high 
vacancy concentrations, many regions with high 𝑛𝑛𝑣𝑣 remained due to ∆𝑧𝑧−𝑏𝑏𝑏𝑏𝑒𝑒 less than the 
starting thickness of material ℎ0. This is consistent with our hypothesis that although milling 
away high 𝑛𝑛𝑣𝑣 areas contributed to an overall lower number of optically active defects at the 
Figure 3.8. SRIM vacancy concentration profile 𝑛𝑛𝑣𝑣(𝑥𝑥, 𝑧𝑧) for FIB exposure level 1 (top left), 2 
(top right), 3 (bottom left) and 4 (bottom right) with milled regions removed according to 
representative AFM line cuts 
69 
  
fourth level, the onset of milling distinct circular holes alone cannot fully justify the PL decrease 
and onset of isolated SPEs we observed. Furthermore, since the level of swelling correlated 
strongly to these observed optical characteristics, we expect that there exists some threshold 
vacancy concentration 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣  that leads to the observed increase in swelling and therefore the 
observed changes to optical characteristics. In addition, by combining 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣  with a 
corresponding threshold vacancy concentration 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣  that leads to the observed onset of directed 
and localized milling, we expect that we can determine different regimes of material damage and 
qualitatively describe how they influence the material’s ability to support SPEs at patterned 
locations.  
To find the value for 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣  from our simulated data, we first calculated the onset exposure 
time 𝑡𝑡0− by fitting the normalized remaining thickness 1 − |∆��𝑧𝑧��−�(𝑡𝑡)|⁄ℎ0 to a complementary 
error function with constant offset, finding the intersection 𝑡𝑡0− of the maximum amplitude with 
the negative slope linear regime, and calculated 𝑛𝑛𝑣𝑣(𝑡𝑡0−) from a linear fit of  𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡). This 
process is illustrated in Figure 3.9, which 
shows the normalized remaining thickness fit 
plotted with 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡). We similarly found 
𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣  by fitting �∆�𝑧𝑧��+�(𝑡𝑡) to an error function 
with constant offset, finding the intersection 
𝑡𝑡0+ of the minimum amplitude with the 
positive slope linear regime, and calculated Figure 3.9. Estimating damage thresholds 
from simulated vacancies Average remaining 
𝑛𝑛𝑣𝑣(𝑡𝑡0+) from the linear fit 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡). From thickness calculated from average mill depth at 
each exposure level used to find 𝑡𝑡0−, denoted by 
these calculations, we found 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠~150 nm-3 the vertical line between 4.0 ms and 5.0 ms, and 𝑣𝑣
calculate 𝑛𝑛𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣 , marked by the star at 
and 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠~243 nm-3𝑣𝑣 . Furthermore, we 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡0−). Error bars on thickness points are standard error on the mean. 
70 
  
similarly fit ∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒 to calculate a simulated concentration that corresponded to the maximum 
observed swelling and obtained 𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒.  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 ~186 nm-3.  
Using these threshold values, we simulated the change in thickness profile as a function of 
position Δ𝑧𝑧(𝑥𝑥) for the third and fourth FIB exposure levels by plotting the maximum vacancy 
concentration as a function of position 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑥𝑥) with positive increasing values if the swell 
threshold was crossed and before the maximum swell, positive decreasing values if the swell 
threshold was crossed and after the maximum swell, and with negative values if the mill 
threshold was crossed, as shown in Figure 3.10. For more details about the simulated thickness 
profile calculations and fits, refer to 3.4 Methods. 
We observed that the simulated height profile was highest towards the center of the exposed 
region for the third exposure level, whereas the simulated height profile for the fourth exposure 
level was highest towards the edges, corresponding to the locations where swelling was 
observed, in good qualitative agreement with our AFM height data [Figure 3.10], as well as 
results from Glushkov et. al.73 that show swelling around the edges of circles patterned with 
high-energy Xe+ FIB. Additionally, results from Glushkov et. al. show that emission is localized 
to the edges of these patterned circles with swelling, thereby demonstrating that optically active 
vacancy defects that exhibit photophysical properties of SPEs remain in regions where swelling 
is observed. Because we only observed swelling in exposed regions with 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 < 𝑛𝑛 < 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣 𝑣𝑣 𝑣𝑣 , 
we expect that vacancy concentrations within this range correspond to a distinct type of 
crystallographic defect, causing higher swelling, reduced region PL, and the onset of non-zero 
SPE yield. 
Supported by the morphological characteristics calculated from AFM height data, the optical 
characteristics calculated from confocal PL and antibunching data, and the vacancy 
71 
  
Figure 3.10. Simulated swelling and milling via vacancy concentration thresholds Simulated 
change in height profile ∆𝑧𝑧(𝑥𝑥) calculated from normalizing 𝑛𝑛 (𝑥𝑥) according to 𝑛𝑛𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒 𝑣𝑣 , 
𝑛𝑛𝑏𝑏𝑏𝑏𝑒𝑒. 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 , and 𝑛𝑛𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣 , with positions where these thresholds are met denotated on 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑥𝑥). 
AFM height profile as plotted in Fig. 1(b) for the corresponding exposure level is reproduced 
here to show qualitative agreement with the simulated height profile. 
concentration characteristics calculated from SRIM and ion exposure simulations, we expect that 
at some critical vacancy concentration 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 < 𝑛𝑛 < 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣 𝑣𝑣 𝑣𝑣  present after FIB exposure at the 
third and fourth levels studied in this work, vacancies coalesce into larger volumetric crystal 
defects, or voids, indicated by the swelling level56,73,81, that are not optically active at an 
excitation wavelength of 532 nm [Figure 3.11]. The process of void nucleation and growth in 
irradiated materials is well studied for materials similar to hBN83,84 and is consistent with recent 
investigations into nanopore growth in irradiated monolayer hBN85–87. Furthermore, the atomic 
density 𝑛𝑛⁄𝑉𝑉𝑐𝑐 of mulitlayer hBN with lattice constants 𝑎𝑎 = 2.5 Å and 𝑒𝑒 = 6.6 Å88 with AA′ 
stacking55 is approximately 166 nm−3, which lies within the simulated range for 𝑛𝑛𝑣𝑣 where 
72 
  
Figure 3.11. Illustration of vacancy distribution, void nucleation and growth, and milling, 
demonstrating the observed PL and SPE yield for each FIB exposure level as dependent on the 
distribution of optically active vacancies and non-optically active voids. Gray areas illustrate 
hBN, white circles with glowing edges illustrate vacancy defects that are color centers, and white 
rectangles illustrate voids 
swelling is observed. This is consistent with our expectation that voids form when vacancies 
dominate the underlying crystal structure, which can be described as 𝑛𝑛𝑣𝑣 > 𝑛𝑛⁄𝑉𝑉𝑐𝑐. 
 While local layer mixing56, Ga+ and recoil ion implantation81, and carbon impurity doping 
from the emitter activation process may contribute in small parts to the observed swelling, they 
cannot fully describe the sudden transition in optical characteristics observed between the second 
and third FIB exposure levels without considering a critical 𝑛𝑛𝑣𝑣 that causes void nucleation. 
Because voids act as sinks for vacancies89, the presence of voids in higher quantities starting at 
the third exposure level further supports the observed decrease in region PL and better isolation 
of SPEs at patterned sites despite similar milling to the first and second exposure levels. Smaller, 
optically active vacancy defects would no longer be evenly distributed within the crystal causing 
background fluorescence during emitter excitation and would instead be concentrated at the void 
locations, coalesced as part of the large volumetric defects that would have different optical 
73 
  
characteristics, and would not change the overall vacancy number, remaining consistent with our 
vacancy concentration calculations.  
 Considering the results, calculations, and analysis presented in this work, we propose the 
following model for the fabrication of isolated SPEs via high-energy Ga+ FIB milling of pristine 
hBN crystals:  
 With no FIB exposure, pristine hBN hosts no SPEs. Un-exposed and carbon-annealed hBN 
hosts randomly located SPEs with a density of 0.0125/μm2 (mean comparison test to zero: p = 
0.0416) likely due to low concentrations of native vacancy defects forming complexes with 
carbon impurities41,90. 
 At the lowest FIB exposure level (1.3 ms), exposed areas host many optically active 
vacancy defects [Figure 3.5(a)-(b), Figure 3.8] with visible-range emission under 532 nm 
excitation, and we do not observe isolated SPEs and ZPLs with traditional diffraction limited 
confocal techniques. This is characterized by the high region PL [Figure 3.5(b)] and zero SPE 
yield [Figure 3.5(c)]. Some surface level milling and sputtering occurs and few areas with 
particularly high vacancy concentrations contain voids. This is characterized by randomly 
located pockets of milling and nonzero swelling within the exposed area [Figure 3.2]. 
 At the next FIB exposure level (2.5 ms), exposed areas host more optically active vacancy 
defects [Figure 3.5(a)-(b), Figure 3.8], and we again did not isolate SPEs with traditional 
confocal techniques. This is characterized by the similar inconsistent sputtering of the exposed 
surface compared to the first level (t-test: p = 0.1566) [Figure 3.2, Figure 3.4], the 90% increase 
in region PL (t-test: p < 0.0001) [Figure 3.5(b)], and zero SPE yield [Figure 3.5(c)]. Areas with 
particularly high vacancy concentrations host voids. This is characterized by small areas of 
74 
  
randomly located nonzero swelling similar to the first level (t-test: p = 0.3245) [Figure 3.2, 
Figure 3.4]. 
 At the third FIB exposure level (3.8 ms), exposed areas have even higher concentrations of 
vacancies [Figure 3.8]. Many more regions have reached a critical concentration of vacancies, 
and therefore contain more voids than optically active vacancies or complexes, but voids have 
not coalesced into structural defects that cause significant material breakage and directed milling. 
This is characterized by surface level milling and sputtering similar to the previous level (t-test: p 
= 0.5058) [Figure 3.2, Figure 3.4], the 44% increase in surface swelling height (t-test: p < 
0.0001) [Figure 3.4] localized towards the center of the exposed area [Figure 3.2, Figure 3.10], 
the 83% decrease in region PL (t-test: p < 0.0001) [Figure 3.5(a)-(b)], and the onset of non-zero 
SPE yield (mean comparison test to zero: p = 0.0027) [Figure 3.5(c)]. 
 At the fourth and highest FIB exposure level tested (5.0 ms), exposed areas have the 
highest concentrations of vacancies, but as voids continue to grow, complete breakage and 
therefore milling of the material occurs, leaving behind voids and vacancies along the edges of 
the milled pattern [Figure 3.8, Figure 3.10]. This is characterized by the 62% increase in mill 
depth (t-test: p < 0.0001) [Figure 3.2, Figure 3.4(a)], the 147% increase in exposed edge surface 
area (t-test: p < 0.0001) [see Supplementary information Figure 3.11], the 287% increase in 
milled volume (t-test: p < 0.0001) [Figure 3.2, Figure 3.4(b)], and the similar level of surface 
swelling height (t-test: p = 0.8108) but 23% increase in swelling volume (t-test: p = 0.0444) 
[Figure 3.2, Figure 3.4]. Furthermore, voids act as sinks for many remaining vacancies causing 
the growth of voids and decrease in isolated vacancies in nearby regions, further driven by 
diffusion during the high temperature annealing used to activate emitters. Therefore, the number 
of optically active defects causing background continues to decrease. This is characterized by the 
75 
  
75% decrease in region PL (t-test: p = 0.0174) [Figure 3.5(a)-(b)] and the 147% increase in SPE 
yield (t-test: p = 0.0073) from the third level [Figure 3.5(c)]. 
 The spatial distribution and relative concentrations of vacancies, complexes, and voids could 
be confirmed by high resolution electron microscopy characterization techniques such as 
transmission electron microscopy (TEM) or scanning transmission electron microscopy (STEM) 
performed on cross-sections sliced from hBN crystals deterministically patterned with high-
energy Ga+ FIB. Studies that explore the large range of Ga+ ion energies and fluences with fine 
control could quantitatively measure damage thresholds that correlate to the observed 
morphological characteristics. Furthermore, combining these characterization techniques with 
computational studies such as Molecular Dynamics (MD) and Density Functional Theory could 
elucidate void nucleation and growth dynamics while also considering stable quantum emitter 
defect candidates that include carbon incorporated during the annealing process and allow for the 
creation of SPEs fabricated with this technique. 
3.3. Conclusions 
 In this work, we demonstrated that morphological characteristics due to FIB exposure play an 
important role in the deterministic placement of SPEs by patterning SPEs in hBN with four 
levels of FIB exposure and characterizing surface topography via AFM. From our analysis, we 
found that the evolution of swelling most strongly correlates to the number of SPEs we observed 
and anticorrelates to PL intensity. We expect that swelling indicates volumetric voids and local 
amorphization of hBN layers, with these conclusions supported by SRIM calculations and 
literature that employ MD modeling25,56,89. Therefore, we conjecture that while fluorescent 
vacancy defects are present at all FIB doses, the increased damage indicated by the swelling 
onset at the third exposure level and the milling onset at the fourth actually decreases the total 
76 
  
number of optically active defects at a given milled site, allowing for the localization of single 
emitters.   
 Further experimental and computational studies are necessary to identify the ideal damage 
threshold to localize emitters in hBN without high-cost and high-precision individual defect 
placement techniques. For suspended structures such as photonic waveguides and 
optomechanical transducers, it may be desirable to identify an ion species that can produce a 
higher yield of localized emitters with minimal sputtering, demonstrated at our third FIB level.70–
72,91,92 Another possible route could be to optimize the milled edge profile for SPE creation and 
design photonic devices that rely on the creation of these edges. By identifying the 
morphological characteristics that can be leveraged to engineer hBN SPEs for quantum 
technologies, we have provided insight into the evolution of these defects-by-design that can be 
applied to the expanding field of solid-state single-photon emitter sources. 
3.4. Methods 
Sample Preparation 
 HPHT hBN Exfoliation: We prepared the samples studied in this work on Si/SiO2 substrates 
with 300 nm of thermal wet oxide in order to optically observe hBN thin films. Substrates were 
patterned with a number array via laser writing and reactive ion etching for labeling purposes 
only. After the labeling process, chips were cleaned in acetone sonication bath for 5 min and IPA 
sonication bath for 5 min. Cleaned chips were exposed to 300 mTorr O2 plasma at 300 W RF 
power in a March plasma cleaner for 1 min to create a hydrophilic surface for material transfer 
adhesion. HPHT hBN crystals were prepared for exfoliation transfer on thermal release tape. The 
thermal release tape with hBN crystals was placed on the substrate surface immediately after 
oxygen plasma exposure. The substrate-hBN-tape stack was placed on a hot plate 120°C for 5 
77 
  
minutes to release the tape from the substrate-hBN stack. The sample was cleaned of any tape 
residue in a quartz tube furnace open to air at 500°C for 3 hours with no gas flow. 
 FIB Patterning: FIB exposure was performed in a FEI Helios Dual Beam with Ga+ source. 
Prior to sample loading, the chamber was plasma cleaned with air for 5 min. The FIB was set to 
acceleration voltage 20 kV and beam current 59 pA for all experiments according to main text 
reference 23. The 500 nm diameter circular pattern exposure mask was drawn directly in the Ion 
Beam quadrant of the instrument user interface with the following pattern parameters: Overlap X 
= 50%, Overlap Y = 50%, Pitch X = 12.64 nm, Pitch Y = 12.64 nm, Loop Time = 1.3 ms, Scan 
Type = Circular, Fill Style = Solid, Total Diameter = 25.3 nm. Once the pattern was drawn in the 
quadrant, we created an array of 500 nm circles with 1 μm pitch. To vary the ion fluence, the 
exposure time was set in pattern execution tool to 1.3 ms, 2.5 ms, 3.8 ms, or 5.0 ms depending on 
the experiment, thereby increasing the number of passes where 1.3 ms = 1 pass, 2.5 ms = 2 
passes, 3.8 ms = 3 passes, and 5.0 ms = 4 passes. The order of exposure times varied for each 
hBN flake patterned. 
 SPE Activation: Carbon-rich annealing was performed according to 2.4 Methods. 
 AFM Characterization and Calculations: AFM height data was collected with a Bruker 
Dimension FastScan AFM in Tapping Mode. AFM scans were collected with a 1:1 aspect ratio 
with 256 points per line and 256 lines per scan. Scans was post-processed with NanoScope 
Analysis Flatten and Plane Fit tools. 
 Morphological characteristics were calculated from AFM height data collected during 16 
μm2 AFM scans on each hBN flake at each exposure level, for a total of 16 AFM scan images. 
Scans were segmented into 𝐴𝐴 = 1 μm2𝑅𝑅  area regions, so each region would contain one exposed 
500 nm circle per region. Segmented regions contained 𝑁𝑁𝑒𝑒 = 64 points per line and 𝑁𝑁𝑦𝑦 = 64 
78 
  
lines per scan for a total of 𝑁𝑁 = 4096 height data points [Figure 3.12(a)]. Exposed circles were 
roughly in the center of each region, so data points in the outer perimeter could be taken as 
unexposed [Figure 3.12(b)]. To determine the unexposed zero height 𝑧𝑧0 for calculations in each 
region, height values 𝑧𝑧𝑚𝑚 in the region perimeter were averaged together by 
∑𝑁𝑁0 𝑧𝑧
𝑧𝑧 = 𝑚𝑚 𝑚𝑚0 (3.1) 𝑁𝑁0
where 𝑁𝑁0 is the number of data points in the designated unmilled perimeter. The average 
unmilled roughness of each region was calculated from the same data points relative to the zero 
height: 
∑𝑁𝑁0𝑚𝑚 (𝑧𝑧𝑚𝑚 − 𝑧𝑧0)
2
𝑅𝑅𝑅𝑅𝑆𝑆0 = (3.2) 𝑁𝑁0 − 1
All height data points in the region were shifted by 𝑧𝑧0 [Figure 3.12(c)] and used to calculate the 
morphological characteristics below, where 𝑧𝑧<𝑅𝑅𝑅𝑅𝑅𝑅0 denotes height values less than −𝑅𝑅𝑅𝑅𝑆𝑆0, 
𝑧𝑧>𝑅𝑅𝑅𝑅𝑅𝑅0 denotes height values greater than 𝑅𝑅𝑅𝑅𝑆𝑆0, 𝑁𝑁𝑧𝑧<𝑅𝑅𝑅𝑅𝑅𝑅  is the number of data points with 0
height values less than −𝑅𝑅𝑅𝑅𝑆𝑆0, 𝑁𝑁𝑧𝑧>𝑅𝑅𝑅𝑅𝑅𝑅  is the number of data points with height values greater 0
than 𝑅𝑅𝑅𝑅𝑆𝑆 , ∆𝑥𝑥 = 1 μm0  is the distance between points in a line, and ∆𝑦𝑦 =
1 μm is the distance 
𝑁𝑁𝑥𝑥 𝑁𝑁𝑦𝑦
between lines in a scan. 
Average mill depth: 
∑𝑁𝑁𝑥𝑥,𝑁𝑁𝑦𝑦 𝑧𝑧
∆��𝑧𝑧��� = 𝑚𝑚,𝑗𝑗
<𝑅𝑅𝑅𝑅𝑅𝑅0�𝑥𝑥𝑚𝑚 ,𝑦𝑦𝑗𝑗�
−  (3.3) 𝑁𝑁𝑧𝑧<𝑅𝑅𝑅𝑅𝑅𝑅0
Maximum mill depth: 
∑100 �𝑧𝑧<𝑅𝑅𝑅𝑅𝑅𝑅 �
∆𝑧𝑧 = 0 𝑏𝑏𝑏𝑏𝑒𝑒 (3.4) −𝑏𝑏𝑏𝑏𝑒𝑒 10
79 
  
Figure 3.12. Analyzing AFM height data to calculate morphological characteristics, (a) AFM 
height scan image (size 4 𝜇𝜇m x 4 𝜇𝜇m) of 500 nm circles patterned at the highest FIB dose 
showing holes milled in the hBN surface. White boxes show examples of 1 𝜇𝜇m by 1 𝜇𝜇m 
segmented regions, labeled by column and row: 01 and 32  (b) AFM height scan images of 
segmented regions 01 (left) and 32 (right) with border edges boxed in purple outline used to 
calculate average unmilled height and unmilled roughness for thresholding (c) Comparison of 
original AFM height scan image (processed with Plane Fit and Flatten functions in NanoScope 
Analysis) with the same AFM height scan image after thresholding. Histograms show the 
distribution of height data points for the full 4 𝜇𝜇m x 4𝜇𝜇m scan (purple), the segmented 1 𝜇𝜇m x 
1𝜇𝜇m regions (pink), and edges used in threshold calculation (tan). Solid black line shows average 
surface height and dotted black lines show the threshold for average unexposed roughness. 
Average swell height: 
∑𝑁𝑁𝑥𝑥,𝑁𝑁𝑦𝑦𝑚𝑚,𝑗𝑗 𝑧𝑧>𝑅𝑅𝑅𝑅𝑅𝑅0�𝑥𝑥𝑚𝑚 ,𝑦𝑦 �∆��
𝑗𝑗
𝑧𝑧��+� =  (3.5) 𝑁𝑁𝑧𝑧>𝑅𝑅𝑅𝑅𝑅𝑅0
Maximum swell height: 
∑100 �𝑧𝑧>𝑅𝑅𝑅𝑅𝑅𝑅 �
∆𝑧𝑧 = 0 𝑏𝑏𝑏𝑏𝑒𝑒+𝑏𝑏𝑏𝑏𝑒𝑒 (3.6)
 
10
 
Exposed edge surface area:  
80 
  
𝑁𝑁𝑥𝑥−1,𝑁𝑁𝑦𝑦−1
1 �(𝑎𝑎1 + 𝑏𝑏1 + 𝑒𝑒) ∗ (−𝑎𝑎1 + 𝑏𝑏1 + 𝑒𝑒) ∗ (𝑎𝑎1 − 𝑏𝑏1 + 𝑒𝑒) ∗ (𝑎𝑎1 + 𝑏𝑏1 − 𝑒𝑒)
𝑆𝑆𝐴𝐴 = � � �
4  
𝑚𝑚,𝑗𝑗 +�(𝑎𝑎2 + 𝑏𝑏2 + 𝑒𝑒) ∗ (−𝑎𝑎2 + 𝑏𝑏2 + 𝑒𝑒) ∗ (𝑎𝑎2 − 𝑏𝑏2 + 𝑒𝑒) ∗ (𝑎𝑎2 + 𝑏𝑏2 − 𝑒𝑒)
−𝐴𝐴𝑅𝑅 (3.7)
where 
 𝑎𝑎1 =  �
2
�𝑧𝑧 2𝑚𝑚,𝑗𝑗 − 𝑧𝑧𝑚𝑚+1,𝑗𝑗� + ∆𝑥𝑥  (3.8) 
� 2𝑏𝑏1 = �𝑧𝑧𝑚𝑚+1,𝑗𝑗 − 𝑧𝑧 2𝑚𝑚+1,𝑗𝑗+1� + ∆𝑦𝑦  (3.9) 
𝑒𝑒 = �
2
�𝑧𝑧 2 2𝑚𝑚,𝑗𝑗 − 𝑧𝑧𝑚𝑚+1,𝑗𝑗+1� + ∆𝑥𝑥 + ∆𝑦𝑦  (3.10) 
𝑎𝑎2 =  �
2
�𝑧𝑧𝑚𝑚,𝑗𝑗 − 𝑧𝑧𝑚𝑚,𝑗𝑗+1� + ∆𝑦𝑦2 (3.11) 
2
𝑏𝑏2 = ��𝑧𝑧 2𝑚𝑚,𝑗𝑗+1 − 𝑧𝑧𝑚𝑚+1,𝑗𝑗+1� + ∆𝑥𝑥  (3.12) 
Mill volume: 
𝑁𝑁𝑥𝑥,𝑁𝑁𝑦𝑦
∆𝑉𝑉+ = � 𝑧𝑧<𝑅𝑅𝑅𝑅𝑅𝑅0�𝑥𝑥𝑚𝑚 ,𝑦𝑦𝑗𝑗� ∗ ∆𝑥𝑥 ∗ ∆𝑦𝑦 (3.13) 
𝑚𝑚,𝑗𝑗
Swell volume: 
𝑁𝑁𝑥𝑥,𝑁𝑁𝑦𝑦
∆𝑉𝑉− = � 𝑧𝑧>𝑅𝑅𝑅𝑅𝑅𝑅0�𝑥𝑥𝑚𝑚 ,𝑦𝑦𝑗𝑗� ∗ ∆𝑥𝑥 ∗ ∆𝑦𝑦 (3.14) 
𝑚𝑚,𝑗𝑗
 Optical Characterization: Confocal PL and time correlated photon counting were 
performed according to 2.4 Methods. 
 Statistical Analysis: All statistical tests were performed in JMP analysis software. For 
analysis of number of SPEs (QEs) and sample PL, samples were divided into 4 𝜇𝜇m2 regions, 
which corresponds to 4 FIB patterned locations for patterned samples because we patterned 
exactly one circle per 1 𝜇𝜇m2. To calculate the mean PL map intensity, we converted all confocal 
81 
  
image points from kcps to cts by the APD detection bin size of 0.1 s and plotted in grayscale 
where 0 corresponded to the minimum counts value and 255 corresponded to the maximum 
counts value for the given PL map. We recorded the mean pixel height in ImageJ for each 
region, converted it to the mean PL intensity in counts, and averaged all regions processed with a 
given combination of experimental input factors. To calculate mean SPE yield, we recorded the 
number of emission sites with 𝑔𝑔2(0) <  0.5 in each 4 𝜇𝜇m2 region and averaged all regions 
processed with a given combination of experimental input factors. We performed mean 
hypothesis tests, pairwise t-tests, and all pairs Tukey-Kramer HSD tests with 𝛼𝛼 = 0.05 
significance level to compare SPE yield and region PL intensity across FIB exposure levels. 
Morphological and optical characteristic correlation was performed in JMP for the average value 
of each property at each FIB exposure level with an 𝛼𝛼 = 0.05 significance level. 
 Vacancy Concentration Simulation: Vacancy profiles were simulated using Stopping 
Range of Ions in Matter (SRIM-2013) Monte-Carlo program with 20 keV Ga+ ion source up to 
2000 ions. The target was set as an hBN-SiO2 stack to mimic experimental conditions [Figure 
3.13].  
 hBN target layer thickness was varied for each run within the range from 1 nm to 15 nm in 1 
nm increments for a total of 15 simulation runs. Damage calculation was set to Monolayer 
Collision Steps/Surface Sputtering according to recommendations for targets of low thickness 
found in the SRIM documentation. hBN layer density was set to 2.1 g/cm3. The hBN 
displacement energies were set according to 93.  SiO2 target layer thickness was varied with the 
hBN thickness to maintain a constant 100 nm total target thickness for the hBN-SiO2 stack. We 
obtained a set of single-point vacancy distribution profiles for each of the specified hBN layer 
thicknesses ℎ = {1, 2, … , 15} as a function of position and depth 𝑁𝑁ℎ𝑣𝑣⁄𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧) in units of 
82 
  
Figure 3.13. Stopping Range of Ions in Matter vacancy concentration (a) SRIM-2013 plot of 
depth showing monolayer collisions calculation for 20 keV Ga+ ions in a 15 nm hBN- 85 nm 
SiO2 target layer stack (b) Vacancy concentration profile in hBN target layer as a function of 
cross-section position and depth 
vacancies/ion-Angstrom and were therefore independent of ion dose. SRIM assumes a beam 
diameter of 0 nm. 
 We modeled the exposure profile of the full pattern 𝐷𝐷(𝑥𝑥,𝑦𝑦) in units of ions/nm2 as the 
superposition of the 2D gaussian beam function of all positions (𝑥𝑥𝑚𝑚 ,𝑦𝑦𝑚𝑚): 
𝑃𝑃
𝐼𝐼𝑡𝑡 1
𝐷𝐷𝐺𝐺𝑏𝑏+(𝑥𝑥,𝑦𝑦) = �
𝑒𝑒𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠
2 exp(−[(𝑥𝑥 − 𝑥𝑥𝑚𝑚)
2 + (𝑦𝑦 − 𝑦𝑦 )2 2𝑚𝑚 ⁄2𝜎𝜎 ]) (3.15) +𝑒𝑒 2𝜋𝜋𝜎𝜎
𝑚𝑚=1
where 𝐼𝐼 is the beam current, 𝑡𝑡𝑒𝑒𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 is the dwell time at each point, +𝑒𝑒 is the elementary charge of 
the Ga+ ion, and 𝜎𝜎 is spread of the gaussian beam. We employed the “apparent beam” model to 
define 𝜎𝜎 as the dose distribution with spot size equal to 108 nm according to 82 due to scattering 
effects, and not the ideal beam diameter referred to in the FIB Patterning subsection of Methods. 
To obtain the ion exposure profile of the full circle pattern, we summed the above 2D gaussian 
over the circle pattern of diameter 𝑑𝑑 equal to 500 nm, updating 𝑥𝑥𝑚𝑚 and 𝑦𝑦𝑚𝑚 for each subsequent 
pixel along the circular pattern trajectory with pitch step 𝑠𝑠 as referenced in the Sample 
Preparation subsection above [Figure 3.14]: 
83 
  
Figure 3.14. Simulating FIB patterning beam path, (a) Plot showing the spiral circle trajectory 
of the focused ion beam during patterning (b) Plot showing overlapping beam positions at each 
point along the trajectory for the specified pitch for Equation 3.15. Both plots evaluated 
according to Equations 3.16-3.18 
𝑑𝑑 𝑠𝑠𝜃𝜃
𝑥𝑥𝑚𝑚 = � −
𝑚𝑚� cos 𝜃𝜃𝑚𝑚 (3.16) 2 2𝜋𝜋
𝑑𝑑 𝑠𝑠𝜃𝜃
𝑦𝑦 = � − 𝑚𝑚𝑚𝑚 � sin 𝜃𝜃𝑚𝑚 (3.17) 2 2𝜋𝜋
Where 𝜃𝜃𝑚𝑚 are given by the pitch constraint: 
𝜃𝜃𝑖𝑖 𝑑𝑑 𝑠𝑠𝜃𝜃
𝑠𝑠 = � 𝑒𝑒(𝜃𝜃)𝑑𝑑𝜃𝜃  where 𝑒𝑒(𝜃𝜃) = − (3.18)
2 2𝜋𝜋  𝜃𝜃𝑖𝑖−1
 To model the vacancy concentration profile 𝑛𝑛𝑣𝑣(𝑥𝑥, 𝑧𝑧) in units of vacancies/nm3 for the 500 
nm circle pattern exposed at the lowest FIB level (1 pass), we convolved the single-point 
vacancy profile for 15 nm hBN 𝑁𝑁15𝑣𝑣⁄𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧)with the dose profile at y = 0 𝐷𝐷𝐺𝐺𝑏𝑏+(𝑥𝑥, 0):  
1000
For 1 pass on un-milled surface: 𝑛𝑛𝑣𝑣(𝑥𝑥, 𝑧𝑧) = � 𝑁𝑁15𝑣𝑣⁄𝐺𝐺𝑏𝑏+(𝑥𝑥 + 𝑥𝑥𝑚𝑚 , 𝑧𝑧) ∗  𝐷𝐷𝐺𝐺𝑏𝑏+(𝑥𝑥𝑚𝑚, 0) (3.19) 
𝑒𝑒𝑖𝑖=0
84 
  
Then, hBN layers that were shown to be milled after 1.3 ms exposure during AFM 
characterization were then removed from the vacancy profile at the corresponding lateral 
position according to height versus lateral position data ℎ1(𝑥𝑥): 
(
𝑛𝑛 (𝑥𝑥 , 𝑧𝑧 ) = �𝑛𝑛𝑣𝑣 𝑥𝑥𝑚𝑚, 𝑧𝑧𝑚𝑚
),  ℎ1(𝑥𝑥𝑚𝑚) < 𝑧𝑧𝑚𝑚
𝑣𝑣,1 𝑚𝑚 𝑚𝑚 0,  ℎ1(𝑥𝑥𝑚𝑚) ≥ 𝑧𝑧
(3.20) 
𝑚𝑚
 To model the vacancy concentration profile 𝑛𝑛𝑣𝑣,𝑏𝑏(𝑥𝑥, 𝑧𝑧) for the next FIB exposure levels (m = 
2, 3, and 4 passes), the single-point vacancy profile obtained from SRIM was again convolved 
with the dose profile at y = 0 𝐷𝐷(𝑥𝑥, 0) for 1 pass with two corrections: (1) 𝑁𝑁ℎ𝑣𝑣⁄𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧) used at 
each 𝑥𝑥𝑚𝑚 of the convolution was dependent on the updated hBN thickness after the previous pass 
ℎ𝑏𝑏−1(𝑥𝑥𝑚𝑚). For example, if the AFM height data showed that ℎ𝑏𝑏−1(𝑥𝑥𝑚𝑚) = −1 nm, then the 
vacancy profile calculated for hBN target layer thickness equal to 14 nm 𝑁𝑁14𝑣𝑣⁄𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧) was used, 
if ℎ𝑏𝑏−1(𝑥𝑥𝑚𝑚) = −2 nm, then the vacancy profile obtained for 13 nm hBN 𝑁𝑁13𝑣𝑣⁄𝐺𝐺𝑏𝑏+(𝑥𝑥, 𝑧𝑧) was used, 
etc.: 
500
For any pass > 1: 𝑛𝑛 (𝑥𝑥, 𝑧𝑧) = �𝑁𝑁15−ℎ𝑚𝑚−1(𝑒𝑒𝑖𝑖)𝑣𝑣 𝑣𝑣⁄𝐺𝐺𝑏𝑏+ (𝑥𝑥 + 𝑥𝑥𝑚𝑚 , 𝑧𝑧) ∗ 𝐷𝐷𝐺𝐺𝑏𝑏+(𝑥𝑥𝑚𝑚 , 0) (3.21) 
𝑒𝑒𝑖𝑖=0
 (2) This new thickness dependent vacancy profile was added on to the complete vacancy 
concentration profile for the previous exposure level before removing layers at positions 
determined by the AFM height data. 
(
( 𝑛𝑛 𝑥𝑥 , 𝑧𝑧
) + 𝑛𝑛 (𝑥𝑥 , 𝑧𝑧 ),  ℎ (𝑥𝑥 ) < 𝑧𝑧
𝑛𝑛𝑣𝑣,𝑏𝑏 𝑥𝑥
𝑣𝑣 𝑚𝑚 𝑚𝑚 𝑣𝑣,𝑏𝑏−1 𝑚𝑚 𝑚𝑚 𝑏𝑏 𝑚𝑚 𝑚𝑚
𝑚𝑚, 𝑧𝑧𝑚𝑚) = � (0,  ℎ (𝑥𝑥 ) ≥ 𝑧𝑧 3.22
) 
𝑏𝑏 𝑚𝑚 𝑚𝑚
 Simulated Height Profile: In order to find threshold 𝑛𝑛𝑣𝑣 values for swelling and milling 
regimes, we first found maximum 𝑛𝑛𝑣𝑣 at each FIB exposure time and fit these values to a linear 
function for 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡): 
𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡) = 𝑚𝑚𝑣𝑣𝑡𝑡 + 𝑏𝑏𝑣𝑣 (3.23) 
85 
  
where 𝑚𝑚𝑣𝑣 is the positive slope in vacancies/(nm3 ms) that characterizes the increase in vacancy 
concentration with increased FIB exposure levels and 𝑏𝑏𝑣𝑣 is the y-axis intercept. 
 To determine the vacancy concentration threshold 𝑛𝑛𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣  after which milling of material 
would be expected, we normalized and converted calculated ∆��𝑧𝑧��−� values to normalized remaining 
thickness at each FIB exposure level: 
∆��𝑧𝑧����(�𝑡𝑡��)
ℎ(𝑡𝑡) = 1 − � − � (3.24) 
ℎ0
where ℎ0 is the starting un-milled thickness of the flake from which ∆���𝑧𝑧�−� was measured. We fit 
ℎ(𝑡𝑡) to an error function with constant offset 𝑒𝑒𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠, which denotes the maximum value of ℎ(𝑡𝑡), 
and extracted best fit parameters: 
𝐴𝐴𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠[1 + 𝑒𝑒𝑒𝑒𝑒𝑒(𝛼𝛼ℎ (𝑡𝑡) = 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠
)]
𝑓𝑓𝑚𝑚𝑡𝑡 + 𝑒𝑒2 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠
(3.25) 
where 𝐴𝐴𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 is the negative amplitude from maximum value to minimum value (𝑒𝑒 > 𝐴𝐴 with 𝐴𝐴 <
0 gives a complementary error function as expected for milling) and 𝛼𝛼𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 = (𝑡𝑡 − 𝜇𝜇𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠)⁄𝜎𝜎𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 
where 𝜇𝜇𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 is the center location of 𝐴𝐴𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠/2 and 𝜎𝜎𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 sets the characteristic width. We then 
found a linear fit for the linear regime of the error function: 
ℎ𝑠𝑠𝑚𝑚𝑝𝑝(𝑡𝑡) = 𝑚𝑚−𝑡𝑡 + 𝑏𝑏𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 (3.26) 
where 𝑚𝑚− is the negative slope in nm/ms that characterizes the milling rate after the exposure 
threshold is reached and 𝑏𝑏𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 is the y-axis intercept. To find the threshold exposure level for 
milling 𝑡𝑡0− we set the linear model equal to the maximum value 𝐴𝐴𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠  + 𝑒𝑒𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠: 
𝑡𝑡0− = (𝐴𝐴𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 + 𝑒𝑒𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠 − 𝑏𝑏𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠)⁄𝑚𝑚− (3.27) 
and then obtained: 
𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠𝑣𝑣 = 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡0−) = 𝑚𝑚𝑣𝑣𝑡𝑡0− + 𝑏𝑏𝑣𝑣 (3.28) 
86 
  
 To determine a vacancy concentration threshold 𝑛𝑛𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣  for the onset of swelling, values for 
∆��𝑧𝑧��+�  at each exposure time were fit to an error function with constant offset of the same form as 
ℎ𝑓𝑓𝑚𝑚𝑡𝑡(𝑡𝑡) above, but with 𝐴𝐴𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 > 0, and best fit parameters were extracted. Then, a linear fit was 
found for the now positively sloped linear regime and threshold exposure level for swelling to 
start 𝑡𝑡0+ was found by setting the linear model equal to the minimum value 𝑒𝑒𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠: 
𝑡𝑡0+ = (𝑒𝑒𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 − 𝑏𝑏𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠)⁄𝑚𝑚+ (3.29) 
and obtaining: 
𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 = 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡0+) = 𝑚𝑚𝑣𝑣𝑡𝑡0+ + 𝑏𝑏𝑣𝑣 (3.30) 
 To determine a vacancy concentration 𝑛𝑛𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣  where maximum swelling would be 
observed, after which the swelling goes down until it reaches 0 prior to milling, values for 
∆𝑧𝑧+𝑏𝑏𝑏𝑏𝑒𝑒 at each exposure time were fit to an error function and best fit parameters were 
extracted. Then, a linear fit was found for the positively sloped linear regime. The exposure level 
for maximum swelling 𝑡𝑡+𝑏𝑏𝑏𝑏𝑒𝑒 was found by setting the linear model to 𝐴𝐴𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 + 𝑒𝑒𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠: 
𝑡𝑡+𝑏𝑏𝑏𝑏𝑒𝑒 = (𝐴𝐴𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 + 𝑒𝑒𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 − 𝑏𝑏max𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠)⁄𝑚𝑚+𝑏𝑏𝑏𝑏𝑒𝑒 (3.31) 
and obtaining: 
𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 = 𝑛𝑛𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒(𝑡𝑡+𝑏𝑏𝑏𝑏𝑒𝑒) = 𝑚𝑚𝑣𝑣𝑡𝑡+𝑏𝑏𝑏𝑏𝑒𝑒 + 𝑏𝑏𝑣𝑣 (3.32) 
 Using the calculated values for 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠, 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 𝑣𝑣 , and 𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠𝑣𝑣 , we simulated the expected 
height profile as ±∆𝑧𝑧 from the reference height 𝑧𝑧0 for FIB level 3 and FIB level 4: 
⎧ 0, 𝑛𝑛𝑣𝑣(𝑥𝑥) < 𝑛𝑛
𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠
𝑣𝑣
⎪ 𝑛𝑛𝑣𝑣(𝑥𝑥) − 𝑛𝑛
𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠
𝑣𝑣 , 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 < 𝑛𝑛 (𝑥𝑥) ≤ 𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠
⎪ 𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 𝑣𝑣 𝑣𝑣 𝑣𝑣𝑣𝑣 − 𝑛𝑛𝑣𝑣
∆𝑧𝑧(𝑥𝑥) = 𝑛𝑛𝑣𝑣(𝑥𝑥) − 𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠1 − 𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠
(3.33) 
⎨ 𝑛𝑛 𝑏𝑏𝑏𝑏𝑒𝑒  𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠 − 𝑛𝑛 𝑠𝑠𝑠𝑠𝑏𝑏𝑠𝑠𝑠𝑠
, 𝑛𝑛𝑣𝑣 < 𝑛𝑛𝑣𝑣(𝑥𝑥) < 𝑛𝑛𝑣𝑣
⎪ 𝑣𝑣 𝑣𝑣
⎪ 𝑛𝑛 (𝑥𝑥) − 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠− 𝑣𝑣 𝑣𝑣 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠
⎩ 𝑛𝑛 𝑏𝑏𝑚𝑚𝑠𝑠𝑠𝑠
, 𝑛𝑛𝑣𝑣(𝑥𝑥) ≥ 𝑛𝑛𝑣𝑣
𝑣𝑣 𝑏𝑏𝑏𝑏𝑒𝑒 − 𝑛𝑛𝑣𝑣
87 
  
3.5. Supplementary information 
Pre-exposure hBN AFM 
 After transfer and cleaning of exfoliated HPHT hBN flakes on Si/SiO2 substrates, we 
characterized the pre-exposure thickness ℎ0 via Tapping-Mode AFM and measured 1um linecut 
across the step from substrate to hBN, shown in Figure 3.15. Furthermore, the AFM height 
scans provide a baseline for the qualitative morphological features of the hBN used in this work. 
We observed none of the randomly located milling and swelling observed in Figure 3.2, and can 
therefore attribute all distinctly visible regions of milling and swelling to be FIB induced. 
Figure 3.15. Representative AFM height scan data of pre-exposure hBN flake studied in this 
work, with corresponding linecut and height measurement. 
Morphological characteristics as a function of exposure time 
 Mean calculated values for average swell height, average mill depth, and exposed edge 
surface area at each FIB exposure level of all samples included in this work are shown below in 
Figure 3.16. The plots for maximum swell height, maximum mill depth, swell volume, and mill 
volume are shown in the main text in Figure 3.4. Average and maximum swell height show 
similar behavior, as do average and maximum mill depth, mill volume, and exposed edge surface 
area. 
88 
  
Figure 3.16. Mean calculated values for average heights ±∆𝒛𝒛 and exposed edge surface area 
𝑺𝑺𝑺𝑺 at each FIB exposure level of all samples included in this work. The lines connecting data 
points in are a guide to the eye. Error bars show standard error on the mean. 
ZPL & 𝒈𝒈𝟎𝟎(𝟎𝟎) range 
 We find ZPL in agreement with expected values for hBN quantum emission and observe a 
distribution of 𝑔𝑔2(0). In Figure 3.12, we show the ranges of ZPL energies and 𝑔𝑔2(0) for all 
emission sites included in the data and analysis in the main text, a subset of data from Chapter 
II. 
Figure 3.17. Measured values for ZPL and 𝒈𝒈𝟎𝟎(𝟎𝟎) found from the three-level system fit of all 
emitters optically characterized in this work, showing different values for FIB exposure level 3 and 
4. On the 𝑔𝑔2(0) distribution, purple diamonds show the average 𝑔𝑔2(0) at the specified exposure 
level and error bars show standard error on the mean standard error on the mean.  
 
89 
  
3.6. Bridge 
 In this chapter, I have reported our analysis of focused ion beam induced morphological and 
structural damage to hexagonal boron nitride crystals during defect engineering. These results 
allowed us to peek “under the surface” so to speak and therefore correlate quantifiable damage 
characteristics to color center creation, which provides a route to better control of single photon 
emitter yields and properties such as the emission energy and quantum purity. While FIB 
exposure damage alone is an important piece of the SPE fabrication process, we also reported in 
Chapter II that the interaction of FIB exposure and carbon annealing time influenced optical 
outputs. Therefore, the next chapter reports our extension of the AFM characterization and 
correlation of observed damage and optical characteristics across multiple anneal times in order 
to refine our model for the physical mechanism of defect creation with the combined FIB-carbon 
annealing technique. Additionally, the simulation results of this chapter are applied to developing 
the diffusion and defect formation parameters used for our calculations to visualize the refined 
model in the following chapter. 
 
 
90 
 
 
CHAPTER IV 
 
MODELING EXTRINSIC IMPURITY DIFFUSION FOR COLOR CENTER  
FORMATION 
 
This chapter contains in-preparation co-authored material; it has been reproduced with 
permission from: R. Klaiss, J. Ziegler, D. Miller, K. Zappitelli, K. Watanabe, T. Tanaguchi, S.K. 
Narayanan, P. Dev, B. J. Alemán, “An adjustable focused ion beam – impurity diffusion 
approach to pattern hexagonal boron nitride quantum emitters.” In this work, I contributed to 
designing the research, performing the research, analyzing the data, and writing the paper. 
 
4.1. Introduction 
 Practical devices applied to quantum information technologies require predictable quantum 
emitter (QE) yields and minimal variation in the emitter properties. However, even the QE 
fabrication techniques that produce the highest yields are extremely material dependent, only 
successful in particular regions of the same sample, while others make it challenging to isolate 
single photons (SPEs).  Methods to localize and align hybrid nanostructures to randomly located 
or unstable emitters do succeed at addressing some of these challenges, but these techniques are 
not scalable to technologies that rely on arrays of SPEs or simplified device geometries. 
Therefore, optimizing a pick-and-place defect engineering technique such as the combined 
focused ion beam (FIB) exposure and carbon annealing technique to pattern SPEs in hexagonal 
boron nitride (hBN) detailed in Chapter II would provide a more controlled route to scalable 
emitter integration.  However, the number of experiments required to test the many possible 
9 1 
  
parameter combinations to optimize this technique scales exponentially, making fabrication 
tedious and impractical. 
 While our results engineering carbon-based QEs in hBN provide an avenue to integrate 
emitters in deterministic locations for photonics and quantum sensing applications, our ability to 
fully optimize this technique is limited by the length of time it takes to fabricate many samples to 
test different fabrication parameters. In order to reduce our reliance on performing many 
fabrication experiments, we develop a model for the proposed physical mechanism of our 
fabrication process by analyzing the interaction of FIB exposure and carbon anneal time through 
material damage characterization. Furthermore, we apply the model to predict and constrain 
optimal annealing time ranges for various material damage levels. Additionally, we discuss 
limitations to the model and how different factors and further calculations could mitigate these 
limitations, leading to a more robust model.  
4.2. Results and discussion 
 Results from our work using atomic force microscopy (AFM) to characterize surface 
morphology [Chapter III] reveal that the four FIB exposure levels on a single sample (no 
variation in anneal time) show distinct differences in mill depth and swell/re-deposition height, 
corresponding to differences in SPE yield and region PL51. Therefore, we extended this 
characterization and analysis to samples fabricated with anneal times 1-3 to elucidate how FIB 
induced damage may be influencing the structural environment for SPE formation and optical 
characteristics at different anneal times [Figure 4.1(a)].   
 In general, regions fabricated with the 3.8 ms and 5.0 ms FIB exposure times had higher 
swelling, deeper milling, and dimmer PL (“high damage”) than regions fabricated with the 1.3 
ms and 2.5 ms exposure times (“low damage”) [Supplementary information Tables 4.1-4.2, 
92 
  
Figure 4.1. Extended atomic force microscopy and photoluminescence characterization (a) 
Top: Array of AFM height scans of 500 nm holes patterned at all 4 FIB exposure and 4 carbon 
anneal times; Bottom: Array of 7 μm x 7 μm representative confocal PL scans for each 
combination of FIB exposure and carbon annealing times (b) Top: Corresponding plots of 
average magnitude swell height and average normalized remaining thickness (calculated from 
mill depth) at each FIB exposure time calculated from AFM scans; Bottom: Interaction plot 
showing average PL map intensity for each combination of FIB exposure and carbon annealing 
times. A, B, C, D, E labels correspond to the same FIB and carbon annealing combinations from 
Figure 2.5 for the conditions that produced the highest probability of QEs. Error bars indicate 
standard error on the mean. 
Chapter II Supplementary information Table 2.9 Tukey-Kramer HSD: p < 0.0001] [Figure 
4.1(b)]. Consistent with our previous results for all but the 20 min anneal, PL got brighter first at 
the 2.5 ms exposure time before it got dimmer [Chapter II Supplementary information Table 
93 
  
2.9 Tukey-Kramer HSD: p<0.0001]. In fact, 2.5 ms FIB exposure and 20 min anneal (D) was the 
only combination that produced SPEs in low damage samples. We conjecture that due to the 
longer exposure time compared to 1.3 ms exposure, the 2.5 ms exposed regions actually 
contained more surface damage not quantifiable as swelling or milling and therefore fewer 
atomic defects available near the surface to form color centers, thereby leading to isolated SPEs 
and dimmer PL overall at the shortest anneal time. 
  Carbon annealing is required for color center formation with this technique. The carbon we 
deposit with this process must diffuse into the crystal and annihilate or bind with point vacancies 
created by FIB exposure. FIB exposure also creates high defect density regions that are 
extremely damaged, demonstrating swelling due to amorphous binding, structural voids, and 
breakage. We expect that this damage changes the migration barrier 𝐸𝐸𝑏𝑏 for impurity in the 
material due to screening, trapping in voids, etc., thereby changing the diffusion constant 𝐷𝐷 ∝
𝑒𝑒−𝐸𝐸𝑏𝑏⁄𝑘𝑘𝐵𝐵𝑇𝑇
2
 and ultimately the time it takes for impurities to diffuse through material (𝑡𝑡~ 𝐹𝐹 �𝐷𝐷).  
 Based on our data that regions patterned with longer FIB exposure (thereby sustaining high 
levels of damage) required longer anneal times to yield SPEs, we propose a multi-layer diffusion 
model for the physical mechanism of defect formation where carbon diffuses through two 
distinct morphologies to form color centers: 1) highly damaged due to continuous, directed 
exposure leading to swelling and milling (𝐷𝐷𝐻𝐻𝐻𝐻) and 2) lower damaged due to the gaussian tail of 
the FIB leading to overall lower dose exposure creating isolated vacancies (𝐷𝐷𝐿𝐿𝐻𝐻) [(Figure 
4.2(a)]. Stopping and Range of Ions in Matter (SRIM) calculations from our previous work51 
show clear distinctions between high defect concentration and low defect concentration regions 
for each dose. We denote the dose dependent characteristic lateral damage length as 𝐹𝐹. [Figure 
4.2(b)].  
94 
  
Figure 4.2. Illustrations demonstrating multilayer diffusion model formulation (a) Sketch 
showing the region limits for multilayer diffusion model based on FIB damage, where labels are 
defined in main text Inset: AFM image of FIB patterned holes showing swelling around edge 
region (b) Sketch demonstrating how FIB induced damage influences the number of potential 
color centers (and therefore dimmer or brighter PL) based on the conditions that impurity carbon 
must diffuse through for each FIB exposure level. 
 We model the carbon concentration for a given FIB exposure level (constant 𝐹𝐹) via diffusion 
in composite media in 1D:  
𝜕𝜕𝐶𝐶 2𝑚𝑚 𝜕𝜕 𝐶𝐶= 𝐷𝐷 𝑚𝑚 (4.1) 
𝜕𝜕𝑡𝑡 𝑚𝑚 𝜕𝜕𝑥𝑥2
𝑖𝑖 = 1,2,3 ;  𝐷𝐷1 = 𝐷𝐷3 = 𝐷𝐷𝐻𝐻𝐻𝐻  ;  𝐷𝐷2 = 𝐷𝐷𝐿𝐿𝐻𝐻  
with perfect matching at interfaces: 
𝜕𝜕𝐶𝐶
𝐷𝐷 𝑚𝑚
𝜕𝜕𝐶𝐶
� = 𝐷𝐷 𝑚𝑚+1𝑚𝑚 𝑚𝑚+1 � (4.2)𝜕𝜕𝑥𝑥  𝑠𝑠− 𝜕𝜕𝑥𝑥 𝑠𝑠+
and constant surface concentration: 
𝐶𝐶(0, 𝑡𝑡) = 𝐶𝐶𝑠𝑠 (4.3) 
𝐶𝐶(1, 𝑡𝑡) = 𝐶𝐶𝑠𝑠 (4.4) 
where 𝐶𝐶𝑅𝑅 is the surface concentration of nanocrystalline graphite deposited on the hBN surface 
during the anneal process.  
95 
  
 For a given color center with formation energy ∆𝐸𝐸𝐻𝐻𝑓𝑓 , the probability that impurity carbon of 
concentration 𝐶𝐶(𝑡𝑡) will annihilate with one vacancy out of 𝑁𝑁𝐻𝐻𝑣𝑣  vacancies to yield a single color 
center (i.e. one isolated SPE) as a function of time would follow a Boltzmann distribution: 
𝑁𝑁𝐻𝐻𝑣𝑣 𝐶𝐶(𝑡𝑡) 𝑒𝑒−𝛽𝛽𝑝𝑝∆𝐸𝐸
𝐷𝐷
𝑓𝑓
𝑝𝑝 𝑛𝑛𝑅𝑅𝑅𝑅𝑃𝑃𝐸𝐸(𝑡𝑡) = (𝐻𝐻 4.5) 𝐷𝐷
1 + ∑𝑁𝑁𝑣𝑣 1 �𝑁𝑁𝑣𝑣 𝐶𝐶
𝑝𝑝
� 𝑒𝑒−𝛽𝛽𝑝𝑝∆𝐸𝐸
𝐷𝐷
𝑓𝑓
𝑝𝑝=1 𝑛𝑛! 𝑛𝑛𝑅𝑅
This model for the SPE yield predicts the experimentally observed behavior of longer carbon 
annealing producing higher SPE yields for high FIB exposure, but lower yields for lower 
exposure [Figure 4.3]. For more information about the calculations and system parameters used 
to evaluate the model, refer to 4.4 Methods and 4.5 Supplementary information. 
Figure 4.3. Evaluated model at tested anneal times Plot of results from the evaluated model 
using a parameter set [Supplementary information Table 4.3] that reproduces the behavior we 
see in data for SPE yield at the tested anneal times 
96 
  
 Furthermore, when we extended this multi-layer diffusion model to anneal times not tested, 
we observed that the anneal time when the SPE yield peaks was highly dependent on the damage 
diffusion constant, the damage length, the number of vacancies available, and the defect 
formation energy, but preserves the same qualitative behavior and saturates at the expected 
Poissonian limit [Figure 4.4]. Full results analyzing these single factor effects on the model 
output as well as interaction effects can be found in 4.5 Supplementary information.  
Figure 4.4. Model parameter dependence of SPE supersaturation characteristics Plot of 
results from the evaluated model labeled by variations to the default parameter set, demonstrating 
the dependence of the single-color center creation probability on the high damage diffusion 
constant, damage length, available vacancy number, and defect formation energy 
 Based on this result, we expect that the ideal annealing time to produce SPEs on samples 
exposed at our first two FIB exposure levels lies below 20 min, as those regions quickly reached 
supersaturation of color centers. Ultimately, we conjecture that SPE fabrication would be most 
successful for a given FIB exposure level with known or calculated damage region length and 
97 
  
vacancy concentration at some onset time prior to the supersaturation limit. Moreover, a more 
complete model including formation energies for multiple defect types could provide a value for 
this time parameter, as the defect formation energy strongly influences the anneal time response 
for optimal SPE probability. 
 To constrain the likely defects formed under these conditions, we refer to previous studies 
that have matched observed ZPLs to calculated ZPLs via Density Functional Theory (DFT) and 
molecular dynamics (MD). Mendelson et. al. put forth the VBC -N  defect center as the source of 
observed visible range single photon emission with ZPL in the range of 585 ± 10 nm (2.12 eV)40, 
while Fischer et. al. favor the VNCB defect center with ZPL observed at 2.2 eV37. However, it 
should be pointed out that, to date, the aforementioned defects remain controversial. For 
instance, our DFT calculations show that VBCN- is a spin-1/2 defect, and not a spin-3/2 defect 
as was calculated for a small few-ring flake of hBN40. Supplementary information Figure 
4.8(b) shows defect formation energy for the VBCN defect in different charge states as a 
function of electronic chemical potential. The three most favorable charge states are 0, -1 and -2, 
out of which only first two are spin active, with spins of 1 and ½, respectively and so, this defect 
might not be the one observed in experiment by Mendelson et. al..  Also, according to our DFT 
calculations, the neutral VNCB is not spin active, with the non-magnetic solution being the lower 
energy state, which was missed by earlier works94. Amongst the energetically most favorable 
charge states (+2, +1 and 0) shown in Supplementary information Figure 4.8(c), only +1 
charge state is spin-active and has a spin-1/2, which was also reported previously95. Additionally, 
our high-energy FIB exposure technique likely produces not only single point vacancies and 
divacancies, but also larger clusters, allowing a class of emitters that contain many carbons in a 
98 
  
cluster96. One such defect candidate in the range of our observed emission is C2CN, modeled by 
Li et. al. with ZPL of 2.13 eV97.  
 We expect that likely defect candidates will vary based on the particular combination of FIB 
exposure and carbon annealing time due to differences in vacancy concentration, cluster size, 
and carbon concentration profile among the samples studied in this work. This is supported by 
the observed ZPL differences reported in Chapter II [Figure 2.6(b)]. Additionally, the structural 
environment within the 2.5 ms exposed regions is markedly different from the environment for 
the 3.8 ms and 5.0 ms exposed regions, altering the initial vacancy concentration and cluster size 
conditions, as well as local strain gradients and substrate interactions as reported in Chapter III 
and in Figure 4.1. Furthermore, the SPEs fabricated with the 2.5 ms exposure time were possible 
only when annealed for the shortest time, while those created in the 3.8 ms and 5.0 ms exposed 
material were created with the highest carbon annealing times, directly impacting the total 
carbon dose available. While we cannot conclusively put forth a particular defect complex 
responsible for the emission we observe from SPEs fabricated with this technique, further 
characterization of the photophysics of these emitters such as polarization and ODMR as well as 
DFT and MD calculations for larger carbon-based defects could narrow the range of defect 
candidates possible for different combinations of FIB exposure and carbon doping levels and 
provide more information to elucidate the anneal time dynamics of color center creation. 
 There are further non-trivial limitations to elucidating a quantitative model for defect 
formation in FIB exposed hBN doped via carbon diffusion. The Monte Carlo method used by 
SRIM to calculate ion induced vacancy concentration do not fully account for sputtering and 
self-annealing, whereas the sputtered volume and stability of large structural voids are important 
considerations in our study of optimal SPE creation. Additionally, evaluating the diffusion 
99 
  
problem in greater than one dimension would provide a more precise concentration model. 
Furthermore, the diffusion constants of carbon in hBN could be even more strongly dependent 
on the surrounding damage environment and vary not just between the high damage and low 
damage regimes but may be a continuous function of position 𝐷𝐷(𝑥𝑥) that is unique for each 
exposure level. However, comparing our results to this qualitative model allows us to constrain 
different testable regimes of FIB exposure and carbon annealing time to optimize and expand 
this fabrication technique. Allowing for more refined control over the design of devices and 
selection of fabrication parameters to align with requirements for different hBN QE reliant 
applications. 
4.3. Conclusions 
 In this work, we used AFM height data and confocal PL images of hBN regions subjected to 
18 different combinations of FIB exposure dose and carbon diffusion time to elucidate that 
carbon impurities likely diffuse through layers of highly damaged material indicated by swelling 
in order to activate color centers. Additionally, by simulating our technique using SRIM and a 
multilayer-diffusion model, we determined that in general, less damaged regions required shorter 
diffusion times with more precision on the order of minutes, while more damaged regions 
required longer diffusion times with less precision on the order of hours, and demonstrated a 
method to constrain the optimal range of diffusion time for a given amount of FIB induced 
damage to produce SPEs. These results, combined with DFT calculations and potential future 
work to characterize polarization and spin characteristics of these patterned emitters could allow 
for both higher resolution spatial control and predictable chemical structure control over defect 
creation, thereby making hBN SPE integration more broadly applicable and scalable for a range 
of quantum information technologies.  
100 
  
4.4. Methods 
 AFM Characterization: AFM height data was collected and analyzed according to 3.4 
Methods. 
 Multi-layer Diffusion Defect Formation Model: We employed the semi-analytical solution 
to the multilayer diffusion problem based on the Laplace transform eigenfunction expansion and 
corresponding MATLAB code published by Carr and Turner to solve for 𝐶𝐶𝑚𝑚(𝑥𝑥, 𝑡𝑡) at each anneal 
time98. Parameter values for diffusivities and damage lengths were determined by order of 
magnitude estimations informed by literature, previous work, and data, and can be found in the 
Supplementary information. 
 We determined an estimate for 𝑁𝑁𝑣𝑣[𝐻𝐻](𝑥𝑥) at each tested exposure dose by convolving SRIM 
calculations with AFM linecuts [refer to morphology paper] and assuming that of all of the 
calculated point vacancies and vacancy complexes, only a small percentage (~1%) would 1) have 
a formation energy low enough to form a defect complex with carbon and 2) when the complex 
is formed, would be optically active under 532 nm excitation. 𝑝𝑝𝑅𝑅𝑃𝑃𝐸𝐸(𝑡𝑡) to produce Figures 4.3 
and 4.4 were calculated in MATLAB using the evaluated solutions for 𝐶𝐶𝑚𝑚(𝑥𝑥, 𝑡𝑡).  
 Ab-initio Calculations: Spin-polarized DFT calculations were performed within the 
generalized gradient approximation99 of Perdew-Burke-Ernzerhof100 as implemented in the 
Quantum Espresso package.101 A 6 × 6 × 1 (72-atoms) supercell of hBN was used in the 
calculations. The Brillouin zone was sampled using a Γ-centered, 4 × 4 × 1 k-grid, which was 
created within the Monkhorst-Pack scheme.102 A 20 Å vacuum layer was added in the direction 
normal to the hBN surface, ensuring minimal interactions between periodic images in the z-
direction. The formation energies of the defects were calculated according to the following 
expression103: 
101 
  
 
∆𝐸𝐸𝐻𝐻𝑏𝑏𝑓𝑓𝑏𝑏𝑐𝑐𝑡𝑡;𝑞𝑞𝑓𝑓𝑒𝑒𝑝𝑝𝑏𝑏 =  𝐸𝐸𝑇𝑇𝑒𝑒𝑡𝑡𝑏𝑏𝑠𝑠[𝑑𝑑𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑡𝑡, 𝑞𝑞] − 𝐸𝐸𝑇𝑇𝑒𝑒𝑡𝑡𝑏𝑏𝑠𝑠[𝑝𝑝𝑒𝑒𝑖𝑖𝑠𝑠𝑡𝑡𝑖𝑖𝑛𝑛𝑒𝑒] −�𝑛𝑛𝑚𝑚𝜇𝜇𝑚𝑚 + 𝑞𝑞(𝐸𝐸𝑉𝑉𝑉𝑉𝑅𝑅[𝑑𝑑𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑡𝑡] + 𝐸𝐸𝐹𝐹) (4.6) 
Here, q is the charge of the defect,  𝑛𝑛𝑚𝑚  (𝜇𝜇𝑚𝑚) represent the number (chemical potential) of the ith 
species of atoms that are removed and/or added, 𝐸𝐸𝑇𝑇𝑒𝑒𝑡𝑡𝑏𝑏𝑠𝑠 represents total energy of the 
defective/pristine supercell, 𝐸𝐸𝐹𝐹 is the Fermi level (electronic chemical potential) relative to the 
valence band maximum (VBM). The 𝐸𝐸𝑉𝑉𝑉𝑉𝑅𝑅 of the defective supercell is calculated using: 
𝐸𝐸𝑉𝑉𝑉𝑉𝑅𝑅[𝑑𝑑𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑡𝑡] = 𝐸𝐸𝑉𝑉𝑉𝑉𝑅𝑅[𝑝𝑝𝑒𝑒𝑖𝑖𝑠𝑠𝑡𝑡𝑖𝑖𝑛𝑛𝑒𝑒] + 𝑉𝑉𝑏𝑏𝑣𝑣[𝑑𝑑𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑒𝑡𝑡] − 𝑉𝑉𝑏𝑏𝑣𝑣[𝑝𝑝𝑒𝑒𝑖𝑖𝑠𝑠𝑡𝑡𝑖𝑖𝑛𝑛𝑒𝑒] (4.7) 
where, the difference in the average potentials (𝑉𝑉𝑏𝑏𝑣𝑣)  is used to line up the band structures of the 
perfect and the defective supercells. Also, for the pristine supercell, one finds the VBM as 
follows:  𝐸𝐸𝑉𝑉𝑉𝑉𝑅𝑅[𝑝𝑝𝑒𝑒𝑖𝑖𝑠𝑠𝑡𝑡𝑖𝑖𝑛𝑛𝑒𝑒] = 𝐸𝐸𝑇𝑇𝑒𝑒𝑡𝑡𝑏𝑏𝑠𝑠[𝑝𝑝𝑒𝑒𝑖𝑖𝑠𝑠𝑡𝑡𝑖𝑖𝑛𝑛𝑒𝑒, 𝑞𝑞 = 0] − 𝐸𝐸𝑇𝑇𝑒𝑒𝑡𝑡𝑏𝑏𝑠𝑠[𝑝𝑝𝑒𝑒𝑖𝑖𝑠𝑠𝑡𝑡𝑖𝑖𝑛𝑛𝑒𝑒, 𝑞𝑞 = +1]. The growth 
conditions are accounted for in the choice of chemical potentials  (𝜇𝜇𝑚𝑚), which represent energies 
of reservoirs with which atoms are being exchanged, along with the requirement that under 
equilibrium conditions: 𝐸𝐸(ℎ𝐵𝐵𝑁𝑁) =   𝜇𝜇𝑉𝑉 +  𝜇𝜇𝑁𝑁. Here, 𝐸𝐸(ℎ𝐵𝐵𝑁𝑁) is total energy of pristine hBN 
unit cell. For example, in the boron rich conditions, the chemical potential for boron, 𝜇𝜇𝑉𝑉, is 
obtained using the total energy of 𝛼𝛼-boron unit cell:  𝜇𝜇𝑉𝑉 = 𝐸𝐸(𝛼𝛼𝐵𝐵)/12, while  𝜇𝜇𝑁𝑁 = 𝐸𝐸(ℎ𝐵𝐵𝑁𝑁) −
 𝜇𝜇𝑉𝑉. In nitrogen-rich condition, N2 is used as a reference to obtain the chemical potential for 
nitrogen,  𝜇𝜇𝑁𝑁, giving  𝜇𝜇𝑁𝑁 = 𝐸𝐸(𝑁𝑁2)/2 and  𝜇𝜇𝑉𝑉 = 𝐸𝐸(ℎ𝐵𝐵𝑁𝑁) −  𝜇𝜇𝑁𝑁. Also, we used graphene as a 
reference to obtain the chemical potential for carbon,  𝜇𝜇𝐶𝐶,  in these calculations. We neglected 
charged defect corrections because they are not expected to affect the qualitative behavior 
comparing reference structure defect formation energies nor spins calculated for this work. 
Therefore, our calculations provide an upper bound on the defect formation energies. 
 
102 
  
4.5. Supplementary information 
AFM Characterization 
 We report the results of Student’s t and Tukey-Kramer HSD Mean Comparison Tests for 
average mill depth and swell height at each exposure time averaged over all anneal times in 
Tables 4.1 and 4.2 respectively, as referenced in 4.2 Results and discussion. 
 
Table 4.1. Average Mill Depth by FIB Exposure time– Mean Comparison Tests where 
averages contain values for samples carbon annealed for 20 min, 1.25 hr, 3 hr, and 5 hr. 𝛼𝛼 = 0.05 
significance level, unequal variance 
Exposure 1 – Exposure 2 Difference Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Difference Each Pair) HSD All Pairs) 
2.5 ms – 5.0 ms 5.22 0.16 < 0.0001 < 0.0001 
1.25 ms – 5.0 ms 4.93 0.16 < 0.0001 < 0.0001 
3.75 ms – 5.0 ms 3.69 0.16 < 0.0001 < 0.0001 
2.5 ms – 3.75 ms 1.52 0.16 < 0.0001 < 0.0001 
1.25 ms – 3.75 ms 1.23 0.16 < 0.0001 < 0.0001 
2.5 ms – 1.25 ms 0.29 0.16 0.0718 0.2729 
 
Table 4.2. Average Swell Height by FIB Exposure time– Mean Comparison Tests where 
averages contain values for samples carbon annealed for 20 min, 1.25 hr, 3 hr, and 5 hr. 𝛼𝛼 = 0.05 
significance level 
Exposure 1 – Exposure 2 Difference Standard Error p-value (Student’s t p-value (Tukey-Kramer 
Difference Each Pair) HSD All Pairs) 
3.75 ms – 2.5 ms 0.570 0.064 < 0.0001 < 0.0001 
3.75 ms – 1.25 ms 0.518 0.064 < 0.0001 < 0.0001 
5.0 ms – 2.5 ms 0.498 0.064 < 0.0001 < 0.0001 
5.0 ms – 1.25 ms 0.446 0.064 < 0.0001 < 0.0001 
3.75 ms – 5.0 ms 0.072 0.064 0.2650 0.6801 
1.25 ms – 2.5 ms 0.052 0.064 0.4163 0.8482 
  
 Data for all calculations from AFM height data collected for this chapter is shown in Figure 
4.5. 
Multi-layer diffusion model 
 When determining the diffusion time for carbon impurities during the CVD/anneal process, 
we do not account for the mass transport time of reactants to the hBN surface due to the high 
concentration of ethylene flow and low-pressure conditions of the CVD process causing growth 
103 
  
Figure 4.5. Full AFM characterization at 4 different anneal times Plots showing average 
values and standard error on mean for all morphological characteristics outlined in Chapter III 
for all FIB exposure-carbon annealing combination (12 hour anneal samples not characterized) 
to occur on shorter time scales than the diffusion. We choose to model in 1D because the in-
plane (between layers) diffusivity of impurities in hBN is orders of magnitude higher than the 
out-of-plane (through layers) diffusivity. Additionally, we do not account for the diffusion of 
vacancy defects because at the process temperature, mono- and di-vacancies quickly reach 
equilibrium concentrations. 
 We determined four different values of l corresponding to expected damage levels for FIB 
exposure based on AFM data, swelling and milling calculations, and SRIM studies performed for 
our previous work. We set the ratio of diffusivities in the high damage and low damage regions 
based on literature values for intercalated carbon in hBN and order of magnitude estimations. All 
model parameters used in the MATLAB calculations are shown in Table 4.3. Furthermore, the 
104 
  
evaluated carbon concentrations for each of the simulated FIB levels intended to replicate our 
experimental conditions are displayed in Figure 4.6. 
Table 4.3. Model Parameters Used for Figure 4.3 
Parameter Variable Value 
Normalized Surface Concentration 𝐶𝐶𝑅𝑅 1 
High-Damage Diffusion Constant (𝜇𝜇𝑚𝑚2/𝑚𝑚𝑖𝑖𝑛𝑛) 𝐷𝐷𝐻𝐻𝐻𝐻 9×10-6 
Low-Damage Diffusion Constant (𝜇𝜇𝑚𝑚2/𝑚𝑚𝑖𝑖𝑛𝑛) 𝐷𝐷𝐿𝐿𝐻𝐻 102 
Lateral Damage Length (nm) {𝐹𝐹1, 𝐹𝐹2, 𝐹𝐹3 𝐹𝐹4} {5, 72, 187, 240} 
Unbound Site Concentration  𝑛𝑛𝑅𝑅 1 
Number of Sites to form defects {𝑁𝑁𝐻𝐻𝑣𝑣 1,𝑁𝑁
𝐻𝐻 𝐻𝐻 𝐻𝐻
𝑣𝑣 2,𝑁𝑁𝑣𝑣 3,𝑁𝑁𝑣𝑣 4} {5, 15, 3, 2} 
Defect Formation Energy (eV) ∆𝐸𝐸𝐻𝐻𝑓𝑓  -0.8 
 
Figure 4.6. Multi-layer diffusion carbon concentration profiles evaluated using the semi-
analytical solution multilayer diffusion problem in MATLAB for parameters [Table 4.3] used to 
calculate the probability of 1 color center at a pattern [Figure 4.3] 
 
105 
  
Color Center Yield Calculation 
 Calculating Boltzmann distribution:  The diagram shown in Figure 4.7 demonstrates the 
process for calculating the single color center probability as a function of time via the Boltzmann 
distribution, where 𝑁𝑁𝐶𝐶 is the number of carbon impurity atoms, 𝑁𝑁𝑅𝑅 is the number of free sites for 
intercalated carbon with energy 𝐸𝐸𝐶𝐶𝑖𝑖, and 𝑁𝑁
𝐻𝐻
𝑣𝑣  is the number of sites where carbon can form a 
stable, optically active defect with energy 𝐸𝐸𝐻𝐻𝑓𝑓 . 
Figure 4.7. Illustration of Boltzmann distribution model for defect formation demonstrating 
the energy and multiplicity calculations for the no color center (n = 0) state, the single color center 
(n = 1) state, and the many color center (n = n) state 
For state with no carbons in stable optically active defect, the energy is 𝑁𝑁𝐶𝐶𝐸𝐸𝐶𝐶𝑖𝑖and the multiplicity 
of that state is 
𝑁𝑁 𝑁𝑁 ! 𝑁𝑁 𝑁𝑁𝐶𝐶
� 𝑅𝑅� = 𝑅𝑅 ≈ 𝑅𝑅 (4.8) 
𝑁𝑁𝐶𝐶 𝑁𝑁𝐶𝐶! (𝑁𝑁𝑅𝑅 − 𝑁𝑁𝐶𝐶)! 𝑁𝑁𝐶𝐶!
106 
  
For state with 1 carbon in stable optically active defect, the energy is (𝑁𝑁𝐶𝐶 − 1)𝐸𝐸 𝐻𝐻𝐶𝐶𝑖𝑖 + 𝐸𝐸𝑓𝑓  and the 
multiplicity of that state is 
𝑁𝑁𝐻𝐻𝑣𝑣 𝑁𝑁𝑅𝑅 𝑁𝑁𝐻𝐻 𝐻𝐻 𝑁𝑁𝐶𝐶−1� � � � = 𝑣𝑣
! 𝑁𝑁𝑅𝑅! 𝑁𝑁𝑣𝑣 𝑁𝑁𝑅𝑅 ( ) 
1 𝑁𝑁𝐶𝐶 − 1 1! (𝑁𝑁𝐻𝐻𝑣𝑣 − 1) (
≈ 4.9
!𝑁𝑁𝐶𝐶! 𝑁𝑁𝑅𝑅 − 𝑁𝑁𝐶𝐶 + 1)! (𝑁𝑁𝐶𝐶 − 1)!
Generalizing to 𝑛𝑛 carbons in stable optically active defect, the energy is (𝑁𝑁𝐶𝐶 − 𝑛𝑛)𝐸𝐸𝐶𝐶𝑖𝑖 + 𝑛𝑛𝐸𝐸
𝐻𝐻
𝑓𝑓  and 
the multiplicity of that state is 
𝑁𝑁𝐻𝐻 𝑁𝑁 𝑁𝑁𝐻𝐻! 𝑁𝑁 ! 𝑁𝑁𝐻𝐻𝑝𝑝𝑁𝑁 𝑁𝑁𝐶𝐶−𝑝𝑝
� 𝑣𝑣 � � 𝑅𝑅 � = 𝑣𝑣 𝑅𝑅 𝑣𝑣 𝑅𝑅𝐻𝐻 ≈ (4.10) 𝑛𝑛 𝑁𝑁𝐶𝐶 − 1 𝑛𝑛! (𝑁𝑁𝑣𝑣 − 𝑛𝑛)!𝑁𝑁𝐶𝐶! (𝑁𝑁𝑅𝑅 − 𝑁𝑁𝐶𝐶 + 𝑛𝑛)! 𝑛𝑛! (𝑁𝑁𝐶𝐶 − 𝑛𝑛)!
The probability of forming exactly 1 defect in the region with 𝑁𝑁𝐶𝐶 carbons and 𝑁𝑁𝑣𝑣[𝐻𝐻] available 
defect sites follows the Boltzmann distribution, generally: 
∑𝑎𝑎𝐹𝐹𝐹𝐹 1 𝑠𝑠𝑡𝑡𝑎𝑎𝑡𝑡𝑒𝑒𝑠𝑠
𝑝𝑝𝑅𝑅𝑃𝑃𝐸𝐸 = (∑ 4.11
) 
𝑎𝑎𝐹𝐹𝐹𝐹 𝑝𝑝𝑝𝑝𝑠𝑠𝑠𝑠𝑖𝑖𝑏𝑏𝐹𝐹𝑒𝑒 𝑠𝑠𝑡𝑡𝑎𝑎𝑡𝑡𝑒𝑒𝑠𝑠
And substituting the state energies and multiplicities from Equations 4.8-4.10: 
𝑁𝑁𝐻𝐻𝑁𝑁 𝑁𝑁𝐶𝐶−1𝑣𝑣 𝑅𝑅 −𝛽𝛽((𝑁𝑁 −1)𝐸𝐸 +𝐸𝐸𝐷𝐷𝐶𝐶 𝐶𝐶 𝑓𝑓 )
( ) 𝑖𝑖
𝑝𝑝 = 𝑁𝑁𝐶𝐶 − 1 !
𝑒𝑒
𝑅𝑅𝑃𝑃𝐸𝐸 (4.12) 
𝑁𝑁𝐷𝐷∑ 𝑣𝑣 𝑁𝑁
𝐻𝐻𝑝𝑝
𝑣𝑣 𝑁𝑁 𝑁𝑁𝐶𝐶−𝑝𝑝𝑅𝑅 𝑒𝑒−𝛽𝛽((𝑁𝑁𝐶𝐶−𝑝𝑝)𝐸𝐸𝐶𝐶𝑖𝑖+𝑝𝑝𝐸𝐸
𝐷𝐷
𝑓𝑓 )
𝑝𝑝=0 𝑛𝑛! (𝑁𝑁𝐶𝐶 − 𝑛𝑛)!
Which we then evaluate: 
𝑁𝑁 𝑁𝑁𝐻𝐻𝑁𝑁 𝑁𝑁𝐶𝐶−1𝐶𝐶 𝑣𝑣 𝑅𝑅 −𝛽𝛽𝑁𝑁 𝐸𝐸 𝛽𝛽𝐸𝐸 −𝛽𝛽𝐸𝐸𝐷𝐷𝑒𝑒 𝐶𝐶 𝐶𝐶𝑁𝑁 ( ) 𝑖𝑖𝑒𝑒
𝐶𝐶𝑖𝑖𝑒𝑒 𝑓𝑓
                                                       = 𝐶𝐶 𝑁𝑁𝐶𝐶 − 1 !
𝑁𝑁 𝑁𝑁𝐶𝐶 𝑝𝑝 𝑁𝑁 −𝑝𝑝
 
𝑅𝑅 𝑒𝑒−𝛽𝛽𝑁𝑁
𝐷𝐷 𝐻𝐻 𝐶𝐶
𝐶𝐶𝐸𝐸𝐶𝐶 + ∑𝑁𝑁𝑖𝑖 𝑣𝑣 𝑁𝑁𝐶𝐶! 𝑁𝑁𝑣𝑣 𝑁𝑁𝑅𝑅 𝑒𝑒−𝛽𝛽𝑁𝑁
𝐷𝐷
𝐶𝐶𝐸𝐸𝐶𝐶𝑖𝑖𝑒𝑒𝛽𝛽𝑝𝑝𝐸𝐸𝐶𝐶 −𝛽𝛽𝑝𝑝𝐸𝐸𝑓𝑓𝑁𝑁𝐶𝐶! 𝑝𝑝=1𝑁𝑁
𝑖𝑖
𝐶𝐶! 𝑛𝑛! (𝑁𝑁𝐶𝐶 − 𝑛𝑛)!
𝑒𝑒
                                               
𝑁𝑁 𝑁𝑁𝐶𝐶 𝐻𝐻𝑅𝑅 𝑒𝑒−𝛽𝛽𝑁𝑁𝐶𝐶𝐸𝐸𝐶𝐶 𝑁𝑁𝑣𝑣 𝑁𝑁𝑖𝑖 𝐶𝐶 𝑒𝑒𝛽𝛽𝐸𝐸
𝐷𝐷
𝐶𝐶
𝑁𝑁 ! 𝑁𝑁 𝑖𝑖𝑒𝑒
−𝛽𝛽𝐸𝐸𝑓𝑓
= 𝐶𝐶 𝑅𝑅
𝑁𝑁 𝑁𝑁
 
𝐶𝐶 𝑁𝑁𝐶𝐶
𝑅𝑅 𝐷𝐷
𝐻𝐻𝑝𝑝 𝐷𝐷
𝑁𝑁 ! 𝑒𝑒
−𝛽𝛽𝑁𝑁𝐶𝐶𝐸𝐸𝐶𝐶 + 𝑁𝑁𝑅𝑅 𝑒𝑒−𝛽𝛽𝑁𝑁𝐶𝐶𝐸𝐸𝐶𝐶 ∑𝑁𝑁𝑣𝑣 𝑁𝑁𝐶𝐶! 𝑁𝑁𝑖𝑖 𝑖𝑖 𝑣𝑣 𝑒𝑒𝛽𝛽𝑝𝑝𝐸𝐸𝐶𝐶 −𝛽𝛽𝑝𝑝𝐸𝐸𝑓𝑓
𝐶𝐶 𝑁𝑁𝐶𝐶! 𝑝𝑝=1 𝑛𝑛! (
𝑖𝑖
𝑁𝑁𝐶𝐶 − 𝑛𝑛)! 𝑁𝑁 𝑝𝑝
𝑒𝑒
𝑅𝑅
𝑛𝑛
Simplifying 𝑁𝑁𝐶𝐶! 𝑁𝑁𝐶𝐶( ) ≈  and setting the defect formation energy ∆𝐸𝐸
𝐻𝐻 𝐻𝐻
𝑝𝑝! 𝑁𝑁 −𝑝𝑝 ! 𝑝𝑝! 𝑓𝑓
= 𝐸𝐸𝑓𝑓 − 𝐸𝐸𝐶𝐶
𝐶𝐶 𝑖𝑖
< 1: 
107 
  
𝑁𝑁𝐻𝐻𝑣𝑣 𝐶𝐶 𝑒𝑒−𝛽𝛽𝑝𝑝∆𝐸𝐸
𝐷𝐷
𝑓𝑓
𝑝𝑝 = 𝑛𝑛𝑅𝑅𝑅𝑅𝑃𝑃𝐸𝐸 ( )𝐷𝐷 𝐻𝐻 𝑝𝑝 4.13  𝐷𝐷
1 + ∑𝑁𝑁𝑣𝑣 1 𝑁𝑁𝑣𝑣 𝐶𝐶 −𝛽𝛽𝑝𝑝∆𝐸𝐸𝑓𝑓𝑝𝑝=1 𝑛𝑛! � 𝑛𝑛 � 𝑒𝑒𝑅𝑅
where the carbon concentration 𝐶𝐶 = 𝑁𝑁𝐶𝐶 and the free site concentration 𝑛𝑛𝑠𝑠 =
𝑁𝑁𝑅𝑅 . 
𝑉𝑉 𝑉𝑉
The same procedure was followed to find the probability of any number of defects, generally: 
∑𝑎𝑎𝐹𝐹𝐹𝐹 𝑛𝑛 > 0 𝑠𝑠𝑡𝑡𝑎𝑎𝑡𝑡𝑒𝑒𝑠𝑠
𝑝𝑝𝐻𝐻 = (∑ 4.14
) 
𝑎𝑎𝐹𝐹𝐹𝐹 𝑝𝑝𝑝𝑝𝑠𝑠𝑠𝑠𝑖𝑖𝑏𝑏𝐹𝐹𝑒𝑒 𝑠𝑠𝑡𝑡𝑎𝑎𝑡𝑡𝑒𝑒𝑠𝑠
Substituting again from Equations 4.8-4.10: 
                                               
𝑁𝑁 𝑁𝑁𝐶𝐶𝑅𝑅 𝑒𝑒−𝛽𝛽𝑁𝑁𝐶𝐶𝐸𝐸
𝐷𝐷
𝐶𝐶 ∑𝑁𝑁𝑣𝑣 𝑁𝑁 ! 𝑁𝑁
𝐻𝐻𝑝𝑝
𝑖𝑖 𝐶𝐶 𝑣𝑣 𝛽𝛽𝑝𝑝𝐸𝐸𝐶𝐶 −𝛽𝛽𝑝𝑝𝐸𝐸
𝐷𝐷
𝑓𝑓
𝑁𝑁𝐶𝐶! 𝑝𝑝=1 𝑛𝑛! (𝑁𝑁𝐶𝐶 − 𝑛𝑛)! 𝑁𝑁 𝑝𝑝
𝑒𝑒 𝑖𝑖𝑒𝑒
= 𝑅𝑅𝑁𝑁  𝑁𝑁 𝐶𝐶 −𝛽𝛽𝑁𝑁 𝐸𝐸 𝑁𝑁 𝑁𝑁𝐶𝐶 𝑝𝑝𝑅𝑅 𝑒𝑒 𝐶𝐶 𝐶𝐶𝑖𝑖 + 𝑅𝑅 𝑒𝑒−𝛽𝛽𝑁𝑁𝐶𝐶𝐸𝐸
𝐷𝐷
𝐶𝐶 𝑁𝑁 𝑁𝑁 ! 𝑁𝑁𝐻𝐻𝑖𝑖 ∑ 𝑣𝑣 𝐶𝐶 𝑣𝑣 𝛽𝛽𝑝𝑝𝐸𝐸𝐶𝐶 −𝛽𝛽𝑝𝑝𝐸𝐸
𝐷𝐷
𝑓𝑓
𝑁𝑁𝐶𝐶! 𝑁𝑁𝐶𝐶! 𝑝𝑝=1
𝑖𝑖
𝑛𝑛! (𝑁𝑁 𝑝𝑝 𝑒𝑒 𝑒𝑒𝐶𝐶 − 𝑛𝑛)! 𝑁𝑁𝑅𝑅
And simplifying to obtain 
𝑁𝑁𝐷𝐷 𝐻𝐻
𝑝𝑝
∑ 𝑣𝑣 1 𝑁𝑁𝑣𝑣 𝐶𝐶 −𝛽𝛽𝑝𝑝∆𝐸𝐸
𝐷𝐷
𝑓𝑓
𝑝𝑝=1 � � 𝑒𝑒
𝑝𝑝 = 𝑛𝑛! 𝑛𝑛𝑅𝑅𝐻𝐻 (4.15) 𝐷𝐷
1 + ∑𝑁𝑁𝑣𝑣 1 �𝑁𝑁
𝐻𝐻 𝑝𝑝
𝑣𝑣 𝐶𝐶 −𝛽𝛽𝑝𝑝∆𝐸𝐸𝐷𝐷𝑓𝑓
𝑝𝑝=1 𝑛𝑛! 𝑛𝑛 � 𝑒𝑒𝑅𝑅
Model Parameters Full Model Effects Test 
 To determine the model parameters that most strongly influence the optimal anneal time to 
maximize the SPE probability, we tested high and low values of the following four factors: 
damage diffusion constant (5E-6 μm2/min, 1E-5 μm2/min), damage length (50 nm, 100 nm), 
number of appropriate vacancy sites (5,10), and defect formation energy (-0.5 eV, -0.1 eV), with 
11 replicates for each of the 16 parameter combinations. Error was introduced for each replicate 
via the diffusion time parameter, where all calculations were performed over the time range 7E0 
min to 3E3 min for 1000 randomly spaced points and the random seed was different for each 
calculation. 
108 
  
 The response variables tested were time of maximum SPE probability (𝑡𝑡𝑏𝑏𝑏𝑏𝑒𝑒 𝑅𝑅𝑃𝑃𝐸𝐸), time of 
SPE onset (𝑡𝑡0 = 𝑡𝑡(𝑝𝑝𝑅𝑅𝑃𝑃𝐸𝐸 ≥ 0.01)), and color center saturation time (𝑡𝑡𝑠𝑠𝑏𝑏𝑡𝑡 = 𝑡𝑡max𝑅𝑅𝑃𝑃𝐸𝐸 − 𝑡𝑡0). 
Average values for each response variable for all 16 combinations are shown in Tables 4.4, 4.5, 
and 4.6 respectively. Single factor and interaction effects with corresponding p-values are shown 
for each of the response variables in Tables 4.7, 4.8, and 4.9, respectively. 
 
Table 4.4. Model Parameter 24 Results –Time of Max SPE Probability 
Experiment (𝐷𝐷𝐻𝐻𝐻𝐻, 𝐹𝐹, 𝑁𝑁𝐻𝐻𝑣𝑣 , ∆𝐸𝐸𝐻𝐻𝑓𝑓 ) 𝑡𝑡𝑏𝑏𝑏𝑏𝑒𝑒 𝑅𝑅𝑃𝑃𝐸𝐸 (Mean±Std Err) [min] 
1E-6 μm2/min, 50 nm, 10, -0.5 eV 26.47±0.22 
1E-6 μm2/min, 50 nm, 5, -0.5 eV 32.96±0.45 
5E-5 μm2/min, 50 nm, 10, -0.5 eV 52.44±0.28 
5E-5 μm2/min, 50 nm, 5, -0.5 eV 64.56±0.40 
1E-6 μm2/min, 100 nm, 10, -0.5 eV 85.94±0.58 
1E-6 μm2/min, 100 nm, 5, -0.5 eV 102.98±0.38 
1E-6 μm2/min, 50 nm, 10, -0.1 eV 136.29±0.25 
5E-5 μm2/min, 100 nm, 10, -0.5 eV 170.61±0.39 
5E-5 μm2/min, 100 nm, 5, -0.5 eV 204.82±0.41 
1E-6 μm2/min, 50 nm, 5, -0.1 eV 237.25±0.56 
5E-5 μm2/min, 50 nm, 10, -0.1 eV 272.41±0.31 
1E-6 μm2/min, 100 nm, 10, -0.1 eV 332.25±0.23 
5E-5 μm2/min, 50 nm, 5, -0.1 eV 473.42±0.45 
1E-6 μm2/min, 100 nm, 5, -0.1 eV 524.87±0.13 
5E-5 μm2/min, 100 nm, 10, -0.1 eV 665.23±0.44 
5E-5 μm2/min, 100 nm, 5, -0.1 eV 11050.08±0.22 
 
Table 4.5. Model Parameter 24 Results –SPE Onset Time 
Experiment (𝐷𝐷𝐻𝐻𝐻𝐻, 𝐹𝐹, 𝑁𝑁𝐻𝐻𝑣𝑣 , ∆𝐸𝐸𝐻𝐻𝑓𝑓 ) 𝑡𝑡0 (Mean±Std Err) [min] 
1E-6 μm2/min, 50 nm, 10, -0.5 eV 12.51±0.45 
1E-6 μm2/min, 50 nm, 5, -0.5 eV 13.87±0.43 
1E-6 μm2/min, 50 nm, 10, -0.1 eV 21.72±0.34 
5E-5 μm2/min, 50 nm, 10, -0.5 eV 22.63±0.29 
5E-5 μm2/min, 50 nm, 5, -0.5 eV 25.13±0.33 
1E-6 μm2/min, 50 nm, 5, -0.1 eV 25.45±0.26 
1E-6 μm2/min, 100 nm, 10, -0.5 eV 39.82±0.44 
5E-5 μm2/min, 50 nm, 10, -0.1 eV 44.05±0.70 
1E-6 μm2/min, 100 nm, 5, -0.5 eV 44.55±0.60 
5E-5 μm2/min, 50 nm, 5, -0.1 eV 50.84±0.46 
1E-6 μm2/min, 100 nm, 10, -0.1 eV 70.40±0.35 
5E-5 μm2/min, 100 nm, 10, -0.5 eV 78.05±0.50 
1E-6 μm2/min, 100 nm, 5, -0.1 eV 81.59±0.41 
5E-5 μm2/min, 100 nm, 5, -0.5 eV 86.42±0.70 
5E-5 μm2/min, 100 nm, 10, -0.1 eV 140.03±0.36 
5E-5 μm2/min, 100 nm, 5, -0.1 eV 162.12±0.33 
109 
  
 
Table 4.6. Model Parameter 24 Results –Color Center Saturation Time 
Experiment (𝐷𝐷 𝐻𝐻𝐻𝐻𝐻𝐻, 𝐹𝐹, 𝑁𝑁𝑣𝑣 , ∆𝐸𝐸𝐻𝐻𝑓𝑓 ) 𝑡𝑡𝑠𝑠𝑏𝑏𝑡𝑡 (Mean±Std Err) [min] 
1E-6 μm2/min, 50 nm, 10, -0.5 eV 13.97±0.45 
1E-6 μm2/min, 50 nm, 5, -0.5 eV 19.09±0.60 
5E-5 μm2/min, 50 nm, 10, -0.5 eV 29.81±0.44 
5E-5 μm2/min, 50 nm, 5, -0.5 eV 39.44±0.60 
1E-6 μm2/min, 100 nm, 10, -0.5 eV 46.12±0.80 
1E-6 μm2/min, 100 nm, 5, -0.5 eV 58.43±0.76 
5E-5 μm2/min, 100 nm, 10, -0.5 eV 92.56±0.73 
1E-6 μm2/min, 50 nm, 10, -0.1 eV 114.56±0.47 
5E-5 μm2/min, 100 nm, 5, -0.5 eV 118.40±0.88 
1E-6 μm2/min, 50 nm, 5, -0.1 eV 211.79±0.46 
5E-5 μm2/min, 50 nm, 10, -0.1 eV 228.36±0.59 
1E-6 μm2/min, 100 nm, 10, -0.1 eV 262.07±0.29 
5E-5 μm2/min, 50 nm, 5, -0.1 eV 422.58±0.56 
1E-6 μm2/min, 100 nm, 5, -0.1 eV 443.29±0.40 
5E-5 μm2/min, 100 nm, 10, -0.1 eV 525.20±0.37 
5E-5 μm2/min, 100 nm, 5, -0.1 eV 887.96±0.36 
 
Table 4.7. Factor Screening Full Model Effects Test for Time of Max SPE Probability 𝛼𝛼 = 
0.05 significance level 
Input Factor or Interaction Effect Size (F) p-value 
Defect formation energy 3848027 < 0.0001 
Damage length 1497734 < 0.0001 
Damage diffusion constant 960328 < 0.0001 
Damage length × Defect formation energy 501450 < 0.0001 
Damage diffusion constant × Defect formation energy 429663 < 0.0001 
Available vacancy sites 397954 < 0.0001 
Available vacancy sites × Defect formation energy 289418 < 0.0001 
Damage diffusion constant × Damage length 166888 < 0.0001 
Damage diffusion constant × Damage length × Defect formation energy 56214 < 0.0001 
Damage diffusion constant × Available vacancy sites 43927 < 0.0001 
Damage length × Available vacancy sites 41887 < 0.0001 
Damage diffusion constant × Available vacancy sites × Defect 32135 < 0.0001 
formation energy 
Damage length × Available vacancy sites × Defect formation energy 26010 < 0.0001 
Damage diffusion constant × Damage length × Available vacancy sites 4772 < 0.0001 
Damage diffusion constant × Damage length × Available vacancy sites 2889 < 0.0001 
× Defect formation energy 
 
 
 
 
 
 
110 
  
Table 4.8. Factor Screening Full Model Effects Test for SPE Onset Time 𝛼𝛼 = 0.05 
significance level 
Input Factor or Interaction Effect Size (F) p-value 
Damage length 72271 < 0.0001 
Damage diffusion constant 27329 < 0.0001 
Defect formation energy 22771 < 0.0001 
Damage diffusion constant × Damage length 7922 < 0.0001 
Damage length × Defect formation energy 5754 < 0.0001 
Damage diffusion constant × Defect formation energy 2834 < 0.0001 
Available vacancy sites 1126 < 0.0001 
Damage diffusion constant × Damage length × Defect formation energy 583 < 0.0001 
Damage length × Available vacancy sites 312 < 0.0001 
Available vacancy sites × Defect formation energy 220 < 0.0001 
Damage diffusion constant × Available vacancy sites 107 < 0.0001 
Damage length × Available vacancy sites × Defect formation energy 56 < 0.0001 
Damage diffusion constant × Damage length × Available vacancy sites 33 < 0.0001 
Damage length × Available vacancy sites × Defect formation energy 26 < 0.0001 
Damage diffusion constant × Damage length × Available vacancy sites 9 0.0036 
× Defect formation energy 
 
Table 4.9. Factor Screening Full Model Effects Test for Color Center Saturation Time 𝛼𝛼 = 
0.05 significance level 
Input Factor or Interaction Effect Size (F) p-value 
Defect formation energy 1363530 < 0.0001 
Damage length 348782 < 0.0001 
Damage diffusion constant 262489 < 0.0001 
Damage length × Defect formation energy 163738 < 0.0001 
Damage diffusion constant × Defect formation energy 150519 < 0.0001 
Available vacancy sites 150034 < 0.0001 
Available vacancy sites × Defect formation energy 116424 < 0.0001 
Damage diffusion constant × Damage length 39092 < 0.0001 
Damage diffusion constant × Damage length × Defect formation energy 18625 < 0.0001 
Damage diffusion constant × Available vacancy sites 16722 < 0.0001 
Damage length × Available vacancy sites 14476 < 0.0001 
Damage diffusion constant × Available vacancy sites × Defect 12903 < 0.0001 
formation energy 
Damage length × Available vacancy sites × Defect formation energy 9980 < 0.0001 
Damage diffusion constant × Damage length × Available vacancy sites 1665 < 0.0001 
Damage diffusion constant × Damage length × Available vacancy sites 1084 < 0.0001 
× Defect formation energy 
 
 
 
 
 
111 
  
Ab-initio Results 
 In order to emulate experimental process of FIB followed by CVD, DFT-based calculations 
were used to determine an upper estimate of defect formation energies for a representative set of 
carbon-based defect complexes, which consist of substitutional carbons adjacent to different 
vacancies. Figure 4.8(a) gives the formation energy of these defect complexes using two 
reference structures – pristine hBN and defective hBN with vacancies already present within the 
matrix. As seen in the histogram, the formation energy of carbon-containing defect complexes is 
much larger in pristine hBN than the energies required to incorporate carbon once vacancies 
have already been created in hBN. Additionally, the reported formation energies for carbon-
containing defect complexes in defective hBN are consistent with the estimations for defect 
formation energy used in the color center yield calculations reported in Figures 4.3-4.4 and 
Table 4.3. 
 Figure 4.8(b) and (c), for VBCN and VNCB defects, respectively, show the ranges of 
electronic chemical potential (Fermi level) over which their different charge states are stable. 
Figure 4.8(b) shows that the three most favorable charge states for VBCN are 0, -1 and -2, out of 
which only first two are spin active, with spins of 1 and ½, respectively and so, this defect might 
not be the one observed by Mendelson et. al in experiment.  Also, according to our DFT 
calculations, the neutral VNCB is not spin active, with the non-magnetic solution being the lower 
energy state, which was missed by earlier works94. Amongst the energetically most favorable 
charge states (+2, +1 and 0) shown in Figure 4.8(c), only +1 charge state is spin-active and has a 
spin-1/2, which was also reported previously95. 
112 
  
Figure 4.8. Density functional theory (DFT) calculations. (a) Defect formation energies for a 
representative set of defect complexes consisting of carbon, with two different starting reference 
structures (abbreviated as Ref. Str.). For sake of comparison, all defect formation energies are 
calculated for boron-rich conditions.  The large reduction in formation energies of carbon-based 
defects if the vacancies are already present in the hBN matrix (say, introduced via FIB method) 
as compared to creating these complexes in pristine hBN corroborate the experimental results. 
(b)  Defect formation energy for VBCN defect as a function of electronic chemical potential 
(Fermi energy). Fermi energy is measured relative to the VBM. The solid line indicates the most 
stable charge states for a given electronic chemical potential, showing that this defect prefers a 
charge state of 0, -1 and -2, depending on the Fermi energy. The ball and stick model in the inset 
shows the equilibrium geometry of the defect in the neutral state. (c) Defect formation energy for 
VNCB defect as a function of Fermi energy, showing the preferred charge states for any given 
Fermi level. Only +1 charge state of VNCB is spin active (spin-1/2). The inset shows the 
equilibrium geometry of the defect in the neutral state, wherein the weak boron-boron bond 
between the boron atoms surrounding the defect result in a structure without any local magnetic 
moment. 
 
 
 
 
113 
 
 
CHAPTER V 
 
CONCLUDING REMARKS 
 
 The primary work in this dissertation to pattern and characterize quantum emitters in 2D 
insulator hBN makes device integration for quantum information technologies more 
straightforward and adjustable. We first achieve this via repeatable success of SPE creation by 
FIB patterning and carbon annealing at multiple parameter combination values with control over 
ZPL and purity as reported in Chapter II. Additionally, by elucidating a physical mechanism for 
the defect formation process via AFM characterization in Chapter III and modeling in Chapter 
IV, we provide an avenue to better control the optimization of color center creation and 
photophysical properties.  
 This technique increases the utility of room temperature carbon-based hBN SPEs in photonic 
nanostructures such as resonant cavities, optomechanical transducers, and magnetic quantum 
sensors by allowing pick-and-place emitter integration rather than relying on multi-step 
alignment or scaffolding processes. Furthermore, our analysis and modeling improve this utility 
by characterizing the effects of FIB exposure and carbon annealing on emitter yields and 
properties so that future applications of this technique can mitigate undesired damage, reduce 
diffusion source byproducts, and control the size and type of defects. The increased level of 
control over SPE creation possible due to this work can be used to improve the stability and 
indistinguishability of hBN emitters, which are requirements for quantum applications. 
 Beyond hBN, many solid-state emission sources based on extrinsic defects benefit from 
adjustable patterning techniques. Color centers in three-dimensional wide bandgap materials 
1 14 
  
such as diamond, GaN, and SiC as well as defect-induced exciton transitions in 2D transition 
metal dichalcogenides (TMDCs) all depend on the controlled introduction of impurities to form 
defect complexes8,45, but current fabrication techniques either rely on niche implantation systems 
or need multi-step photolithography processes to localize emitters, and sometimes both. The 
general technique presented in this work depends only on the ability to first open a diffusion 
window and create vacancies via FIB, then introduce the impurity dopant via CVD or annealing, 
using equipment and resources more commonly found in standard fabrication facilities. Although 
other single photon sources are based on different extrinsic defect complexes (e.g. NV in 
diamond, V V  in 4H-SiC104Si C , Crw in WS 1052 ), adjusting the dopant feedstock gas or evaporating 
a precursor upstream can provide control over impurity sources and concentrations, while 
adjusting FIB patterning parameters can control the available vacancy concentration and the 
introduction of impurities via a high damage diffusion barrier. Applying this technique broadly 
to pattern solid-state single photon emitters can enable widespread emitter integration for 
quantum photonics and sensors that have unique advantages. 
 These results provide a scalable technique to create solid-state single photon sources in 
materials beyond hBN with low native defect densities, allowing exposure dose, impurity 
species, and diffusion time to be chosen based on system requirements and fabrication 
constraints. This fabrication and modeling framework can be used to design and integrate 
monolithic and hybrid quantum systems that rely on coupling single photons in quantum 
information technologies for transduction, communication, and sensing with more readily 
accessible resources and more time efficient engineering and optimization. 
115 
 
 
APPENDICIES 
 
APPENDIX A 
 
CORRELATIVE SCANNING NEAR-FIELD OPTICAL FLUORESCENCE  
AND TOPOGRAPHY OF HBN COLOR CENTERS 
 
This appendix summarizes preliminary results obtained during my work on the WITec alpha 
3000 Correlative Scanning Near-Field Optical Microscope (SNOM) housed in the Center for 
Advanced Materials Characterization in Oregon (CAMCOR) Surface Analytic Facility in order 
to adapt it for fluorescence measurements. These results were previously reported at the Oregon 
Center for Optical, Molecular, and Quantum Science (OMQ) Fall 2019 Symposium via poster 
presentation. 
 
 Techniques to engineer quantum emitters (QEs) in deterministic locations in solid-state 
hosts such as the 2D material hexagonal boron nitride (hBN) rely on ion irradiation or focused 
ion beam (FIB) patterning that can damage the material surface while creating the desired 
vacancy or implantation defects. Damage to the material and substrate can produce local strain 
gradients and charge traps, which can ultimately lead to variation and instabilities in the zero-
phonon line (ZPL) energy, fluorescence intensity, and quantum purity. Current techniques to 
characterize emitters and correlate material damage are limited by optical localization 
capabilities (i.e. the diffraction limit) and alignment calibrations when transferring samples 
between optical and topographical characterization tools. It would therefore be advantageous to 
1 16 
  
develop a technique to simultaneously characterize optical properties of emitters and the local 
topography of the host material at the defect location with resolutions beyond the diffraction 
limit.  
We leveraged the high-resolution imaging technique of scanning near-field optical 
fluorescence microscopy (SNOM-f) to investigate the location and correlated topography of 
activated emitters to elucidate intrinsic strain dependences on QE formation. In order to only 
collect off-resonant fluorescent color centers, we adapted the WITec alpha 3000 Correlative 
SNOM transmission collection optics by adding a 532 nm notch filter and 550 nm longpass filter 
[Figure A.1]. In this particular instrument, a 532 nm excitation laser is focused into a hollow 
cantilever with 100 nm aperture. The cantilever utilizes AFM feedback principles to maintain a 
sub-wavelength contact distance from the sample during a scan. The evanescent waves 
transmitted through the cantilever aperture are within near-field distance of defects to excite, 
thereby beating the diffraction limit. 
We tested the resolution of SNOM-f by characterizing the same region of an exfoliated 
hBN flake first in a home-built confocal microscope and then in the SNOM. We found that our 
SNOM-f technique improved the resolution by a factor of 2 [Figure A.2]. For the region labeled 
Spot 1, we were able to resolve two distinct color centers from the SNOM-f scan, while there 
was only one present in the corresponding confocal scan. We expect that this indicated there 
were multiple defect-based color centers present in the same localized region that were unable to 
be resolved via traditional confocal characterization. 
 
 
117 
  
Figure A.1. Schematic of WITec SNOM adapted for fluorescence (Inset) SNOM principles 
 
Figure A.2. Improving confocal resolution with SNOM (a) Confocal PL scan of exfoliated hBN 
flake and corresponding intensity line cut for spot 1 (b) SNOM-f scan of the boxed area and 
corresponding intensity line cut for spot 1 
 
118 
  
During optical collection, we obtained simultaneous AFM height topography data as 
shown in Figure A.3. We observed surface features in locations that corresponded to the 
fluorescent locations, which could indicate that the color centers are either strain activated or 
were due to some other fluorescent feature on the surface of hBN and not from defects within the 
hBN itself. 
Figure A.3. Simultaneous correlation of hBN color center SNOM fluorescence to AFM 
height topography (a) AFM height scan and corresponding line cut of marked area (b) SNOM-f 
scan of the same area taken simultaneously and corresponding intensity line cut on the same line. 
Highlighted areas and areas on all images correspond to correlated points 
To determine if the observed correlations between AFM topography and SNOM 
fluorescence were due to material dependent defects or were some other imaging artefact, we 
characterized another sample and identified one of the color centers as a QE [Figure A.4]. In 
particular, the highlighted defect appeared on a flake edge in a folded region, supporting 
previous work that observed that hBN SPEs prefer to form at edges. Furthermore, the folded area 
119 
  
Figure A.4. Correlative imaging of an hBN quantum emitter (a) SNOM fluorescence image of 
exfoliated hBN with SPE circled (b) AFM topography image taken simultaneously (Inset) 
2
Antibunching dip of the emitter circled in (a) with g (0) < 0.5 characteristic of SPEs 
may have introduced further local strain or confinement that allowed this color center to be an 
isolated SPE. 
Acquiring SNOM data on many SPEs could help answer questions about why SPEs tend 
to form at edges or near grain boundaries, which may be a result of intrinsic strain gradients 
created by deformations in the material. This technique could be used in future work to identify 
where SPEs tend to form and potential intrinsic strain, as the tunable range of hBN SPEs 
depends on intrinsic strain.106 By understanding under which conditions SPEs prefer to form, we 
can further tailor our milling parameters. 
 
 
120 
 
 
APPENDIX B 
 
SAMPLE IMAGES 
 
This appendix contains optical microscope, scanning electron microscope, and confocal 
photoluminescence images of the samples characterized for this dissertation.  
 
Figure B.1. Optical and SEM of factor screening samples 
1 21 
  
Figure B.2. Confocal PL scans of factor screening samples 
 
122 
  
 123 
Figure B.3. Optical and SEM of optimization samples 
  
 
 
 
 
124 
Figure B.4. Confocal PL of optimization samples 
 
 
APPENDIX C 
 
HOME-BUILT CONFOCAL MICROSCOPE USE INSTRUCTIONS 
 
This appendix outlines the steps to operate the home-built confocal microscope in the Alemán 
lab used to collect data for this dissertation. 
 
Confocal Scan 
• Load Sample 
o Use a thin strip of double sided capton tape to adhere sample chip to glass 
microscope slide 
o Place microscope slide on sample stage by sliding under clips 
• Make sure flip mirror is flipped down so that light is collected by camera 
• Turn on camera lamp (red box) 
• Turn on Laser (black box) 
o Flip switch on back of box 
o Wait until red light turns on, the turn key to horizontal orientation 
o Make sure ND filter is in the range 4.1-4.6 
• Adjust position and focus knobs until you can see the light reflect off the sample by eye 
• Open ThorCam Live recording on computer 
1 25 
  
o  
o If sample is not visible (typically blue even if unfocused), Auto-Scale Display 
Settings 
• Course focus by hand until laser spot is small and round 
• Find desired measurement location, adjusting focus along the way 
• To fine focus:  
o Make sure APD power supply strip is OFF 
o On Confocal 2 LabVIEW program click Start Counter 
o        
o Use Picometer on step size + or – 1 to “Step Z pico” (move the stage by small 
amounts), while watching the ThorCam for best focus 
126 
  
• When ready to start Confocal Scan (must be done in order): 
o Turn off camera lamp (red box) 
o Push button to flip mirror up (so light goes to collection optics) 
o Change ND filter to 1.6 
 If taking many measurements over multiple days, measure power output 
and adjust filter to keep constant excitation power 
o Turn off room lights  
o Decrease computer brightness 
 Right click LED display on bottom task bar (it looks like a blue computer 
screen) 
 Click Configure 
 Slide scale bar to 95% 
o Turn on APD power strip 
 Open up TimeHarp software 
 Watch sync/cps and input/cps count numbers. The APDs are ready once 
the count numbers flip back to 0 and then show counts again 
• Change file location to correct project folder with user initials  
o   
• Set scan size by scrolling through pre-set XY sizes (in microns) 
127 
  
o  
• Start Scan 
Photon Counting [𝒈𝒈𝟎𝟎(𝒕𝒕)] measurements 
• When Confocal Scan is completed, turn off monitor and move to Remote Desktop 
• click Start Counter 
• Drag center cursor over bright spots 
o  
 
• Hover for ~30 seconds to see if signal attenuates or if there is consistent fluorescence  
128 
  
Consistent fluorescence, 
bleaching good to try g2 
o  
• Open TimeHarp software 
• Click GO 
o  
 
129 
 
 
REFERENCES CITED 
 
1. Kimble, H. J. The quantum internet. Nature 453, 1023–1030 (2008). 
2. Wang, Y., Li, J., Zhang, S., Su, K., Zhou, Y., Liao, K., Du, S., Yan, H. & Zhu, S. L. 
Efficient quantum memory for single-photon polarization qubits. Nat. Photonics 13, 346–
351 (2019). 
3. Schell, A. W., Takashima, H., Tran, T. T., Aharonovich, I. & Takeuchi, S. Coupling 
Quantum Emitters in 2D Materials with Tapered Fibers. ACS Photonics 4, 761–767 
(2017). 
4. Gao, X., Yin, Z.-Q. & Li, T. High-Speed Quantum Transducer with a Single-Photon 
Emitter in a 2D Resonator. Ann. Phys. 532, 2000233 (2020). 
5. Hollenbach, M., Berencén, Y., Kentsch, U., Helm, M. & Astakhov, G. V. Engineering 
telecom single-photon emitters in silicon for scalable quantum photonics. Opt. Express 28, 
26111 (2020). 
6. Durand, A., Baron, Y., Redjem, W., Herzig, T., Benali, A., Pezzagna, S., Meijer, J., 
Kuznetsov, A. Y., Gérard, J. M., Robert-Philip, I., Abbarchi, M., Jacques, V., Cassabois, 
G. & Dréau, A. Broad Diversity of Near-Infrared Single-Photon Emitters in Silicon. Phys. 
Rev. Lett. 126, 083602 (2021). 
7. Prabhu, M., Errando-Herranz, C., De Santis, L., Christen, I., Chen, C., Gerlach, C. & 
Englund, D. Individually addressable and spectrally programmable artificial atoms in 
silicon photonics. Nat. Commun. 14, 2380 (2023). 
8. Moody, G., Sorger, V. J., Blumenthal, D. J., Juodawlkis, P. W., Loh, W., Sorace-agaskar, 
C., Jones, A. E., Balram, K. C., Matthews, J. C. F., Laing, A., Davanco, M., Chang, L., 
Bowers, J. E., Quack, N., Galland, C., Aharonovich, I., Wolff, M. A., Schuck, C., Sinclair, 
N., Lon, M., Komljenovic, T., Weld, D., Mookherjea, S., Buckley, S., Radulaski, M., 
Reitzenstein, S., Pingault, B., Machielse, B. & Mukhopadhyay, D. 2022 Roadmap on 
integrated quantum photonics. J. Phys. Photonics 4, 012501 (2022). 
9. Kaur, P., Boes, A., Ren, G., Nguyen, T. G., Roelkens, G. & Mitchell, A. Hybrid and 
heterogeneous photonic integration. APL Photonics 6, 061102 (2021). 
10. Yeo, I., De Assis, P. L., Gloppe, A., Dupont-Ferrier, E., Verlot, P., Malik, N. S., Dupuy, 
E., Claudon, J., Gérard, J. M., Auffèves, A., Nogues, G., Seidelin, S., Poizat, J. P., Arcizet, 
O. & Richard, M. Strain-mediated coupling in a quantum dot-mechanical oscillator hybrid 
system. Nat. Nanotechnol. 9, 106–110 (2014). 
11. Golter, D. A., Oo, T., Amezcua, M., Stewart, K. A. & Wang, H. Optomechanical 
Quantum Control of a Nitrogen-Vacancy Center in Diamond. Phys. Rev. Lett. 116, 
143602 (2016). 
1 30 
  
12. Lee, K. W., Lee, D., Ovartchaiyapong, P., Minguzzi, J., Maze, J. R. & Bleszynski Jayich, 
A. C. Strain Coupling of a Mechanical Resonator to a Single Quantum Emitter in 
Diamond. Phys. Rev. Appl. 6, 034005 (2016). 
13. Castellanos-Gomez, A., Singh, V., Van Der Zant, H. S. J. & Steele, G. A. Mechanics of 
freely-suspended ultrathin layered materials. Ann. Phys. 527, 27–44 (2015). 
14. Tran, T. T., Bray, K., Ford, M. J., Toth, M. & Aharonovich, I. Quantum emission from 
hexagonal boron nitride monolayers. Nat. Nanotechnol. 11, 37–41 (2015). 
15. Tran, T. T., Elbadawi, C., Totonjian, D., Lobo, C. J., Grosso, G., Moon, H., Englund, D. 
R., Ford, M. J., Aharonovich, I. & Toth, M. Robust Multicolor Single Photon Emission 
from Point Defects in Hexagonal Boron Nitride. ACS Nano 10, 7331–7338 (2016). 
16. Kianinia, M., Regan, B., Tawfik, S. A., Tran, T. T., Ford, M. J., Aharonovich, I. & Toth, 
M. Robust Solid-State Quantum System Operating at 800 K. ACS Photonics 4, 768–773 
(2017). 
17. Proscia, N. V, Hotan, Z. A. V. S., Ayakumar, H. A. J., Eddy, P. R. R., Ohen, C. H. C., 
Ollar, M. I. D., Lkauskas, A. U. A., Oherty, M. A. D., Eriles, C. A. A. M. & Enon, V. I. 
M. M. Near-deterministic activation of room-temperature quantum emitters in hexagonal 
boron nitride. Optica 5, 1128–1134 (2018). 
18. Mendelson, N., Xu, Z. Q., Tran, T. T., Kianinia, M., Scott, J., Bradac, C., Aharonovich, I. 
& Toth, M. Engineering and Tuning of Quantum Emitters in Few-Layer Hexagonal Boron 
Nitride. ACS Nano 13, 3132–3140 (2019). 
19. Tawfik, S. A., Ali, S., Fronzi, M., Kianinia, M., Tran, T. T., Stampfl, C., Aharonovich, I., 
Toth, M. & Ford, M. J. First-principles investigation of quantum emission from hBN 
defects. Nanoscale 9, 13575–13582 (2017). 
20. Kozawa, D., Li, S. X., Ichihara, T., Rajan, A. G., Gong, X., He, G., Koman, V. B., Zeng, 
Y., Kuehne, M., Silmore, K. S., Parviz, D., Liu, P., Liu, A. T., Faucher, S., Yuan, Z., 
Warner, J., Blankschtein, D. & Strano, M. S. Discretized hexagonal boron nitride quantum 
emitters and their chemical interconversion. Nanotechnology 34, 115702 (2023). 
21. Narayanan, S. K. & Dev, P. Substrate-Induced Modulation of Quantum Emitter Properties 
in 2D Hexagonal Boron Nitride: Implications for Defect-Based Single Photon Sources in 
2D Layers. ACS Appl. Nano Mater. 6, 3446–3452 (2023). 
22. Pelini, T., Waselowski, V., Maze, J. R., Gil, B., Cassabois, G. & Jacques, V. Efficient 
single photon emission from a high-purity hexagonal boron nitride crystal. 121405, 1–5 
(2016). 
23. Ziegler, J., Klaiss, R., Blaikie, A., Miller, D., Horowitz, V. R. & Alemán, B. J. 
Deterministic Quantum Emitter Formation in Hexagonal Boron Nitride via Controlled 
Edge Creation. Nano Lett. 19, 2121–2127 (2019). 
131 
  
24. Stern, H. L., Wang, R., Fan, Y., Mizuta, R., Stewart, J. C., Needham, L., Roberts, T. D., 
Wai, R., Ginsberg, N. S., Klenerman, D., Hofmann, S. & Lee, S. F. Spectrally Resolved 
Photodynamics of Individual Emitters in Large-Area Monolayers of Hexagonal Boron 
Nitride. ACS Nano 13, 4538–4547 (2019). 
25. Duong, N. M. H., Nguyen, M. A. P., Kianinia, M., Ohshima, T., Abe, H., Watanabe, K., 
Taniguchi, T., Edgar, J. H., Aharonovich, I. & Toth, M. Effects of High-Energy Electron 
Irradiation on Quantum Emitters in Hexagonal Boron Nitride. ACS Appl. Mater. 
Interfaces 10, 24886–24891 (2018). 
26. Frantsuzov, P., Kuno, M., Jankó, B. & Marcus, R. A. Universal emission intermittency in 
quantum dots, nanorods and nanowires. Nat. Phys. 4, 519–522 (2008). 
27. Kuno, M., Fromm, D. P., Hamann, H. F., Gallagher, A. & Nesbitt, D. J. Nonexponential 
“blinking” kinetics of single CdSe quantum dots: A universal power law behavior. J. 
Chem. Phys. 112, 3117–3120 (2000). 
28. Shotan, Z., Jayakumar, H., Considine, C. R., Mackoit, M., Fedder, H., Wrachtrup, J., 
Alkauskas, A., Doherty, M. W., Menon, V. M. & Meriles, C. A. Photoinduced 
Modification of Single-Photon Emitters in Hexagonal Boron Nitride. ACS Photonics 3, 
2490–2496 (2016). 
29. Gottscholl, A., Kianinia, M., Soltamov, V., Orlinskii, S., Mamin, G., Bradac, C., Kasper, 
C., Krambrock, K., Sperlich, A., Toth, M., Aharonovich, I. & Dyakonov, V. Initialization 
and read-out of intrinsic spin defects in a van der Waals crystal at room temperature. Nat. 
Mater. 19, 540–545 (2020). 
30. Healey, A. J., Scholten, S. C., Yang, T., Scott, J. A., Abrahams, G. J., Robertson, I. O., 
Hou, X. F., Guo, Y. F., Rahman, S., Lu, Y., Kianinia, M., Aharonovich, I. & Tetienne, J. 
P. Quantum microscopy with van der Waals heterostructures. Nat. Phys. 19, 87–91 
(2023). 
31. Mendelson, N., Doherty, M., Toth, M., Aharonovich, I. & Tran, T. T. Strain-Induced 
Modification of the Optical Characteristics of Quantum Emitters in Hexagonal Boron 
Nitride. Adv. Mater. 32, 1908316 (2020). 
32. Nikolay, N., Mendelson, N., Sadzak, N., Böhm, F., Tran, T. T., Sontheimer, B., 
Aharonovich, I. & Benson, O. Very Large and Reversible Stark-Shift Tuning of Single 
Emitters in Layered Hexagonal Boron Nitride. Phys. Rev. Appl. 10, 041001 (2019). 
33. Auburger, P. & Gali, A. Towards ab initio identification of paramagnetic substitutional 
carbon defects in hexagonal boron nitride acting as quantum bits. Phys. Rev. B 104, 
075410 (2021). 
34. Chejanovsky, N., Mukherjee, A., Geng, J., Chen, Y. C., Kim, Y., Denisenko, A., Finkler, 
A., Taniguchi, T., Watanabe, K., Dasari, D. B. R., Auburger, P., Gali, A., Smet, J. H. & 
Wrachtrup, J. Single-spin resonance in a van der Waals embedded paramagnetic defect. 
Nat. Mater. 20, 1079–1084 (2021). 
132 
  
35. Lee, H. Y., Al Ezzi, M. M., Raghuvanshi, N., Chung, J. Y., Watanabe, K., Taniguchi, T., 
Garaj, S., Adam, S. & Gradečak, S. Tunable Optical Properties of Thin Films Controlled 
by the Interface Twist Angle. Nano Lett. 21, 2832–2839 (2021). 
36. Su, C., Zhang, F., Kahn, S., Shevitski, B., Jiang, J., Dai, C., Ungar, A., Park, J. H., 
Watanabe, K., Taniguchi, T., Kong, J., Tang, Z., Zhang, W., Wang, F., Crommie, M., 
Louie, S. G., Aloni, S. & Zettl, A. Tuning colour centres at a twisted hexagonal boron 
nitride interface. Nat. Mater. 21, 896–902 (2022). 
37. Fischer, M., Caridad, J. M., Sajid, A., Ghaderzadeh, S., Ghorbani-Asl, M., Gammelgaard, 
L., Bøggild, P., Thygesen, K. S., Krasheninnikov, A. V., Xiao, S., Wubs, M. & Stenger, 
N. Controlled generation of luminescent centers in hexagonal boron nitride by irradiation 
engineering. Sci. Adv. 7, 7138 (2021). 
38. Kianinia, M., White, S., Fröch, J. E., Bradac, C. & Aharonovich, I. Generation of Spin 
Defects in Hexagonal Boron Nitride. ACS Photonics 7, 2147–2152 (2020). 
39. Xu, Z. Q., Elbadawi, C., Tran, T. T., Kianinia, M., Li, X., Liu, D., Hoffman, T. B., 
Nguyen, M., Kim, S., Edgar, J. H., Wu, X., Song, L., Ali, S., Ford, M., Toth, M. & 
Aharonovich, I. Single photon emission from plasma treated 2D hexagonal boron nitride. 
Nanoscale 10, 7957–7965 (2018). 
40. Mendelson, N., Chugh, D., Reimers, J. R., Cheng, T. S., Gottscholl, A., Long, H., Mellor, 
C. J., Zettl, A., Dyakonov, V., Beton, P. H., Novikov, S. V., Jagadish, C., Tan, H. H., 
Ford, M. J., Toth, M., Bradac, C. & Aharonovich, I. Identifying carbon as the source of 
visible single-photon emission from hexagonal boron nitride. Nat. Mater. 20, 321–328 
(2021). 
41. Lyu, C., Zhu, Y., Gu, P., Qiao, J., Watanabe, K., Taniguchi, T. & Ye, Y. Single-photon 
emission from two-dimensional hexagonal boron nitride annealed in a carbon-rich 
environment. Appl. Phys. Lett. 117, 244002 (2020). 
42. Horder, J., White, S., Gale, A., Li, C., Watanabe, K., Taniguchi, T., Kianinia, M., 
Aharonovich, I. & Toth, M. Coherence properties of electron beam activated emitters in 
hexagonal boron nitride under resonant excitation. Phys. Rev. Appl. 18, 064021 (2022). 
43. Xu, X., Martin, Z. O., Sychev, D., Lagutchev, A. S., Chen, Y. P., Taniguchi, T., 
Watanabe, K., Shalaev, V. M. & Boltasseva, A. Creating Quantum Emitters in Hexagonal 
Boron Nitride Deterministically on Chip-Compatible Substrates. Nano Lett. 21, 8182–
8189 (2021). 
44. Zhang, G., Cheng, Y., Chou, J. P. & Gali, A. Material platforms for defect qubits and 
single-photon emitters. Appl. Phys. Rev. 7, 031308 (2020). 
45. Aharonovich, I., Englund, D. & Toth, M. Solid-state single-photon emitters. Nat. 
Photonics 10, 631–641 (2016). 
 
133 
  
46. Kianinia, M., Xu, Z. Q., Toth, M. & Aharonovich, I. Quantum emitters in 2D materials: 
Emitter engineering, photophysics, and integration in photonic nanostructures. Appl. Phys. 
Rev. 9, 011306 (2022). 
47. Tran, T. T., Zachreson, C., Berhane, A. M., Bray, K., Sandstrom, R. G., Li, L. H., 
Taniguchi, T., Watanabe, K., Aharonovich, I. & Toth, M. Quantum Emission from 
Defects in Single-Crystalline Hexagonal Boron Nitride. Phys. Rev. Appl. 5, 034005 
(2016). 
48. Aharonovich, I., Tetienne, J. P. & Toth, M. Quantum Emitters in Hexagonal Boron 
Nitride. Nano Lett. 22, 9227–9235 (2022). 
49. Yang, T., Mendelson, N., Li, C., Gottscholl, A., Scott, J., Kianinia, M., Dyakonov, V., 
Toth, M. & Aharonovich, I. Spin defects in hexagonal boron nitride for strain sensing on 
nanopillar arrays. Nanoscale 14, 5239–5244 (2022). 
50. Gu, R., Wang, L., Zhu, H., Han, S., Bai, Y., Zhang, X., Li, B., Qin, C., Liu, J., Guo, G., 
Shan, X., Xiong, G., Gao, J., He, C., Han, Z., Liu, X. & Zhao, F. Engineering and 
Microscopic Mechanism of Quantum Emitters Induced by Heavy Ions in hBN. ACS 
Photonics 8, 2912–2922 (2021). 
51. Klaiss, R., Ziegler, J., Miller, D., Zappitelli, K., Watanabe, K., Taniguchi, T. & Alemán, 
B. Uncovering the morphological effects of high-energy Ga+focused ion beam milling on 
hBN single-photon emitter fabrication. J. Chem. Phys. 157, 074703 (2022). 
52. Weston, L., Wickramaratne, D., Mackoit, M., Alkauskas, A. & Van De Walle, C. G. 
Native point defects and impurities in hexagonal boron nitride. Phys. Rev. B 97, 214104 
(2018). 
53. Chu, P. K. & Li, L. Characterization of amorphous and nanocrystalline carbon films. 
Mater. Chem. Phys. 96, 253–277 (2006). 
54. Blundo, E., Surrente, A., Spirito, D., Pettinari, G., Yildirim, T., Chavarin, C. A., 
Baldassarre, L., Felici, M. & Polimeni, A. Vibrational Properties in Highly Strained 
Hexagonal Boron Nitride Bubbles. Nano Lett. 22, 1525–1533 (2022). 
55. Chen, L., Elibol, K., Cai, H., Jiang, C., Shi, W., Chen, C., Wang, H. S., Wang, X., Mu, X., 
Li, C., Watanabe, K., Taniguchi, T., Guo, Y., Meyer, J. C. & Wang, H. Direct observation 
of layer-stacking and oriented wrinkles in multilayer hexagonal boron nitride. 2D Mater. 
8, 024001 (2021). 
56. Ghaderzadeh, S., Kretschmer, S., Ghorbani-Asl, M., Hlawacek, G. & Krasheninnikov, A. 
V. Atomistic simulations of defect production in monolayer and bulk hexagonal boron 
nitride under low-and high-fluence ion irradiation. Nanomaterials 11, 1214 (2021). 
 
 
134 
  
57. Zeng, W., Nie, W., Li, L. & Chen, A. Ground-state cooling of a mechanical oscillator in a 
hybrid optomechanical system including an atomic ensemble. Sci. Rep. 7, 17258 (2017). 
58. Midolo, L., Schliesser, A. & Fiore, A. Nano-opto-electro-mechanical systems. Nat. 
Nanotechnol. 13, 11–18 (2018). 
59. Noh, G., Choi, D., Kim, J. H., Im, D. G., Kim, Y. H., Seo, H. & Lee, J. Stark Tuning of 
Single-Photon Emitters in Hexagonal Boron Nitride. Nano Lett. 18, 4710–4715 (2018). 
60. Chakraborty, C., Goodfellow, K. M., Dhara, S., Yoshimura, A., Meunier, V. & 
Vamivakas, A. N. Quantum-Confined Stark Effect of Individual Defects in a van der 
Waals Heterostructure. Nano Lett. 17, 2253–2258 (2017). 
61. Teissier, J., Barfuss, A., Appel, P., Neu, E. & Maletinsky, P. Strain Coupling of a 
Nitrogen-Vacancy Center Spin to a Diamond Mechanical Oscillator. Phys. Rev. Lett. 113, 
020503 (2014). 
62. Wang, L., Zihlmann, S., Baumgartner, A., Overbeck, J., Watanabe, K., Taniguchi, T., 
Makk, P. & Schönenberger, C. In Situ Strain Tuning in hBN-Encapsulated Graphene 
Electronic Devices. Nano Lett. 19, 4097–4102 (2019). 
63. Dev, P. Fingerprinting quantum emitters in hexagonal boron nitride using strain. Phys. 
Rev. Res. 2, 22050 (2020). 
64. Gao, X., Jiang, B., Llacsahuanga Allcca, A. E., Shen, K., Sadi, M. A., Solanki, A. B., Ju, 
P., Xu, Z., Upadhyaya, P., Chen, Y. P., Bhave, S. A. & Li, T. High-Contrast Plasmonic-
Enhanced Shallow Spin Defects in Hexagonal Boron Nitride for Quantum Sensing. Nano 
Lett. 21, 7708–7714 (2021). 
65. Stern, H. L., Gu, Q., Jarman, J., Eizagirre Barker, S., Mendelson, N., Chugh, D., Schott, 
S., Tan, H. H., Sirringhaus, H., Aharonovich, I. & Atatüre, M. Room-temperature 
optically detected magnetic resonance of single defects in hexagonal boron nitride. Nat. 
Commun. 13, 618 (2022). 
66. Caldwell, J. D., Aharonovich, I., Cassabois, G., Edgar, J. H., Gil, B. & Basov, D. N. 
Photonics with hexagonal boron nitride. Nat. Rev. Mater. 4, 552–567 (2019). 
67. Alaloul, M. & Khurgin, J. B. High-Performance All-Optical Modulator Based on 
Graphene-hBN Heterostructures. IEEE J. Sel. Top. Quantum Electron. 28, 3500108 
(2022). 
68. Lien, D. H., Uddin, S. Z., Yeh, M., Amani, M., Kim, H., Ager, J. W., Yablonovitch, E. & 
Javey, A. Electrical suppression of all nonradiative recombination pathways in monolayer 
semiconductors. Science 364, 468–471 (2019). 
 
 
135 
  
69. Klein, J., Lorke, M., Florian, M., Sigger, F., Sigl, L., Rey, S., Wierzbowski, J., Cerne, J., 
Müller, K., Mitterreiter, E., Zimmermann, P., Taniguchi, T., Watanabe, K., Wurstbauer, 
U., Kaniber, M., Knap, M., Schmidt, R., Finley, J. J. & Holleitner, A. W. Site-selectively 
generated photon emitters in monolayer MoS2 via local helium ion irradiation. Nat. 
Commun. 10, 2755 (2019). 
70. Kim, S., Fröch, J. E., Christian, J., Straw, M., Bishop, J., Totonjian, D., Watanabe, K., 
Taniguchi, T., Toth, M. & Aharonovich, I. Photonic Crystal Cavities from Hexagonal 
Boron Nitride. Nat. Commun. 9, 2623 (2018). 
71. Li, C., Fröch, J. E., Nonahal, M., Tran, T. N., Toth, M., Kim, S. & Aharonovich, I. 
Integration of hBN Quantum Emitters in Monolithically Fabricated Waveguides. ACS 
Photonics 8, 2966–2972 (2021). 
72. Fröch, J. E., Spencer, L. P., Kianinia, M., Totonjian, D. D., Nguyen, M., Gottscholl, A., 
Dyakonov, V., Toth, M., Kim, S. & Aharonovich, I. Coupling Spin Defects in Hexagonal 
Boron Nitride to Monolithic Bullseye Cavities. Nano Lett. 21, 6549–6555 (2021). 
73. Glushkov, E., Macha, M., Räth, E., Navikas, V., Ronceray, N., Cheon, C. Y., Ahmed, A., 
Avsar, A., Watanabe, K., Taniguchi, T., Shorubalko, I., Kis, A., Fantner, G. & Radenovic, 
A. Engineering Optically Active Defects in Hexagonal Boron Nitride Using Focused Ion 
Beam and Water. ACS Nano 16, 3695–3703 (2022). 
74. Choi, S., Tran, T. T., Elbadawi, C., Lobo, C., Wang, X., Juodkazis, S., Seniutinas, G., 
Toth, M. & Aharonovich, I. Engineering and Localization of Quantum Emitters in Large 
Hexagonal Boron Nitride Layers. ACS Appl. Mater. Interfaces 8, 29642–29648 (2016). 
75. Gale, A., Li, C., Chen, Y., Watanabe, K., Taniguchi, T., Aharonovich, I. & Toth, M. Site-
Specific Fabrication of Blue Quantum Emitters in Hexagonal Boron Nitride. ACS 
Photonics 9, 2170–2177 (2022). 
76. Li, C., Mendelson, N., Ritika, R., Chen, Y., Xu, Z. Q., Toth, M. & Aharonovich, I. 
Scalable and Deterministic Fabrication of Quantum Emitter Arrays from Hexagonal 
Boron Nitride. Nano Lett. 21, 3626–3632 (2021). 
77. Fournier, C., Plaud, A., Roux, S., Pierret, A., Rosticher, M., Watanabe, K., Taniguchi, T., 
Buil, S., Quélin, X., Barjon, J., Hermier, J. P. & Delteil, A. Position-controlled quantum 
emitters with reproducible emission wavelength in hexagonal boron nitride. Nat. 
Commun. 12, 3779 (2021). 
78. Gilbert, S. M., Liu, S., Schumm, G. & Zettl, A. Nanopatterning Hexagonal Boron Nitride 
with Helium Ion Milling: Towards Atomically-thin, Nanostructured Insulators. MRS Adv. 
3, 327–331 (2018). 
79. Proscia, N. Control of Light-Matter Interactions via Nanostructured Photonic Materials. 
CUNY Acad. Work. (2019). 
 
136 
  
80. Lassaline, N., Thureja, D., Chervy, T., Petter, D., Murthy, P. A., Knoll, A. W. & Norris, 
D. J. Freeform Electronic and Photonic Landscapes in Hexagonal Boron Nitride. Nano 
Lett. 21, 8175–8181 (2021). 
81. Ziegler, J. F., Ziegler, M. D. & Biersack, J. P. SRIM - The stopping and range of ions in 
matter. Nucl. Instruments Methods Phys. Res. Sect. B Beam Interact. with Mater. Atoms 
268, 1818–1823 (2010). 
82. Vladov, N., Segal, J. & Ratchev, S. Apparent beam size definition of focused ion beams 
based on scanning electron microscopy images of nanodots. J. Vac. Sci. Technol. B, 
Nanotechnol. Microelectron. Mater. Process. Meas. Phenom. 33, 041803 (2015). 
83. Gupta, S., Periasamy, P. & Narayanan, B. Defect dynamics in two-dimensional black 
phosphorus under argon ion irradiation. Nanoscale 13, 8575–8590 (2021). 
84. Gao, Y., Zhang, Y., Schwen, D., Jiang, C., Sun, C. & Gan, J. Formation and self-
organization of void superlattices under irradiation: A phase field study. Materialia 1, 78–
88 (2018). 
85. Gilbert, S. M., Dunn, G., Azizi, A., Pham, T., Shevitski, B., Dimitrov, E., Liu, S., Aloni, 
S. & Zettl, A. Fabrication of Subnanometer-Precision Nanopores in Hexagonal Boron 
Nitride. Sci. Rep. 7, 15096 (2017). 
86. Mouhoub, O., Martinez-Gordillo, R., Nelayah, J., Wang, G., Park, J. H., Kim, K. K., Lee, 
Y. H., Bichara, C., Loiseau, A., Ricolleau, C., Amara, H. & Alloyeau, D. Quantitative 
insights into the growth mechanisms of nanopores in hexagonal boron nitride. Phys. Rev. 
Mater. 4, 014005 (2020). 
87. Pham, T., Gibb, A. L., Li, Z., Gilbert, S. M., Song, C., Louie, S. G. & Zettl, A. Formation 
and Dynamics of Electron-Irradiation-Induced Defects in Hexagonal Boron Nitride at 
Elevated Temperatures. Nano Lett. 16, 7142–7147 (2016). 
88. Paszkowicz, W., Pelka, J. B., Knapp, M., Szyszko, T. & Podsiadlo, S. Lattice parameters 
and anisotropic thermal expansion of hexagonal boron nitride in the 10-297.5 K 
temperature range. Appl. Phys. A Mater. Sci. Process. 75, 431–435 (2002). 
89. Millett, P. C., Rokkam, S., El-Azab, A., Tonks, M. & Wolf, D. Void nucleation and 
growth in irradiated polycrystalline metals: A phase-field model. Model. Simul. Mater. 
Sci. Eng. 17, 064003 (2009). 
90. Weston, L., Wickramaratne, D., Mackoit, M., Alkauskas, A. & Van De Walle, C. G. 
Native point defects and impurities in hexagonal boron nitride. Phys. Rev. B 97, 214104 
(2018). 
91. Abdi, M., Hwang, M., Aghtar, M. & Plenio, M. B. Spin-Mechanical Scheme with Color 
Centers in Hexagonal Boron Nitride Membranes. Phys. Rev. Lett. 119, 233602 (2017). 
 
137 
  
92. Abdi, M. & Plenio, M. B. Quantum Effects in a Mechanically Modulated Single-Photon 
Emitter. Phys. Rev. Lett. 122, 23602 (2019). 
93. Aradi, E., Naidoo, S. R., Cummings, F., Motochi, I. & Derry, T. E. Cross-sectional 
transmission electron microscopy studies of boron ion implantation in hexagonal boron 
nitride. Diam. Relat. Mater. 92, 168–173 (2019). 
94. Tawfik, S. A., Ali, S., Fronzi, M., Kianinia, M., Tran, T. T., Stampfl, C., Aharonovich, I., 
Toth, M. & Ford, M. J. First-principles investigation of quantum emission from hBN 
defects. Nanoscale 9, 13575–13582 (2017). 
95. Abdi, M., Chou, J. P., Gali, A. & Plenio, M. B. Color Centers in Hexagonal Boron Nitride 
Monolayers: A Group Theory and Ab Initio Analysis. ACS Photonics 5, 1967–1976 
(2018). 
96. Li, S., Pershin, A., Thiering, G., Udvarhelyi, P. & Gali, A. Ultraviolet quantum emitters in 
hexagonal boron nitride from carbon clusters. J. Phys. Chem. Lett. 13, 3150–3157 (2022). 
97. Li, K., Smart, T. J. & Ping, Y. Carbon trimer as a 2 eV single-photon emitter candidate in 
hexagonal boron nitride: A first-principles study. Phys. Rev. Mater. 6, L042201 (2022). 
98. Carr, E. J. & Turner, I. W. A semi-analytical solution for multilayer diffusion in a 
composite medium consisting of a large number of layers. Appl. Math. Model. 40, 7034–
7050 (2016). 
99. Perdew, J. P. & Yue, W. Accurate and simple density functional for the electronic 
exchange energy: Generalized gradient approximation. Phys. Rev. B 33, 8800–8802 
(1986). 
100. Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized Gradient Approximation Made 
Simple. Phys. Rev. Lett. 77, 3865–3868 (1996). 
101. Giannozzi, P., Baroni, S., Bonini, N., Calandra, M., Car, R., Cavazzoni, C., Ceresoli, D., 
Chiarotti, G. L., Cococcioni, M., Dabo, I., Corso, A. D., Gironcoli, S. De, Fabris, S., 
Fratesi, G., Gebauer, R., Gerstmann, U., Gougoussis, C., Kokalj, A., Lazzeri, M., Martin-
Samos, L., Marzari, N., Mauri, F., Mazzarello, R., Paolini, S., Pasquarello, A., Paulatto, 
L., Sbraccia, C., Scandolo, S., Sclauzero, G., Seitsonen, A. P., Smogunov, A., Umari, P. & 
Wentzcovitch, R. M. QUANTUM ESPRESSO: a modular and open-source software 
project for quantum simulations of materials. J. Phys. Condens. Matter 21, 395502 (2009). 
102. Monkhorst, H. J. & Pack, J. D. Special points for Brillouin-zone integrations. Phys. Rev. B 
13, 5188–5192 (1976). 
103. Tanaka, T., Matsunaga, K., Ikuhara, Y. & Yamamoto, T. First-principles study on 
structures and energetics of intrinsic vacancies in SrTiO 3. Phys. Rev. B 68, 205213 
(2003). 
 
138 
  
104. Li, Q., Wang, J.-F., Yan, F.-F., Zhou, J.-Y., Wang, H.-F., Liu, H., Guo, L.-P., Zhou, X., 
Gali, A., Liu, Z.-H., Wang, Z.-Q., Sun, K., Guo, G.-P., Tang, J.-S., Li, H., You, L.-X., Xu, 
J.-S., Li, C.-F. & Guo, G.-C. Room-temperature coherent manipulation of single-spin 
qubits in silicon carbide with a high readout contrast. Natl. Sci. Rev. 9, (2022). 
105. Schuler, B., Lee, J. H., Kastl, C., Cochrane, K. A., Chen, C. T., Refaely-Abramson, S., 
Yuan, S., Van Veen, E., Roldán, R., Borys, N. J., Koch, R. J., Aloni, S., Schwartzberg, A. 
M., Ogletree, D. F., Neaton, J. B. & Weber-Bargioni, A. How Substitutional Point Defects 
in Two-Dimensional WS2 Induce Charge Localization, Spin-Orbit Splitting, and Strain. 
ACS Nano 13, 10520–10534 (2019). 
106. Grosso, G., Moon, H., Lienhard, B., Ali, S., Efetov, D. K., Furchi, M. M., Jarillo-Herrero, 
P., Ford, M. J., Aharonovich, I. & Englund, D. Tunable and high-purity room temperature 
single-photon emission from atomic defects in hexagonal boron nitride. Nat. Commun. 8, 
705 (2017). 
 
139