Source code for temfield_mpylab.waveform_validation

"""Model-based validation of an in-loop amplitude-modulated field waveform."""

from dataclasses import dataclass
import math

import numpy as np
from scipy.optimize import curve_fit

from mpylab.tools.sin_fit import fit_sin

__all__ = [
    "AMWaveformLimits",
    "InLoopAMValidationResult",
    "analyse_am_waveform",
    "validate_am_fit",
    "validate_carrier_level",
    "validate_waveform_stability",
]


[docs] @dataclass(frozen=True) class AMWaveformLimits: """Engineering limits for TEMField's in-loop AM validation. The 2 dB difference between the ideal 5.1 dB AM peak ratio and the historical 3.1 dB saturation boundary motivates an outer compression limit. TEMField uses the more conservative 1 dB default because the EUT is already present during this validation. """ maximum_peak_compression_db: float = 1.0 minimum_modulation_depth_percent: float = 75.0 maximum_modulation_depth_percent: float = 85.0 minimum_modulation_frequency_khz: float = 0.95 maximum_modulation_frequency_khz: float = 1.05 minimum_r_squared: float = 0.995 maximum_normalized_rmse: float = 0.03 carrier_relative_tolerance: float = 0.02 maximum_peak_relative_tolerance: float = 0.02 maximum_stability_relative_spread: float = 0.02 required_stable_waveforms: int = 3 def __post_init__(self): positive = ( "maximum_peak_compression_db", "minimum_modulation_frequency_khz", "maximum_modulation_frequency_khz", "minimum_r_squared", "maximum_normalized_rmse", "carrier_relative_tolerance", "maximum_peak_relative_tolerance", "maximum_stability_relative_spread", ) for name in positive: value = float(getattr(self, name)) if not math.isfinite(value) or value < 0.0: raise ValueError("%s must be finite and non-negative" % name) if not ( 0.0 < self.minimum_modulation_depth_percent <= self.maximum_modulation_depth_percent < 100.0 ): raise ValueError("modulation-depth limits must lie between 0 and 100") if self.carrier_relative_tolerance >= 1.0: raise ValueError("carrier_relative_tolerance must be less than one") if self.minimum_modulation_frequency_khz > self.maximum_modulation_frequency_khz: raise ValueError("minimum modulation frequency exceeds maximum") if not 0.0 <= self.minimum_r_squared <= 1.0: raise ValueError("minimum_r_squared must lie between zero and one") if ( not isinstance(self.required_stable_waveforms, int) or self.required_stable_waveforms < 1 ): raise ValueError("required_stable_waveforms must be a positive integer")
[docs] def as_dict(self): """Return a serialization-friendly mapping.""" return dict(vars(self))
[docs] @dataclass(frozen=True) class InLoopAMValidationResult: """Result of TEMField's non-normative waveform validation. Parameters ---------- status : str Machine-readable decision status. passed : bool Whether the final waveform set satisfies all configured limits. reason : str Human-readable decision reason. factor : float Positive AM peak factor ``1 + m``. target_field_v_per_m : float Requested fitted carrier field in volts per metre. start_field_v_per_m : float Safe CW starting field in volts per metre. limits : AMWaveformLimits Engineering limits used for the decision. ramp_points : tuple Fitted waveform diagnostics acquired during the upward ramp. stable_waveforms : tuple Final waveform fits used for the reported decision. carrier_correction_attempts : int, optional Number of downward final-carrier corrections applied. carrier_correction_waveforms : tuple, optional Otherwise valid final fits superseded by the downward correction. method : str, optional Validation method identifier. normative : bool, optional Whether the method claims normative IEC status. ramp_strategy : str, optional Actor-level ramp strategy identifier. """ status: str passed: bool reason: str factor: float target_field_v_per_m: float start_field_v_per_m: float limits: AMWaveformLimits ramp_points: tuple stable_waveforms: tuple carrier_correction_attempts: int = 0 carrier_correction_waveforms: tuple = () method: str = "in_loop_waveform_rapp" normative: bool = False ramp_strategy: str = "adaptive_square_law"
[docs] def as_dict(self): """Return a serialization-friendly mapping.""" return { "status": self.status, "passed": self.passed, "reason": self.reason, "method": self.method, "normative": self.normative, "ramp_strategy": self.ramp_strategy, "factor": self.factor, "target_field_v_per_m": self.target_field_v_per_m, "start_field_v_per_m": self.start_field_v_per_m, "limits": self.limits.as_dict(), "ramp_points": list(self.ramp_points), "stable_waveforms": list(self.stable_waveforms), "carrier_correction_attempts": self.carrier_correction_attempts, "carrier_correction_waveforms": list( self.carrier_correction_waveforms ), "final_fit": ( None if not self.stable_waveforms else self.stable_waveforms[-1] ), }
def _rapp_envelope(normalized_input, scale, knee, smoothness): """Return the output envelope of a memoryless Rapp AM/AM model.""" drive = np.maximum(np.asarray(normalized_input, dtype=float), 0.0) exponent = 2.0 * smoothness return scale * drive / ( 1.0 + np.power(drive / knee, exponent) ) ** (1.0 / exponent)
[docs] def analyse_am_waveform(times_ms, values, modulation_depth_percent): """Fit sinusoidal and Rapp models to one measured field waveform. ``peak_compression_db`` is the differential AM/AM gain loss between the carrier input and the positive AM peak. It is independent of the fitted absolute field scale and is therefore directly comparable across frequencies and target field strengths. """ times = np.asarray(times_ms, dtype=float) field = np.asarray(values, dtype=float) if times.ndim != 1 or len(times) < 5: raise ValueError("waveform needs at least five one-dimensional samples") if field.shape != times.shape: raise ValueError("waveform time and field arrays must have equal lengths") if not np.all(np.isfinite(times)) or not np.all(np.isfinite(field)): raise ValueError("waveform contains non-finite values") if np.any(field < 0.0): raise ValueError("AM envelope contains negative field magnitudes") depth = float(modulation_depth_percent) / 100.0 if not math.isfinite(depth) or not 0.0 < depth < 1.0: raise ValueError("modulation depth must be between 0 and 100 percent") sine = fit_sin(times, field) amplitude = float(sine["amp"]) phase = float(sine["phase"]) if amplitude < 0.0: amplitude = -amplitude phase += math.pi fitted = amplitude * np.sin(float(sine["omega"]) * times + phase) + float( sine["offset"] ) residual = field - fitted residual_sum = float(np.sum(residual * residual)) total_sum = float(np.sum((field - np.mean(field)) ** 2)) offset = float(sine["offset"]) if offset <= 0.0: raise ValueError("fitted AM carrier field must be positive") normalized_input = 1.0 + depth * np.sin( float(sine["omega"]) * times + phase ) scale_guess = offset rapp_parameters, _covariance = curve_fit( _rapp_envelope, normalized_input, field, p0=(scale_guess, 3.0, 2.0), bounds=( (np.finfo(float).eps, 1.0, 0.5), (max(10.0 * float(np.max(field)), scale_guess * 10.0), 1000.0, 20.0), ), maxfev=20000, ) scale, knee, smoothness = (float(value) for value in rapp_parameters) rapp_fitted = _rapp_envelope(normalized_input, scale, knee, smoothness) rapp_residual = field - rapp_fitted rapp_residual_sum = float(np.sum(rapp_residual * rapp_residual)) gain_carrier = float(_rapp_envelope([1.0], 1.0, knee, smoothness)[0]) peak_input = 1.0 + depth gain_peak = float( _rapp_envelope([peak_input], 1.0, knee, smoothness)[0] / peak_input ) peak_compression_db = max( 0.0, 20.0 * math.log10(gain_carrier / gain_peak) ) normalized_rmse = math.sqrt(residual_sum / len(field)) / amplitude return { "offset": offset, "amplitude": amplitude, "modulation_depth_percent": amplitude / offset * 100.0, "frequency": abs(float(sine["freq"])), "minimum": float(np.min(field)), "maximum": float(np.max(field)), "fitted_peak": offset + amplitude, "rmse": math.sqrt(residual_sum / len(field)), "normalized_rmse": normalized_rmse, "r_squared": ( math.nan if total_sum == 0.0 else 1.0 - residual_sum / total_sum ), "rapp_peak_compression_db": peak_compression_db, "rapp_knee": knee, "rapp_smoothness": smoothness, "rapp_rmse": math.sqrt(rapp_residual_sum / len(field)), "sample_count": len(times), "sample_interval_ms": float(np.mean(np.diff(times))), "duration_ms": float(times[-1] - times[0]), }
[docs] def validate_am_fit( fit, target_field_v_per_m, modulation_depth_percent, limits, *, check_carrier=True, ): """Return ``(passed, reason)`` for one final-target waveform fit. ``check_carrier=False`` is used when several stable final waveforms are available and their mean carrier level is evaluated separately with :func:`validate_carrier_level`. All per-waveform signal-quality, modulation, compression, and peak-safety limits remain active. """ target = float(target_field_v_per_m) expected_peak = target * (1.0 + float(modulation_depth_percent) / 100.0) checks = [ ( limits.minimum_modulation_depth_percent <= fit["modulation_depth_percent"] <= limits.maximum_modulation_depth_percent, "measured modulation depth is outside the accepted interval", ), ( limits.minimum_modulation_frequency_khz <= fit["frequency"] <= limits.maximum_modulation_frequency_khz, "measured modulation frequency is outside the accepted interval", ), ( fit["r_squared"] >= limits.minimum_r_squared, "sinusoidal fit quality is below the accepted R-squared", ), ( fit["normalized_rmse"] <= limits.maximum_normalized_rmse, "normalized sinusoidal residual exceeds the accepted limit", ), ( fit["rapp_peak_compression_db"] <= limits.maximum_peak_compression_db, "modelled amplifier peak compression exceeds the accepted limit", ), ( fit["fitted_peak"] <= expected_peak * (1.0 + limits.maximum_peak_relative_tolerance), "fitted AM peak exceeds the safety limit", ), ] if check_carrier: lower = target * (1.0 - limits.carrier_relative_tolerance) upper = target * (1.0 + limits.carrier_relative_tolerance) checks[:0] = [ ( fit["offset"] >= lower, "fitted carrier field %.9g V/m is below %.9g V/m" % (fit["offset"], lower), ), ( fit["offset"] <= upper, "fitted carrier field %.9g V/m exceeds %.9g V/m" % (fit["offset"], upper), ), ] for passed, reason in checks: if not passed: return False, reason return True, "waveform satisfies the in-loop AM validation limits"
[docs] def validate_carrier_level(fits, target_field_v_per_m, limits): """Validate the mean carrier field of stable final-target waveforms.""" if not fits: return False, "no final-target carrier fields are available" target = float(target_field_v_per_m) mean_offset = float(np.mean([float(fit["offset"]) for fit in fits])) lower = target * (1.0 - limits.carrier_relative_tolerance) upper = target * (1.0 + limits.carrier_relative_tolerance) if mean_offset < lower: return ( False, "mean fitted carrier field %.9g V/m is below %.9g V/m" % (mean_offset, lower), ) if mean_offset > upper: return ( False, "mean fitted carrier field %.9g V/m exceeds %.9g V/m" % (mean_offset, upper), ) return ( True, "mean fitted carrier field %.9g V/m is within %.9g...%.9g V/m" % (mean_offset, lower, upper), )
[docs] def validate_waveform_stability(fits, limits): """Return whether consecutive final-target waveforms are sufficiently stable.""" if len(fits) < limits.required_stable_waveforms: return False, "not enough consecutive final-target waveforms" selected = fits[-limits.required_stable_waveforms:] # Every waveform already has to satisfy the compression limit on its own. # Near the linear limit, the Rapp knee and smoothness are weakly # identifiable, so requiring their derived near-zero compression estimates # to agree would primarily measure probe noise and optimizer variation. for key in ("offset", "modulation_depth_percent"): values = [float(fit[key]) for fit in selected] scale = abs(float(np.mean(values))) scale = max(scale, np.finfo(float).eps) relative_spread = (max(values) - min(values)) / scale if relative_spread > limits.maximum_stability_relative_spread: return False, "%s is not stable across consecutive waveforms" % key return True, "consecutive final-target waveforms are stable"