高位表中热接触和电线的随机扰动带来的影响（YUNSWJ建模）

Assumed '30 nV/√Hz (good)': en=30.0 nV/√Hz, fc=1.0 Hz, k_rw=0.15 µV/√s
  τ= 0.02 s: Allan σ_A ≈   0.112 µV
  τ= 0.20 s: Allan σ_A ≈   0.080 µV
  τ= 2.00 s: Allan σ_A ≈   0.214 µV
  τ=20.00 s: Allan σ_A ≈   0.671 µV

Assumptions (representative):
  en=30.0 nV/√Hz, 1/f corner fc=1.0 Hz, autozero timescale T_az=5.0s (f_L=0.200Hz)
  random-walk k_rw=0.15 µV/√s, floor=0.20 µV
  temperature OU: sigma_T=0.050 °C RMS, tau_T=10.0 s
  temp coefficients: alpha=0.20 ppm/°C, beta=1.00 µV/°C
  Range=100 mV => 1 count = 0.050 µV

Breakdown for Vin=80 mV:
  NPLC=  0.1: white= 0.335µV, 1/f= 0.080µV, RW+floor= 0.200µV, temp= 0.051µV => total= 0.402µV ( 8.04 counts)
  NPLC=  1.0: white= 0.106µV, 1/f= 0.066µV, RW+floor= 0.201µV, temp= 0.051µV => total= 0.242µV ( 4.84 counts)
  NPLC= 10.0: white= 0.034µV, 1/f= 0.048µV, RW+floor= 0.211µV, temp= 0.051µV => total= 0.225µV ( 4.49 counts)
  NPLC=100.0: white= 0.011µV, 1/f= 0.014µV, RW+floor= 0.292µV, temp= 0.049µV => total= 0.296µV ( 5.92 counts)

Part 1) β sweep numeric check (100 mV range, Vin=80 mV, 1 count=0.050 µV):

β=0.1 µV/°C
  NPLC=  0.1: total= 0.399 µV ( 7.97 counts) | temp= 0.006 µV
  NPLC=  1.0: total= 0.237 µV ( 4.74 counts) | temp= 0.006 µV
  NPLC= 10.0: total= 0.219 µV ( 4.38 counts) | temp= 0.006 µV
  NPLC=100.0: total= 0.292 µV ( 5.84 counts) | temp= 0.006 µV

β=1 µV/°C
  NPLC=  0.1: total= 0.402 µV ( 8.04 counts) | temp= 0.051 µV
  NPLC=  1.0: total= 0.242 µV ( 4.84 counts) | temp= 0.051 µV
  NPLC= 10.0: total= 0.225 µV ( 4.49 counts) | temp= 0.051 µV
  NPLC=100.0: total= 0.296 µV ( 5.92 counts) | temp= 0.049 µV

β=10 µV/°C
  NPLC=  0.1: total= 0.640 µV (12.80 counts) | temp= 0.501 µV
  NPLC=  1.0: total= 0.554 µV (11.08 counts) | temp= 0.501 µV
  NPLC= 10.0: total= 0.545 µV (10.90 counts) | temp= 0.499 µV
  NPLC=100.0: total= 0.566 µV (11.32 counts) | temp= 0.485 µV

import numpy as np
import matplotlib.pyplot as plt

# -----------------------------
# Extend: show 10 V / 1 V / 100 mV ranges and Allan deviation view
# -----------------------------

f_line = 50.0
nplc = np.logspace(-1, 2, 600)
Ti = nplc / f_line

# Ranges to compare
ranges = [10.0, 1.0, 0.1]  # V_FS
counts_6p5 = 2_000_000

