A Python library for simulating optically detected magnetic resonance (ODMR) experiments with spin-defect centres, herein referred to as SpinDefectSim. Includes a general spin Hamiltonian for arbitrary defect species, inhomogeneous ensemble modelling, AC sensing protocols, screened Coulomb electrostatics, magnetometry and electrometry from 2-D field distributions.
Physical models: all formulae, assumptions, and literature references are documented in
docs/theory.md.
pip install numpy scipy matplotlib # clone or copy SpinDefectSim/ into your project, then: import sys; sys.path.insert(0, "/path/to/parent") import SpinDefectSim
The library ships with the following built-in presets, selectable by name string. All can be overridden per-parameter.
| Name | Host | Spin | D0 (GHz) | d⊥ (Hz V−1m) | Notes |
|---|---|---|---|---|---|
vb_minus |
hBN | 1 | 3.46 | 0.35 | Default. VB−, Gottscholl 2020 |
nv_minus |
diamond | 1 | 2.87 | 0.17 | NV− centre, Dolde 2011 |
v_sic |
4H-SiC | 1 | 1.28 | 0.10 | V2 silicon vacancy |
p1 |
diamond | 1/2 | 0 | 0 | P1 N centre, Zeeman-only |
cr_gaN |
GaN | 3/2 | 1.80 | 0 | Cr4+ impurity, approximate |
Custom defect types can be created directly:
from SpinDefectSim.spin.defects import DefectType my_defect = DefectType(name="my_centre", spin=1, D0_Hz=1.5e9, d_perp=0.2)
For a spin-S defect with quantization axis z′, the Hamiltonian (H/h, in Hz) is:
H = D0·(Sz′2 − S(S+1)/3·I) ← ZFS axial (S ≥ 1)
+ E0·(Sx′2 − Sy′2) ← ZFS strain (S ≥ 1)
+ d∥·Ez′·(Sz′2 − S(S+1)/3·I) ← E axial (S ≥ 1)
+ d⊥·[Ey′·(Sx′2−Sy′2) + Ex′·{Sx′,Sy′}] ← E transverse (S ≥ 1)
+ γe·(Bx′·Sx′ + By′·Sy′ + Bz′·Sz′) ← electron Zeeman (all S)
+ Σk S·Ak·Ik ← hyperfine (nuclear spins)
+ Σk γn,k·(Bx′·Ix,k + By′·Iy,k + Bz′·Iz,k) ← nuclear Zeeman
+ Σk Pk·(Iz,k2 − Ik(Ik+1)/3·I) ← quadrupole (Ik ≥ 1)
The electron-only Hamiltonian (top five terms) is built by odmr_hamiltonian_Hz.
The full electron + nuclear Hamiltonian (all terms) is built by
full_hyperfine_hamiltonian_Hz, which assembles the tensor-product Hilbert space
H_e ⊗ H_n1 ⊗ H_n2 ⊗ ... automatically.
Applied B and E vectors are supplied in the lab frame and automatically rotated into the defect's local frame via the quantization axis before H is built.
ODMR observables come from transitions out of the |ms=0⟩-like eigenstate (the
one with highest overlap with the Sz eigenstate at index ms0_index).
NuclearSpin encapsulates a single nuclear spin with its hyperfine tensor,
gyromagnetic ratio, and (optionally) quadrupole coupling. All built-in defect
presets already carry their natural nuclear environments.
from SpinDefectSim.spin.nuclear import ( NuclearSpin, axial_A_tensor, isotropic_A_tensor, GAMMA_14N, GAMMA_15N, GAMMA_13C, GAMMA_11B, ) from SpinDefectSim.spin.hamiltonian import ( SpinParams, full_hyperfine_hamiltonian_Hz, odmr_transitions_Hz, ) # ── NV− on-site 14N (I = 1) ────────────────────────────────────────────── N14_NV = NuclearSpin( spin=1, A_tensor_Hz=axial_A_tensor(A_zz_Hz=-2.14e6, A_perp_Hz=-2.70e6), gamma_Hz_T=GAMMA_14N, label="14N", quadrupole_Hz=-4.95e6, ) # ── VB− three equivalent 14N neighbours (in-plane, axial A ≈ 0) ────────── N14_VB = NuclearSpin( spin=1, A_tensor_Hz=axial_A_tensor(A_zz_Hz=0, A_perp_Hz=47.8e6), gamma_Hz_T=GAMMA_14N, label="14N_VB", ) # ── Build the full electron + nuclear Hamiltonian ───────────────────────── sp = SpinParams(spin=1, D0_Hz=2.87e9) # NV− ground state H = full_hyperfine_hamiltonian_Hz(sp, E_vec_lab=[0, 0, 0], nuclear_spins=[N14_NV]) print(H.shape) # → (9, 9) [3 electron ×ばつ 3 nuclear] # Extract ODMR transitions (hyperfine-split lines) freqs = odmr_transitions_Hz(H, electron_dim=3, ms0_basis_index=1) print([f"{f/1e6:.2f} MHz" for f in freqs]) # ── 13C nearest-neighbour in diamond ───────────────────────────────────── C13 = NuclearSpin( spin=0.5, A_tensor_Hz=isotropic_A_tensor(13.0e6), # Fermi contact ~ 13 MHz gamma_Hz_T=GAMMA_13C, label="13C", ) # ── Use the preset NV− defect (includes 14N automatically) ─────────────── from SpinDefectSim.spin.hamiltonian import SpinDefect sd_nv = SpinDefect("nv_minus", B_mT=3.0) H_full = sd_nv.full_hamiltonian() # includes on-site 14N freqs_hf = sd_nv.hyperfine_transitions() # hyperfine-split ODMR lines
Isotope gyromagnetic ratios shipped with the library:
| Constant | Isotope | I | γn (MHz/T) |
|---|---|---|---|
GAMMA_14N |
14N | 1 | +3.077 |
GAMMA_15N |
15N | 1/2 | −4.316 |
GAMMA_11B |
11B | 3/2 | +13.660 |
GAMMA_10B |
10B | 3 | +4.575 |
GAMMA_13C |
13C | 1/2 | +10.708 |
GAMMA_29Si |
29Si | 1/2 | −5.319 |
RateModel solves the optical / ISC rate equations in steady state to predict
the CW ODMR contrast $C = (\text{PL}\text{off} - \text{PL}\text{on}) / \text{PL}_\text{off}$
for any spin-$S$ defect.
from SpinDefectSim.spin.rates import ( RateModel, RateParams, NV_RATES, VB_RATES, VSIC_RATES, P1_RATES, ) # ── Inspect pre-built contrasts ────────────────────────────────────────── for name, rp, ms0 in [("NV−", NV_RATES, 1), ("VB−", VB_RATES, 1), ("V_SiC",VSIC_RATES, 1), ("P1", P1_RATES, 0)]: model = RateModel(rp, ms0_index=ms0) print(f"{name:6s} C = {model.contrast()*100:.1f} %") # NV− C = 22.9 % # VB− C = 2.0 % # V_SiC C = 15.1 % # P1 C = 0.0 % ← no ISC path # ── Custom defect ──────────────────────────────────────────────────────── my_rates = RateParams( spin=1, k_optical=20e6, k_rad=100e6, k_isc_excited=[50e6, 5e6, 50e6], # ms = +1, 0, -1 k_from_shelving=[0.0, 10e6, 0.0], # spin-polarising return to ms = 0 ) model = RateModel(my_rates, ms0_index=1) print(f"contrast = {model.contrast()*100:.1f} %") # ── Photoluminescence vs. MW power ────────────────────────────────────── import numpy as np k_mw_range = np.logspace(3, 9, 200) # MW Rabi rate (Hz) pl_vals = [model.pl(k_mw) for k_mw in k_mw_range] # normalised PL
The rate model is integrated into Defaults.get_contrast(): if contrast=None
(the default), it auto-computes rate_params.
import sys; sys.path.insert(0, "/app") import numpy as np from SpinDefectSim.spin.hamiltonian import SpinDefect # VB− in hBN (default) — transition frequencies sd = SpinDefect(B_mT=1.5) f1, f2 = sd.transition_frequencies() # (Hz) # NV− in diamond, B along z-axis sd_nv = SpinDefect("nv_minus", B_vec_mT=[0, 0, 3.0]) # --- Run a sensing protocol on a single defect --- # E-field sensing: signal has E turned on, reference has E = 0 exp = sd_nv.to_experiment(E_vec_Vpm=[1e7, 0, 0], sensing="E") tau, *_, dS_peak = exp.echo_static() # B-field sensing: signal has extra stray B, reference has none exp = sd_nv.to_experiment(B_extra_T=[0, 0, 5e-5], sensing="B") f_axis = np.linspace(2.8e9, 2.95e9, 2000) pl_w, pl_no, dpl = exp.cw_odmr(f_axis) # Both E and B contribute to the contrast simultaneously exp = sd_nv.to_experiment( E_vec_Vpm=[1e7, 0, 0], B_extra_T=[0, 0, 5e-5], sensing="both", ) # Spin-1/2 P1 centre sd_p1 = SpinDefect("p1", B_mT=10.0) (f_single,) = sd_p1.transition_frequencies() # one transition # Fully custom spin-3/2 sd_custom = SpinDefect(spin=1.5, D0_Hz=2.0e9, B_mT=2.0) # Specific quantization axis (e.g. NV axis tilted 45° in xz-plane) sd_tilt = SpinDefect("nv_minus", B_mT=3.0, quantization_axis=[np.sin(np.pi/4), 0, np.cos(np.pi/4)])
from SpinDefectSim.base.params import Defaults # Use Defaults.for_defect() to load Hamiltonian params from a preset d = Defaults.for_defect("nv_minus", B_mT=5.0, T2echo=10e-6) sp = d.to_spin_params() # SpinParams with NV− D0, d⊥, spin=1, ms0_index=1 # Or keep VB− defaults and just change a field d = Defaults(B_mT=3.0, T2echo=20e-6)
import numpy as np from SpinDefectSim.spin.spectra import ensemble_transitions_from_Efields, ensemble_odmr_spectrum d = Defaults() sp = d.to_spin_params() f = np.linspace(3.35e9, 3.60e9, 2000) # Inhomogeneous ensemble — random E-field disorder E_fields = np.random.default_rng(0).normal(0, 5e4, (500, 3)) E_fields[:, 2] = 0 # in-plane disorder transitions = ensemble_transitions_from_Efields(E_fields, sp) fwhm = 1.0 / (np.pi * d.T2star) pl = ensemble_odmr_spectrum(f, transitions, fwhm, d.contrast)
from SpinDefectSim.spin.echo import lock_in_difference_echo, lock_in_difference_ramsey tr_signal = ensemble_transitions_from_Efields(E_fields, sp) tr_ref = ensemble_transitions_from_Efields(np.zeros_like(E_fields), sp) # Hahn-echo tau = np.linspace(0, 3 * d.T2echo, 600) S_w, S_n, dS = lock_in_difference_echo(tr_signal, tr_ref, tau, d.T2echo) # Ramsey FID tau_r = np.linspace(0, 3 * d.T2star, 400) S_w_r, S_n_r, dS_r, _ = lock_in_difference_ramsey(tr_signal, tr_ref, tau_r, d.T2star)
The sensing parameter of to_experiment() selects which fields contribute
to the contrast between signal and reference branches:
sensing |
Signal branch uses | Reference branch uses |
|---|---|---|
"E" |
E_fields + bias B | bias B only, E = 0 |
"B" |
stray B, E = 0 | bias B only, B_extra = 0 |
"both" |
E_fields + stray B | bias B only, both zero |
The bias B (from Defaults.B_mT) is always active in both branches.
from SpinDefectSim.analysis.ensemble import DefectEnsemble from scipy.constants import e as e_charge d = Defaults(T2echo=10e-6, B_mT=1.5) ens = DefectEnsemble(N_def=300, defaults=d) ens.generate_defects(seed=42) # ── E-field sensing ─────────────────────────────────────────────── rng = np.random.default_rng(7) xyzq = np.column_stack([ rng.uniform(-d.R_patch, d.R_patch, (20, 2)), np.zeros(20), rng.choice([-1, 1], 20) * e_charge, ]) ens.compute_efields(E0_gate=(0, 0, 1e4), disorder_xyzq=xyzq) exp_E = ens.to_experiment(sensing="E") tau, Sw, Sn, dS, tau_opt, dS_peak = exp_E.echo_static() print(f"E-sensing: Peak ΔS = {dS_peak:.4f} at τ = {tau_opt*1e6:.2f} μs") # ── B-field sensing (stray field from a magnetic sample) ────────── from SpinDefectSim.magnetometry import SquareGeometry geom = SquareGeometry(side=500e-9, n_boundary_pts=200) ens.compute_bfields(magnetization=lambda x, y: 1e-4, geometry=geom) exp_B = ens.to_experiment(sensing="B") print("B-sensing dS_peak:", exp_B.echo_static()[-1]) # ── Both fields contribute simultaneously ───────────────────────── exp_both = ens.to_experiment(sensing="both") print("Combined dS_peak:", exp_both.echo_static()[-1])
When defects have random orientations (e.g. in powdered samples or disordered films), generate uniformly distributed quantization axes:
ens = DefectEnsemble(N_def=500, defaults=d) ens.generate_defects(seed=0, quantization_axis="random") # isotropic powder ens.compute_efields(disorder_xyzq=xyzq) exp = ens.to_experiment() # Or a specific tilt shared by all defects: ens.generate_defects(seed=0, quantization_axis=[0, 0, 1]) # all along z ens.generate_defects(seed=0, quantization_axis=[1, 0, 0]) # all along x
Axes can also be set after placement:
ens.set_quantization_axis("random", seed=99)
from SpinDefectSim.sensing.sequences import RamseySequence, HahnEchoSequence, XY8Sequence from SpinDefectSim.sensing.snr import snr, n_avg_for_threshold tau_ac = 5e-6 # sensing τ for a 100 kHz AC field for name, seq in [("Ramsey", RamseySequence()), ("Hahn-echo", HahnEchoSequence()), ("XY8", XY8Sequence())]: r = seq.repetition_rate(tau_ac) print(f"{name}: {r/1e3:.1f} kHz rep rate") N_thresh = n_avg_for_threshold(0.05, snr_target=5.0, contrast=d.contrast, n_photons=d.n_photons) print(f"Averages for SNR=5: {N_thresh:.0f}")
Three screening models for the Coulomb interaction between disorder charges and defects:
screening_model |
Description |
|---|---|
None |
Bare Coulomb 1/r2 |
"yukawa" |
Yukawa: exp(−r/λ)/r2 with lambda_screen |
"dual_gate" |
Image-charge sum enforcing V=0 at ±d_gate |
from SpinDefectSim.electrometry.efield import E_disorder_point_charges E_vec = E_disorder_point_charges( obs_xyz=[0, 0, 0.34e-9], charges_xyzq=xyzq, epsilon_eff=7.0, screening_model="dual_gate", d_gate=15e-9, n_images=30, )
Any E-field source — COMSOL, Lumerical, meep, FEniCS, or a custom simulation — can be fed directly into the ensemble. Choose whichever path matches your data format.
The most common case: your solver exports field components on a Cartesian mesh.
# 2-D slice (e.g. COMSOL "Export > Data" at fixed z) # xs : (Nx,) array of x-coordinates in metres # ys : (Ny,) array of y-coordinates in metres # Ex2d, Ey2d, Ez2d : (Nx, Ny) arrays in V/m ens.efields_from_grid(Ex2d, Ey2d, Ez2d, xs, ys) # 3-D volumetric (e.g. FDTD output) — evaluated at z = z_defect ens.efields_from_grid(Ex3d, Ey3d, Ez3d, xs, ys, zs, z_defect=0.34e-9) # Nearest-neighbour instead of bilinear — useful for very coarse grids ens.efields_from_grid(Ex2d, Ey2d, Ez2d, xs, ys, method="nearest")
The add=True flag superimposes the grid field on top of an already-computed
disorder/gate contribution:
ens.compute_efields(disorder_xyzq=charges) # analytic disorder ens.efields_from_grid(Ex_fdtd, Ey_fdtd, Ez_fdtd, xs, ys, add=True) # + FDTD gate
Out-of-bounds defects are filled with zero by default (fill_value=0.0).
Set bounds_error=True to raise an error instead.
Wrap any interpolator, closed-form expression, or file-reading function:
# Uniform field — trivial example ens.efields_from_callable(lambda xyz: [1e4, 0, 0]) # Wrapping scipy RegularGridInterpolator built from solver output from scipy.interpolate import RegularGridInterpolator interp_x = RegularGridInterpolator((xs, ys), Ex2d) interp_y = RegularGridInterpolator((xs, ys), Ey2d) interp_z = RegularGridInterpolator((xs, ys), Ez2d) def my_field(xyz): pt = [[xyz[0], xyz[1]]] return [float(interp_x(pt)), float(interp_y(pt)), float(interp_z(pt))] ens.efields_from_callable(my_field) # add=True works here too ens.compute_efields(disorder_xyzq=charges) ens.efields_from_callable(my_field, add=True)
The callable receives xyz as a length-3 numpy array [x, y, z_defect] in metres
and must return [Ex, Ey, Ez] in V/m.
from SpinDefectSim.analysis.sweep import ParameterSweep import pandas as pd ps = ParameterSweep(N_def=200, seed=0) def run(n_charges): ens = ps.make_ensemble() rng2 = np.random.default_rng(1) n = int(n_charges) xyzq2 = np.column_stack([ rng2.uniform(-d.R_patch, d.R_patch, (n, 2)), np.zeros(n), rng2.choice([-1, 1], n) * e_charge, ]) ens.compute_efields(disorder_xyzq=xyzq2) _, _, _, _, tau_opt, dS_peak = ens.to_experiment().echo_static() return dict(n_charges=n, dS_peak=dS_peak, tau_opt_us=tau_opt * 1e6) results = ps.sweep(run, n_charges=[5, 10, 20, 40]) df = pd.DataFrame(results)
Two complementary workflows for E-field sensing:
Use ElectrometryExperiment to compute ODMR frequency maps across a sample
surface — the E-field analogue of MagnetometryExperiment.
import numpy as np from SpinDefectSim.electrometry import ElectrometryExperiment from SpinDefectSim.base.params import Defaults from scipy.constants import e as e_charge # Disorder charges at known positions charges = np.array([ [ 50e-9, 0.0, 0.0, e_charge], [-50e-9, 0.0, 0.0, -e_charge], ]) exp = ElectrometryExperiment(charges, Defaults(), z_defect=0.34e-9) # Point-by-point E-field E_vec = exp.E_field(0.0, 0.0) # (3,) V/m f1, f2 = exp.transition_frequencies(0.0, 0.0) # Hz # 2-D maps x_arr = np.linspace(-300e-9, 300e-9, 80) Ez_map = exp.E_z_map(x_arr, x_arr) # (80, 80) V/m df_map = exp.frequency_shift_map(x_arr, x_arr) # (80, 80) Hz # Gate-only experiment (no discrete charges) exp_gate = ElectrometryExperiment( charges_xyzq=None, defaults=Defaults(), z_defect=0.34e-9, E0_gate=(0, 0, 1e4), # uniform z-gate (V/m) epsilon_eff=4.0, screening_model="dual_gate", d_gate=15e-9, )
Use DefectEnsemble.compute_efields() to compute the screened-Coulomb E-field
at each defect position and feed it into the echo / Ramsey contrast:
ens.compute_efields( E0_gate=(0, 0, 1e4), # background gate bias disorder_xyzq=charges, # discrete disorder charges ) exp = ens.to_experiment(sensing="E") tau, *_, dS_peak = exp.echo_static() # Or import E fields directly from a FEM simulation: ens.efields_from_grid(Ex2d, Ey2d, Ez2d, xs, ys) exp = ens.to_experiment(sensing="E") # Or from a callable: ens.efields_from_callable(lambda xyz: my_comsol_interp(xyz[0], xyz[1]))
Two complementary workflows for B-field sensing:
Use MagnetometryExperiment to compute ODMR frequency maps across a sample
surface — the standard NV scanning-magnetometry workflow.
from SpinDefectSim.magnetometry import SquareGeometry, MagnetometryExperiment geom = SquareGeometry(side=1e-6, n_boundary_pts=300) mag_exp = MagnetometryExperiment( geom, magnetization=lambda x, y: 1e3 * np.exp(-(x**2 + y**2) / (250e-9)**2), defaults=Defaults(), z_defect=50e-9, bias_B_T=[0, 0, 1e-3], ) x_arr = np.linspace(-800e-9, 800e-9, 50) Bz_map = mag_exp.B_z_map(x_arr, x_arr) # (50, 50) array f1, f2 = mag_exp.transition_frequencies(0.0, 0.0)
Use DefectEnsemble.compute_bfields() to compute the Biot-Savart stray B at
each defect position and feed it into the echo / Ramsey contrast:
ens.compute_bfields( magnetization=lambda x, y: 1e-3, # uniform 1 mA geometry=geom, n_pts=80, ) exp = ens.to_experiment(sensing="B") tau, *_, dS_peak = exp.echo_static() # Or import B fields directly from a micromagnetic simulation: ens.bfields_from_grid(Bx2d, By2d, Bz2d, xs, ys) exp = ens.to_experiment(sensing="B") # Or from a callable (mumax3 wrap, OOMMF output, etc.): ens.bfields_from_callable(lambda xyz: my_mumax_interp(xyz[0], xyz[1]))
SensingExperiment is the central object for computing observables. It is
normally obtained from DefectEnsemble.to_experiment(), but can be constructed
directly.
| Method | Returns |
|---|---|
cw_odmr(f_axis) |
pl_with, pl_no, dpl |
ramsey() |
tau, S_w, S_n, dS, tau_opt, dS_peak |
echo_static() |
tau, S_w, S_n, dS, tau_opt, dS_peak |
echo_odmr_lockIn() |
f_axis, dpl |
snr(tau_s, N_avg) |
scalar SNR |
n_avg_to_detect() |
N averages for specified SNR at optimal τ |
Defaults( defect_type="vb_minus", # spin species D0_Hz=3.46e9, E0_Hz=50e6, d_perp=0.35, B_mT=1.5, T2star=50e-9, T2echo=10e-6, contrast=0.02, n_photons=500, N_def=1200, R_patch=200e-9, z_defect=0.34e-9, epsilon_eff=7.0, screening_model="dual_gate", d_gate=15e-9, ) Defaults.for_defect("nv_minus", B_mT=3.0) # classmethod
SpinDefect(defect_type="vb_minus", *, spin=None, B_mT=None, B_vec_mT=None, quantization_axis=None, D0_Hz=None, E0_Hz=None, d_perp=None, d_parallel=None, defaults=None) # Single-defect sensing — returns a full SensingExperiment exp = defect.to_experiment( E_vec_Vpm=(0, 0, 0), # lab-frame E-field for signal branch (V/m) B_extra_T=(0, 0, 0), # lab-frame stray B for signal branch (T) sensing="E" | "B" | "both", )
ens = DefectEnsemble(N_def, R_patch, defaults) ens.generate_defects(seed, quantization_axis=None|"random"|[x,y,z]|(N,3)) ens.generate_defects_gaussian(beam_waist_m, seed) ens.set_defects(positions) ens.set_quantization_axis(spec, seed) # E-field sources ens.compute_efields(E0_gate, gate_grad, disorder_xyzq, verbose) ens.efields_from_grid(Ex, Ey, Ez, x_coords, y_coords, [z_coords], *, add=False) ens.efields_from_callable(E_func, *, add=False) ens.set_efields(E_fields) # inject (N, 3) V/m directly # B-field sources (stray field signal) ens.compute_bfields(magnetization, geometry, *, n_pts, add=False) ens.bfields_from_grid(Bx, By, Bz, x_coords, y_coords, [z_coords], *, add=False) ens.bfields_from_callable(B_func, *, add=False) ens.set_bfields(B_fields) # inject (N, 3) T directly # Run experiment — sensing selects which fields create the contrast exp = ens.to_experiment(B_mT=None, sensing="E" | "B" | "both")
from SpinDefectSim.spin.defects import DefectType, get_defect, list_defects dt = get_defect("nv_minus") # retrieve preset list_defects() # print table of all presets dt = DefectType(name="custom", spin=1, D0_Hz=2e9, d_perp=0.1)
SpinDefectSim/
├── base/
│ ├── params.py Defaults dataclass, PhysicalParams base, mixins
│ └── mixins.py PlottingMixin, SerializationMixin, SweepMixin
├── spin/
│ ├── defects.py DefectType class and built-in presets
│ ├── matrices.py Spin operators for any spin S
│ ├── hamiltonian.py odmr_hamiltonian_Hz(), full_hyperfine_hamiltonian_Hz(), SpinDefect
│ ├── nuclear.py NuclearSpin dataclass, isotope γ constants, A-tensor helpers
│ ├── rates.py RateModel, RateParams, CW ODMR contrast from ISC equations
│ ├── spectra.py CW lineshapes, ensemble spectrum builders
│ └── echo.py Hahn-echo, Ramsey lock-in signals
├── coulomb/
│ └── kernels.py Screened Coulomb kernels (bare, Yukawa, dual-gate)
├── sensing/
│ ├── protocols.py SensingExperiment — CW, Ramsey, Hahn-echo
│ ├── sequences.py RamseySequence, HahnEchoSequence, XY8Sequence
│ └── snr.py Shot-noise SNR utilities
├── electrometry/ ← E-field sources and E-field imaging
│ ├── efield.py ElectricFieldBuilder, E_gate_bias, E_disorder_point_charges
│ └── electrometry.py ElectrometryExperiment — charges → E → ODMR maps
├── magnetometry/ ← B-field sources and B-field imaging
│ ├── geometry.py SquareGeometry, DiskGeometry sample shapes
│ ├── bfield.py Biot-Savart stray B-field calculations
│ └── magnetometry.py MagnetometryExperiment — M(x,y) → B → ODMR maps
└── analysis/
├── ensemble.py DefectEnsemble — positions, E-fields, B-fields, axes
└── sweep.py ParameterSweep wrapper