Coverage for src / tracekit / component / reactive.py: 93%
211 statements
« prev ^ index » next coverage.py v7.13.1, created at 2026-01-11 23:04 +0000
1"""Reactive component extraction for TraceKit.
3This module provides capacitance and inductance measurement from
4waveform data, including parasitic extraction.
7Example:
8 >>> from tracekit.component import measure_capacitance, measure_inductance
9 >>> C = measure_capacitance(voltage_trace, current_trace)
10 >>> L = measure_inductance(voltage_trace, current_trace)
12References:
13 IEEE 181-2011: Standard for Transitional Waveform Definitions
14"""
16from __future__ import annotations
18from dataclasses import dataclass, field
19from typing import TYPE_CHECKING, Any, Literal, cast
21import numpy as np
23from tracekit.core.exceptions import AnalysisError, InsufficientDataError
25if TYPE_CHECKING:
26 from numpy.typing import NDArray
28 from tracekit.core.types import WaveformTrace
31@dataclass
32class CapacitanceMeasurement:
33 """Result of capacitance measurement.
35 Attributes:
36 capacitance: Measured capacitance in Farads.
37 esr: Equivalent Series Resistance in ohms.
38 method: Measurement method used.
39 confidence: Confidence in measurement (0-1).
40 statistics: Additional measurement statistics.
41 """
43 capacitance: float
44 esr: float = 0.0
45 method: str = ""
46 confidence: float = 1.0
47 statistics: dict = field(default_factory=dict) # type: ignore[type-arg]
50@dataclass
51class InductanceMeasurement:
52 """Result of inductance measurement.
54 Attributes:
55 inductance: Measured inductance in Henrys.
56 dcr: DC Resistance in ohms.
57 q_factor: Quality factor at measurement frequency.
58 method: Measurement method used.
59 confidence: Confidence in measurement (0-1).
60 statistics: Additional measurement statistics.
61 """
63 inductance: float
64 dcr: float = 0.0
65 q_factor: float | None = None
66 method: str = ""
67 confidence: float = 1.0
68 statistics: dict = field(default_factory=dict) # type: ignore[type-arg]
71@dataclass
72class ParasiticExtraction:
73 """Result of parasitic parameter extraction.
75 Attributes:
76 capacitance: Parasitic capacitance in Farads.
77 inductance: Parasitic inductance in Henrys.
78 resistance: Parasitic resistance in ohms.
79 model_type: Equivalent circuit model type.
80 resonant_freq: Self-resonant frequency (if applicable).
81 fit_quality: Quality of model fit (R-squared).
82 """
84 capacitance: float
85 inductance: float
86 resistance: float
87 model_type: Literal["series_RLC", "parallel_RLC", "pi", "tee"]
88 resonant_freq: float | None = None
89 fit_quality: float = 0.0
92def measure_capacitance(
93 voltage_trace: WaveformTrace,
94 current_trace: WaveformTrace | None = None,
95 *,
96 method: Literal["charge", "slope", "frequency"] = "charge",
97 resistance: float | None = None,
98) -> CapacitanceMeasurement:
99 """Measure capacitance from voltage/current waveforms.
101 Calculates capacitance using the relationship C = Q/V or C = I/(dV/dt).
103 Args:
104 voltage_trace: Voltage waveform across capacitor.
105 current_trace: Current waveform through capacitor (optional for
106 some methods).
107 method: Measurement method:
108 - "charge": C = integral(I*dt) / delta_V
109 - "slope": C = I / (dV/dt)
110 - "frequency": Extract from RC time constant
111 resistance: Known resistance for frequency method.
113 Returns:
114 CapacitanceMeasurement with capacitance value.
116 Raises:
117 AnalysisError: If measurement conditions are not met.
118 InsufficientDataError: If insufficient samples are provided.
120 Example:
121 >>> C = measure_capacitance(voltage, current)
122 >>> print(f"C = {C.capacitance * 1e12:.1f} pF")
124 References:
125 COMP-002
126 """
127 voltage = voltage_trace.data.astype(np.float64)
128 sample_rate = voltage_trace.metadata.sample_rate
129 dt = 1.0 / sample_rate
131 if len(voltage) < 10:
132 raise InsufficientDataError(
133 "Capacitance measurement requires at least 10 samples",
134 required=10,
135 available=len(voltage),
136 analysis_type="capacitance",
137 )
139 if method == "charge" and current_trace is not None:
140 # C = Q / V = integral(I*dt) / delta_V
141 current = current_trace.data.astype(np.float64)
142 min_len = min(len(voltage), len(current))
143 voltage = voltage[:min_len]
144 current = current[:min_len]
146 # Integrate current to get charge
147 charge = np.cumsum(current) * dt
148 delta_v = np.max(voltage) - np.min(voltage)
150 if delta_v > 1e-10:
151 delta_q = np.max(charge) - np.min(charge)
152 capacitance = delta_q / delta_v
153 else:
154 raise AnalysisError("Voltage change too small for capacitance measurement")
156 # Estimate ESR from phase relationship
157 esr = _estimate_esr(voltage, current, sample_rate)
159 return CapacitanceMeasurement(
160 capacitance=float(abs(capacitance)),
161 esr=esr,
162 method="charge_integration",
163 confidence=0.9,
164 statistics={
165 "delta_v": delta_v,
166 "delta_q": delta_q,
167 "num_samples": min_len,
168 },
169 )
171 elif method == "slope" and current_trace is not None:
172 # C = I / (dV/dt)
173 current = current_trace.data.astype(np.float64)
174 min_len = min(len(voltage), len(current))
175 voltage = voltage[:min_len]
176 current = current[:min_len]
178 # Calculate dV/dt
179 dv_dt = np.diff(voltage) / dt
181 # Find region where dV/dt is significant
182 significant_mask = np.abs(dv_dt) > np.max(np.abs(dv_dt)) * 0.1
183 if np.sum(significant_mask) < 5: 183 ↛ 184line 183 didn't jump to line 184 because the condition on line 183 was never true
184 raise AnalysisError("Insufficient voltage slope for capacitance measurement")
186 # Use corresponding current values
187 current_for_slope = current[:-1][significant_mask]
188 dv_dt_significant = dv_dt[significant_mask]
190 # C = I / (dV/dt)
191 capacitance_values = current_for_slope / dv_dt_significant
192 capacitance = float(np.median(np.abs(capacitance_values)))
194 return CapacitanceMeasurement(
195 capacitance=capacitance,
196 method="slope",
197 confidence=0.85,
198 statistics={
199 "num_valid_points": int(np.sum(significant_mask)),
200 "capacitance_std": float(np.std(np.abs(capacitance_values))),
201 },
202 )
204 elif method == "frequency":
205 # Extract from RC time constant
206 if resistance is None:
207 raise AnalysisError("Resistance value required for frequency method")
209 # Find time constant from step response
210 tau = _extract_time_constant(voltage, sample_rate)
212 # C = tau / R
213 capacitance = tau / resistance
215 return CapacitanceMeasurement(
216 capacitance=float(capacitance),
217 method="time_constant",
218 confidence=0.8,
219 statistics={
220 "time_constant": tau,
221 "resistance": resistance,
222 },
223 )
225 else:
226 raise AnalysisError(f"Method '{method}' requires current_trace or resistance parameter")
229def measure_inductance(
230 voltage_trace: WaveformTrace,
231 current_trace: WaveformTrace | None = None,
232 *,
233 method: Literal["flux", "slope", "frequency"] = "slope",
234 resistance: float | None = None,
235) -> InductanceMeasurement:
236 """Measure inductance from voltage/current waveforms.
238 Calculates inductance using the relationship V = L * dI/dt.
240 Args:
241 voltage_trace: Voltage waveform across inductor.
242 current_trace: Current waveform through inductor (optional for
243 some methods).
244 method: Measurement method:
245 - "flux": L = integral(V*dt) / delta_I
246 - "slope": L = V / (dI/dt)
247 - "frequency": Extract from RL time constant
248 resistance: Known resistance for frequency method.
250 Returns:
251 InductanceMeasurement with inductance value.
253 Raises:
254 AnalysisError: If measurement conditions are not met.
255 InsufficientDataError: If insufficient samples are provided.
257 Example:
258 >>> L = measure_inductance(voltage, current)
259 >>> print(f"L = {L.inductance * 1e6:.1f} uH")
261 References:
262 COMP-003
263 """
264 voltage = voltage_trace.data.astype(np.float64)
265 sample_rate = voltage_trace.metadata.sample_rate
266 dt = 1.0 / sample_rate
268 if len(voltage) < 10:
269 raise InsufficientDataError(
270 "Inductance measurement requires at least 10 samples",
271 required=10,
272 available=len(voltage),
273 analysis_type="inductance",
274 )
276 if method == "flux" and current_trace is not None:
277 # L = flux / I = integral(V*dt) / delta_I
278 current = current_trace.data.astype(np.float64)
279 min_len = min(len(voltage), len(current))
280 voltage = voltage[:min_len]
281 current = current[:min_len]
283 # Integrate voltage to get flux linkage
284 flux = np.cumsum(voltage) * dt
285 delta_i = np.max(current) - np.min(current)
287 if delta_i > 1e-10:
288 delta_flux = np.max(flux) - np.min(flux)
289 inductance = delta_flux / delta_i
290 else:
291 raise AnalysisError("Current change too small for inductance measurement")
293 # Estimate DCR from steady-state
294 dcr = _estimate_dcr(voltage, current)
296 return InductanceMeasurement(
297 inductance=float(abs(inductance)),
298 dcr=dcr,
299 method="flux_integration",
300 confidence=0.9,
301 statistics={
302 "delta_i": delta_i,
303 "delta_flux": delta_flux,
304 "num_samples": min_len,
305 },
306 )
308 elif method == "slope" and current_trace is not None:
309 # L = V / (dI/dt)
310 current = current_trace.data.astype(np.float64)
311 min_len = min(len(voltage), len(current))
312 voltage = voltage[:min_len]
313 current = current[:min_len]
315 # Calculate dI/dt
316 di_dt = np.diff(current) / dt
318 # Find region where dI/dt is significant
319 significant_mask = np.abs(di_dt) > np.max(np.abs(di_dt)) * 0.1
320 if np.sum(significant_mask) < 5: 320 ↛ 321line 320 didn't jump to line 321 because the condition on line 320 was never true
321 raise AnalysisError("Insufficient current slope for inductance measurement")
323 # Use corresponding voltage values
324 voltage_for_slope = voltage[:-1][significant_mask]
325 di_dt_significant = di_dt[significant_mask]
327 # L = V / (dI/dt)
328 inductance_values = voltage_for_slope / di_dt_significant
329 inductance = float(np.median(np.abs(inductance_values)))
331 return InductanceMeasurement(
332 inductance=inductance,
333 method="slope",
334 confidence=0.85,
335 statistics={
336 "num_valid_points": int(np.sum(significant_mask)),
337 "inductance_std": float(np.std(np.abs(inductance_values))),
338 },
339 )
341 elif method == "frequency":
342 # Extract from RL time constant
343 if resistance is None: 343 ↛ 344line 343 didn't jump to line 344 because the condition on line 343 was never true
344 raise AnalysisError("Resistance value required for frequency method")
346 # Find time constant from step response
347 tau = _extract_time_constant(voltage, sample_rate)
349 # L = tau * R
350 inductance = tau * resistance
352 return InductanceMeasurement(
353 inductance=float(inductance),
354 method="time_constant",
355 confidence=0.8,
356 statistics={
357 "time_constant": tau,
358 "resistance": resistance,
359 },
360 )
362 else:
363 raise AnalysisError(f"Method '{method}' requires current_trace or resistance parameter")
366def extract_parasitics(
367 voltage_trace: WaveformTrace,
368 current_trace: WaveformTrace,
369 *,
370 model: Literal["series_RLC", "parallel_RLC"] = "series_RLC",
371 frequency_range: tuple[float, float] | None = None,
372) -> ParasiticExtraction:
373 """Extract parasitic R, L, C parameters from impedance measurement.
375 Fits an equivalent circuit model to measured voltage/current data
376 to extract parasitic component values.
378 Args:
379 voltage_trace: Voltage waveform.
380 current_trace: Current waveform.
381 model: Equivalent circuit model type.
382 frequency_range: Frequency range for analysis (Hz).
384 Returns:
385 ParasiticExtraction with R, L, C values.
387 Raises:
388 AnalysisError: If measurement conditions are not met.
389 InsufficientDataError: If insufficient samples are provided.
391 Example:
392 >>> params = extract_parasitics(voltage, current)
393 >>> print(f"C = {params.capacitance*1e12:.1f}pF, L = {params.inductance*1e9:.1f}nH")
395 References:
396 COMP-004
397 """
398 voltage = voltage_trace.data.astype(np.float64)
399 current = current_trace.data.astype(np.float64)
400 sample_rate = voltage_trace.metadata.sample_rate
402 min_len = min(len(voltage), len(current))
403 voltage = voltage[:min_len]
404 current = current[:min_len]
406 if min_len < 100:
407 raise InsufficientDataError(
408 "Parasitic extraction requires at least 100 samples",
409 required=100,
410 available=min_len,
411 analysis_type="parasitic_extraction",
412 )
414 # Compute impedance in frequency domain
415 from scipy.fft import fft, fftfreq
417 V_fft = fft(voltage)
418 I_fft = fft(current)
420 freqs = fftfreq(min_len, 1 / sample_rate)
422 # Use only positive frequencies
423 pos_mask = freqs > 0
424 freqs = freqs[pos_mask]
425 Z_fft = V_fft[pos_mask] / (I_fft[pos_mask] + 1e-20)
427 # Apply frequency range filter
428 if frequency_range is not None:
429 freq_mask = (freqs >= frequency_range[0]) & (freqs <= frequency_range[1])
430 freqs = freqs[freq_mask]
431 Z_fft = Z_fft[freq_mask]
433 if len(freqs) < 10: 433 ↛ 434line 433 didn't jump to line 434 because the condition on line 433 was never true
434 raise AnalysisError("Insufficient frequency points for parasitic extraction")
436 # Fit RLC model to impedance data
437 if model == "series_RLC":
438 R, L, C = _fit_series_rlc(freqs, Z_fft)
439 else: # parallel_RLC
440 R, L, C = _fit_parallel_rlc(freqs, Z_fft)
442 # Calculate resonant frequency
443 f_res = 1 / (2 * np.pi * np.sqrt(L * C)) if L > 0 and C > 0 else None
445 # Calculate fit quality
446 Z_model = _calculate_rlc_impedance(freqs, R, L, C, model)
447 ss_res = np.sum(np.abs(Z_fft - Z_model) ** 2)
448 ss_tot = np.sum(np.abs(Z_fft - np.mean(Z_fft)) ** 2)
449 r_squared = 1 - ss_res / (ss_tot + 1e-20)
451 return ParasiticExtraction(
452 capacitance=float(C),
453 inductance=float(L),
454 resistance=float(R),
455 model_type=model,
456 resonant_freq=f_res,
457 fit_quality=float(max(0, r_squared)),
458 )
461def _extract_time_constant(data: NDArray[np.float64], sample_rate: float) -> float:
462 """Extract time constant from step response."""
463 # Normalize data
464 data_min = np.min(data)
465 data_max = np.max(data)
466 if data_max - data_min < 1e-10: 466 ↛ 467line 466 didn't jump to line 467 because the condition on line 466 was never true
467 return 1 / sample_rate # Return minimum time constant
469 normalized = (data - data_min) / (data_max - data_min)
471 # Find 63.2% point (time constant)
472 target = 0.632
473 idx_raw = np.argmax(normalized >= target)
474 idx = int(idx_raw)
475 if idx == 0:
476 idx = int(len(data) // 2)
478 tau = float(idx) / sample_rate
479 min_tau = 1 / sample_rate
480 return float(tau if tau > min_tau else min_tau)
483def _estimate_esr(
484 voltage: NDArray[np.float64],
485 current: NDArray[np.float64],
486 sample_rate: float,
487) -> float:
488 """Estimate ESR from voltage/current phase relationship."""
489 # Use correlation to estimate resistive component
490 # ESR causes in-phase voltage drop
491 correlation = np.correlate(voltage - np.mean(voltage), current - np.mean(current))
492 if np.std(current) > 1e-10: 492 ↛ 495line 492 didn't jump to line 495 because the condition on line 492 was always true
493 esr = np.max(correlation) / (len(current) * np.var(current))
494 return float(abs(esr))
495 return 0.0
498def _estimate_dcr(voltage: NDArray[np.float64], current: NDArray[np.float64]) -> float:
499 """Estimate DC resistance from steady-state V/I."""
500 # Use last 10% of data as steady-state
501 n = len(voltage)
502 start = int(0.9 * n)
503 v_ss = np.mean(voltage[start:])
504 i_ss = np.mean(current[start:])
505 if abs(i_ss) > 1e-10: 505 ↛ 507line 505 didn't jump to line 507 because the condition on line 505 was always true
506 return float(v_ss / i_ss)
507 return 0.0
510def _fit_series_rlc(
511 freqs: NDArray[np.float64], Z: NDArray[np.complex128]
512) -> tuple[float, float, float]:
513 """Fit series RLC model to impedance data."""
514 omega = 2 * np.pi * freqs
516 # Initial estimates
517 R_init = float(np.real(np.mean(Z)))
518 # From imaginary part at low/high frequency
519 L_init = float(np.abs(np.imag(Z[-1])) / omega[-1]) if len(omega) > 0 else 1e-9
520 C_init = 1e-12
522 def model(omega: NDArray[np.float64], R: float, L: float, C: float) -> NDArray[np.complex128]:
523 return R + 1j * (omega * L - 1 / (omega * C + 1e-20))
525 try:
526 # Fit real and imaginary parts separately
527 np.real(Z)
528 np.imag(Z)
530 # Simple optimization
531 from scipy.optimize import minimize
533 def objective(params: NDArray[np.float64]) -> np.floating[Any]: # type: ignore[name-defined]
534 R, L, C = params
535 Z_model = model(omega, R, L, C)
536 return float(np.sum(np.abs(Z - Z_model) ** 2)) # type: ignore[return-value]
538 result = minimize(
539 objective,
540 [R_init, L_init, C_init],
541 bounds=[(1e-6, 1e6), (1e-15, 1e-3), (1e-15, 1e-3)],
542 method="L-BFGS-B",
543 )
544 return tuple(result.x)
545 except Exception:
546 return (R_init, L_init, C_init)
549def _fit_parallel_rlc(
550 freqs: NDArray[np.float64], Z: NDArray[np.complex128]
551) -> tuple[float, float, float]:
552 """Fit parallel RLC model to impedance data."""
553 # Convert to admittance
554 Y = 1 / (Z + 1e-20)
555 omega = 2 * np.pi * freqs
557 # Initial estimates
558 G_init = float(np.real(np.mean(Y)))
559 R_init = 1 / G_init if G_init > 1e-10 else 1e3
560 C_init = float(np.abs(np.imag(Y[-1])) / omega[-1]) if len(omega) > 0 else 1e-12
561 L_init = 1e-9
563 try:
564 from scipy.optimize import minimize
566 def objective(params: NDArray[np.float64]) -> np.floating[Any]: # type: ignore[name-defined]
567 R, L, C = params
568 Y_model = 1 / R + 1j * omega * C + 1 / (1j * omega * L + 1e-20)
569 Z_model = 1 / (Y_model + 1e-20)
570 return float(np.sum(np.abs(Z - Z_model) ** 2)) # type: ignore[return-value]
572 result = minimize(
573 objective,
574 [R_init, L_init, C_init],
575 bounds=[(1e-6, 1e9), (1e-15, 1e-3), (1e-15, 1e-3)],
576 method="L-BFGS-B",
577 )
578 return tuple(result.x)
579 except Exception:
580 return (R_init, L_init, C_init)
583def _calculate_rlc_impedance(
584 freqs: NDArray[np.float64],
585 R: float,
586 L: float,
587 C: float,
588 model: str,
589) -> NDArray[np.complex128]:
590 """Calculate impedance of RLC circuit."""
591 omega = 2 * np.pi * freqs
593 if model == "series_RLC":
594 Z: NDArray[np.complex128] = R + 1j * (omega * L - 1 / (omega * C + 1e-20))
595 return Z
596 else: # parallel_RLC
597 Y = 1 / R + 1j * omega * C + 1 / (1j * omega * L + 1e-20)
598 return cast("NDArray[np.complex128]", 1 / (Y + 1e-20))