# Noise model params (same as previous)
en_list = [
    ("10 nV/√Hz (very good)", 10e-9),
    ("30 nV/√Hz (good)", 30e-9),
    ("100 nV/√Hz (okay)", 100e-9),
]
T_az = 5.0
f_L = 1.0 / T_az
f_c_map = {
    "10 nV/√Hz (very good)": 0.2,
    "30 nV/√Hz (good)": 1.0,
    "100 nV/√Hz (okay)": 5.0,
}
k_rw_map = {
    "10 nV/√Hz (very good)": 0.05e-6,
    "30 nV/√Hz (good)": 0.15e-6,
    "100 nV/√Hz (okay)": 0.40e-6,
}
sigma_floor = 0.2e-6

def sigma_white(en, Ti):
    return en / (2.0 * np.sqrt(Ti))

def sigma_1f(en, f_c, Ti):
    f_H = 1.0 / (2.0 * Ti)
    ratio = np.maximum(f_H / f_L, 1.0)
    var = (en**2 * f_c) * np.log(ratio)
    return np.sqrt(var)

def sigma_drift(k_rw, Ti):
    return np.sqrt((k_rw**2) * Ti + sigma_floor**2)

def sigma_total(label, en, Ti):
    f_c = f_c_map[label]
    k_rw = k_rw_map[label]
    s_w = sigma_white(en, Ti)
    s_1f = sigma_1f(en, f_c, Ti)
    s_d = sigma_drift(k_rw, Ti)
    return np.sqrt(s_w**2 + s_1f**2 + s_d**2), s_w, s_1f, s_d

# Choose a representative quality for multi-range plot to avoid clutter
label_sel = "30 nV/√Hz (good)"
en_sel = dict(en_list)[label_sel]
sig_tot, sig_w, sig_1f, sig_d = sigma_total(label_sel, en_sel, Ti)

# --- Plot A: counts RMS vs NPLC for different ranges (6½-digit counts fixed) ---
plt.figure(figsize=(8.6, 5.4))
for V_FS in ranges:
    V_LSB = V_FS / counts_6p5
    plt.loglog(nplc, sig_tot / V_LSB, label=f"Range {V_FS:g} V (1 count={V_LSB*1e6:.3f} µV)")

plt.axhline(1.0, linestyle="--", linewidth=1, label="1 count RMS")
plt.axhline(0.3, linestyle=":", linewidth=1, label="0.3 count RMS")
plt.xlabel("NPLC")
plt.ylabel("RMS noise σ (counts)")
plt.title(f"Counts noise vs NPLC for different ranges (assumed '{label_sel}', 50 Hz)")
plt.grid(True, which="both")
plt.legend()
plt.tight_layout()
plt.show()

# --- Plot B: absolute noise σV vs NPLC with component breakdown (selected quality) ---
plt.figure(figsize=(8.6, 5.4))
plt.loglog(nplc, sig_tot, label="Total")
plt.loglog(nplc, sig_w, linestyle="--", label="White")
plt.loglog(nplc, sig_1f, linestyle="--", label="1/f")
plt.loglog(nplc, sig_d, linestyle="--", label="Drift/thermal")
plt.xlabel("NPLC")
plt.ylabel("RMS noise σ_V (V)")
plt.title(f"Noise components vs NPLC (assumed '{label_sel}', 50 Hz)")
plt.grid(True, which="both")
plt.legend()
plt.tight_layout()
plt.show()

# --- Allan deviation (heuristic mapping) ---
# For a measurement averaged over time tau, Allan deviation contributions:
# - White noise ~ A / sqrt(tau)
# - Flicker (1/f) ~ B (flat)
# - Random walk ~ C * sqrt(tau)
# We'll fit A,B,C from our assumed params for a rough "DMM-like" curve.
tau = np.logspace(-2, 2, 600)  # 0.01 .. 100 s averaging time

# Map A from en: for averaging, sigma_white(tau) ~ en/sqrt(2*tau) (order) ; keep consistent scale with our model
A = en_sel / 2.0  # so white ~ A / sqrt(tau) matches earlier sigma_white with Ti=tau
# Flicker floor: approximate from 1/f at tau=1 s (order)
B = sigma_1f(en_sel, f_c_map[label_sel], np.array([1.0]))[0]
# Random-walk coefficient from k_rw
C = k_rw_map[label_sel]

allan = np.sqrt((A**2)/tau + B**2 + (C**2)*tau)

plt.figure(figsize=(8.6, 5.4))
plt.loglog(tau, allan, label="Allan deviation (heuristic)")
plt.loglog(tau, A/np.sqrt(tau), linestyle="--", label="White ~ 1/√τ")
plt.loglog(tau, np.full_like(tau, B), linestyle="--", label="Flicker ~ constant")
plt.loglog(tau, C*np.sqrt(tau), linestyle="--", label="Random walk ~ √τ")
plt.xlabel("Averaging time τ (s)")
plt.ylabel("Allan deviation σ_A(τ) (V)")
plt.title(f"Heuristic Allan deviation (assumed '{label_sel}')")
plt.grid(True, which="both")
plt.legend()
plt.tight_layout()
plt.show()

# Print a quick table: at τ = 0.02, 0.2, 2, 20 s
tau_pts = np.array([0.02, 0.2, 2.0, 20.0])
allan_pts = np.sqrt((A**2)/tau_pts + B**2 + (C**2)*tau_pts)
print(f"Assumed '{label_sel}': en={en_sel*1e9:.1f} nV/√Hz, fc={f_c_map[label_sel]} Hz, k_rw={C*1e6:.2f} µV/√s")
for t, a in zip(tau_pts, allan_pts):
    print(f"  τ={t:>5.2f} s: Allan σ_A ≈ {a*1e6:>7.3f} µV")

import numpy as np
import matplotlib.pyplot as plt

# -----------------------------
# Extend drift model:
#   drift/thermal = random walk + linear tempco × temperature random process (OU)
# This captures: "touch the leads" / "air draft" -> slow wander.
# -----------------------------

# Base measurement setup
f_line = 50.0
nplc = np.logspace(-1, 2, 600)              # 0.1..100
Ti = nplc / f_line                           # seconds (integration/run-up time proxy)
counts_6p5 = 2_000_000
ranges = [10.0, 1.0, 0.1]                    # 10V / 1V / 100mV

# White + 1/f front-end assumptions (same as before; representative "good" DMM front-end)
label_sel = "30 nV/√Hz (good)"
en_sel = 30e-9                               # V/sqrt(Hz)
f_c = 1.0                                    # Hz (1/f corner)
T_az = 5.0                                   # s (autozero / reversal effective reset time)
f_L = 1.0 / T_az

# Random-walk drift (offset wander independent of input, e.g., thermals / contact)
k_rw = 0.15e-6                               # V/sqrt(s)
sigma_floor = 0.2e-6                         # V RMS (micro-thermal floor)

# Temperature random process (ambient + local gradients):
# Model ΔT(t) as OU process with std sigma_T and correlation time tau_T.
sigma_T = 0.05                               # °C RMS (gentle lab air movement level)
tau_T = 10.0                                 # s correlation time (draft / convection)

# Linear temperature coefficients
# Gain tempco (ppm/°C) multiplies Vin; offset thermo-EMF coefficient (µV/°C) adds regardless of Vin.
alpha_ppm_per_C = 0.2                        # ppm/°C (good precision front-end + resistors)
beta_uV_per_C = 1.0                          # µV/°C (thermoelectric at junctions/leads; can be 0.1~10)

alpha = alpha_ppm_per_C * 1e-6               # 1/°C
beta = beta_uV_per_C * 1e-6                  # V/°C

def sigma_white(en, Ti):
    # dual-slope-like: sigma ≈ en/(2*sqrt(Ti))
    return en / (2.0 * np.sqrt(Ti))

def sigma_1f(en, f_c, Ti):
    # 1/f PSD: en^2 * f_c / f ; integrate with averaging-like effective passband [f_L, f_H]
    f_H = 1.0 / (2.0 * Ti)
    ratio = np.maximum(f_H / f_L, 1.0)
    var = (en**2 * f_c) * np.log(ratio)
    return np.sqrt(var)

def sigma_rw(k_rw, Ti):
    # random-walk offset: grows ~ sqrt(T)
    return k_rw * np.sqrt(Ti)

def sigma_drift_base(k_rw, Ti):
    # random-walk + floor
    return np.sqrt((k_rw**2) * Ti + sigma_floor**2)

def var_mean_OU(sigma_T, tau, T):
    """
    OU process with stationary variance sigma_T^2 and correlation time tau:
    autocov R(Δ)=sigma_T^2 * exp(-|Δ|/tau)
    Variance of time-average over window length T:
      Var(mean) = (2/T^2) * ∫_0^T (T-Δ) R(Δ) dΔ
                = sigma_T^2 * (2*tau/T^2) * (T - tau*(1 - exp(-T/tau)))
    """
    T = np.asarray(T)
    return (sigma_T**2) * (2.0 * tau / (T**2)) * (T - tau * (1.0 - np.exp(-T / tau)))

def sigma_temp_linear(Vin, Ti):
    # Temperature-induced voltage error = (beta + alpha*Vin) * mean(ΔT over Ti)
    # RMS = |beta + alpha*Vin| * sqrt(Var(mean))
    gain = np.abs(beta + alpha * Vin)  # V/°C
    return gain * np.sqrt(var_mean_OU(sigma_T, tau_T, Ti))

def sigma_total(Vin, Ti):
    s_w = sigma_white(en_sel, Ti)
    s_1f = sigma_1f(en_sel, f_c, Ti)
    s_d0 = sigma_drift_base(k_rw, Ti)
    s_T = sigma_temp_linear(Vin, Ti)
    s_tot = np.sqrt(s_w**2 + s_1f**2 + s_d0**2 + s_T**2)
    return s_tot, s_w, s_1f, s_d0, s_T

# ---- Plot 1: counts RMS vs NPLC for different ranges (including temp-induced drift) ----
plt.figure(figsize=(8.8, 5.6))
for V_FS in ranges:
    Vin = 0.8 * V_FS  # assume measuring near 80% of range (worst-ish for gain tempco term)
    V_LSB = V_FS / counts_6p5
    s_tot, *_ = sigma_total(Vin, Ti)
    plt.loglog(nplc, s_tot / V_LSB, label=f"Range {V_FS:g} V (Vin≈{Vin:g}V, 1 count={V_LSB*1e6:.3f} µV)")

plt.axhline(1.0, linestyle="--", linewidth=1, label="1 count RMS")
plt.axhline(0.3, linestyle=":", linewidth=1, label="0.3 count RMS")
plt.xlabel("NPLC")
plt.ylabel("RMS noise σ (counts)")
plt.title("6½-digit counts noise vs NPLC (white + 1/f + RW drift + temperature process)")
plt.grid(True, which="both")
plt.legend(fontsize=8)
plt.tight_layout()
plt.show()

# ---- Plot 2: component breakdown for a representative range (100 mV) ----
V_FS_rep = 0.1
Vin_rep = 0.08
s_tot, s_w, s_1f, s_d0, s_T = sigma_total(Vin_rep, Ti)

plt.figure(figsize=(8.8, 5.6))
plt.loglog(nplc, s_tot, label="Total")
plt.loglog(nplc, s_w, linestyle="--", label="White")
plt.loglog(nplc, s_1f, linestyle="--", label="1/f")
plt.loglog(nplc, s_d0, linestyle="--", label="RW drift + floor")
plt.loglog(nplc, s_T, linestyle="--", label="Tempco × OU(ΔT)")
plt.xlabel("NPLC")
plt.ylabel("RMS noise σ_V (V)")
plt.title("Noise components vs NPLC (example: 100 mV range, Vin≈80 mV)")
plt.grid(True, which="both")
plt.legend()
plt.tight_layout()
plt.show()

# ---- Plot 3: time trace showing "touch the leads" (temperature step + decay) ----
# Simulate temperature: OU noise + touch event modeled as exponential bump ΔT_touch * exp(-(t-t0)/tau_touch)
fs = 10.0          # Hz sampling for visualization
dt = 1/fs
t = np.arange(0, 120.0, dt)
tau_touch = 25.0
dT_touch = 1.0     # °C (finger contact / warm air locally can easily be ~0.5-2°C transient)

t0 = 30.0
touch = np.where(t >= t0, dT_touch * np.exp(-(t - t0)/tau_touch), 0.0)

# OU noise generation (discrete-time exact update)
np.random.seed(0)
a = np.exp(-dt/tau_T)
sigma_eps = sigma_T * np.sqrt(1 - a*a)
T_ou = np.zeros_like(t)
for i in range(1, len(t)):
    T_ou[i] = a*T_ou[i-1] + sigma_eps*np.random.randn()

dT = T_ou + touch

# Convert temperature to voltage error (offset + gain*Vin)
Vin_trace = 0.08  # 80 mV signal
v_temp = (beta + alpha*Vin_trace) * dT  # V

# Add random-walk offset drift (discrete)
k_step = k_rw * np.sqrt(dt)   # V per step RMS for random walk
rw = np.cumsum(k_step * np.random.randn(len(t)))
v_rw = rw

v_total = v_temp + v_rw  # show slow wander components only (not including white noise)

plt.figure(figsize=(9.2, 5.6))
plt.plot(t, dT)
plt.xlabel("Time (s)")
plt.ylabel("ΔT (°C)")
plt.title("Temperature process: OU noise + touch event (ΔT bump with decay)")
plt.grid(True)
plt.tight_layout()
plt.show()

plt.figure(figsize=(9.2, 5.6))
plt.plot(t, v_total*1e6)
plt.xlabel("Time (s)")
plt.ylabel("Equivalent error (µV)")
plt.title("Equivalent slow error from temperature + random-walk drift (example: 100 mV range, Vin=80 mV)")
plt.grid(True)
plt.tight_layout()
plt.show()

# ---- Quick numeric summary at a few NPLC points for 100 mV range ----
points = np.array([0.1, 1, 10, 100.0])
Ti_p = points / f_line
Vin_p = 0.08
V_LSB_p = V_FS_rep / counts_6p5

s_tot_p, s_w_p, s_1f_p, s_d0_p, s_T_p = sigma_total(Vin_p, Ti_p)

print("Assumptions (representative):")
print(f"  en={en_sel*1e9:.1f} nV/√Hz, 1/f corner fc={f_c} Hz, autozero timescale T_az={T_az}s (f_L={f_L:.3f}Hz)")
print(f"  random-walk k_rw={k_rw*1e6:.2f} µV/√s, floor={sigma_floor*1e6:.2f} µV")
print(f"  temperature OU: sigma_T={sigma_T:.3f} °C RMS, tau_T={tau_T:.1f} s")
print(f"  temp coefficients: alpha={alpha_ppm_per_C:.2f} ppm/°C, beta={beta_uV_per_C:.2f} µV/°C")
print(f"  Range=100 mV => 1 count = {V_LSB_p*1e6:.3f} µV")
print("\nBreakdown for Vin=80 mV:")
for n, w, onef, d0, st, tot in zip(points, s_w_p, s_1f_p, s_d0_p, s_T_p, s_tot_p):
    print(f"  NPLC={n:>5}: white={w*1e6:>6.3f}µV, 1/f={onef*1e6:>6.3f}µV, RW+floor={d0*1e6:>6.3f}µV, temp={st*1e6:>6.3f}µV"
          f" => total={tot*1e6:>6.3f}µV ({tot/V_LSB_p:>5.2f} counts)")