Adriano
7249cb8a8c
Prima Cp/Cpk usava stdev campionaria e Pp/Ppk la population stdev: stesso indice due volte, proporzionale. Ora Cp/Cpk (capability a breve termine) usa sigma_within dalla moving range dei valori consecutivi in ordine cronologico; Pp/Ppk (performance complessiva) resta sulla stdev campionaria (n-1). MRbar == 0 -> Cp/Cpk None. Stessa sigma_within per UCL/LCL della carta di controllo a valori individuali (mean +/- 2.66*MRbar), non piu la stdev globale che allargava i limiti e nascondeva la deriva. Aggiunto campo opzionale std_dev_within a CapabilityData. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
290 lines
8.6 KiB
Python
290 lines
8.6 KiB
Python
"""SPC (Statistical Process Control) computation service.
|
|
|
|
Pure functions — no DB dependencies. Receives lists of floats and tolerance limits,
|
|
returns structured data matching schemas in schemas/statistics.py.
|
|
|
|
Uses only Python stdlib (math, statistics). No numpy/scipy needed.
|
|
"""
|
|
import math
|
|
import statistics as stats
|
|
from datetime import datetime
|
|
|
|
from src.backend.models.api.statistics import (
|
|
CapabilityData,
|
|
ControlChartData,
|
|
HistogramData,
|
|
SummaryData,
|
|
)
|
|
|
|
|
|
def compute_summary(
|
|
pass_fail_values: list[str],
|
|
) -> SummaryData:
|
|
"""Compute pass/fail/warning summary from a list of pass_fail strings.
|
|
|
|
Args:
|
|
pass_fail_values: List of "pass", "warning", or "fail" strings.
|
|
|
|
Returns:
|
|
SummaryData with counts and rates.
|
|
"""
|
|
total = len(pass_fail_values)
|
|
if total == 0:
|
|
return SummaryData(
|
|
total=0,
|
|
pass_count=0,
|
|
warning_count=0,
|
|
fail_count=0,
|
|
pass_rate=0.0,
|
|
warning_rate=0.0,
|
|
fail_rate=0.0,
|
|
)
|
|
|
|
pass_count = pass_fail_values.count("pass")
|
|
warning_count = pass_fail_values.count("warning")
|
|
fail_count = pass_fail_values.count("fail")
|
|
|
|
return SummaryData(
|
|
total=total,
|
|
pass_count=pass_count,
|
|
warning_count=warning_count,
|
|
fail_count=fail_count,
|
|
pass_rate=round(pass_count / total * 100, 2),
|
|
warning_rate=round(warning_count / total * 100, 2),
|
|
fail_rate=round(fail_count / total * 100, 2),
|
|
)
|
|
|
|
|
|
def _sigma_within(values: list[float]) -> float | None:
|
|
"""Estimate within-subgroup sigma from the moving range of individuals.
|
|
|
|
sigma_within = MRbar / d2, with d2 = 1.128 (moving range of consecutive
|
|
pairs, i.e. subgroup size 2). `values` must be in chronological order
|
|
(the order in which measurements arrived) — the moving range is not
|
|
meaningful on a re-sorted series.
|
|
|
|
Returns None when fewer than 2 values are available, or when MRbar == 0
|
|
(no variation between consecutive points: the index would be undefined).
|
|
"""
|
|
n = len(values)
|
|
if n < 2:
|
|
return None
|
|
moving_ranges = [abs(values[i] - values[i - 1]) for i in range(1, n)]
|
|
mr_bar = stats.mean(moving_ranges)
|
|
if mr_bar == 0:
|
|
return None
|
|
return mr_bar / 1.128
|
|
|
|
|
|
def compute_capability(
|
|
values: list[float],
|
|
utl: float | None,
|
|
ltl: float | None,
|
|
nominal: float | None,
|
|
) -> CapabilityData:
|
|
"""Compute capability indices Cp, Cpk, Pp, Ppk.
|
|
|
|
`values` must be in chronological order (arrival order), since Cp/Cpk are
|
|
derived from the moving range between consecutive points:
|
|
|
|
- sigma_within = MRbar / 1.128 (MRbar = mean of |x_i - x_(i-1)|)
|
|
- Cp = (UTL - LTL) / (6 * sigma_within)
|
|
- Cpk = min((UTL - mean) / (3 * sigma_within), (mean - LTL) / (3 * sigma_within))
|
|
- Pp/Ppk: same formulas using the overall sample std dev (n-1, `stats.stdev`)
|
|
of the whole series — this is the "process performance" pair and is
|
|
allowed to be wider than Cp/Cpk when there is drift between subgroups.
|
|
|
|
Args:
|
|
values: Measured values, in chronological order.
|
|
utl: Upper Tolerance Limit (None if not defined).
|
|
ltl: Lower Tolerance Limit (None if not defined).
|
|
nominal: Nominal value (None if not defined).
|
|
|
|
Returns:
|
|
CapabilityData with indices and statistics.
|
|
"""
|
|
n = len(values)
|
|
if n < 2:
|
|
return CapabilityData(
|
|
cp=None, cpk=None, pp=None, ppk=None,
|
|
mean=values[0] if n == 1 else 0.0,
|
|
std_dev=0.0, std_dev_within=None, n=n,
|
|
utl=utl, ltl=ltl, nominal=nominal,
|
|
)
|
|
|
|
mean = stats.mean(values)
|
|
# Overall sample std dev (n-1) for Pp/Ppk
|
|
std_dev_sample = stats.stdev(values)
|
|
# Within-subgroup sigma (moving range) for Cp/Cpk
|
|
sigma_within = _sigma_within(values)
|
|
|
|
cp = cpk = pp = ppk = None
|
|
|
|
if sigma_within is not None:
|
|
if utl is not None and ltl is not None:
|
|
cp = round((utl - ltl) / (6 * sigma_within), 4)
|
|
|
|
cpk_values = []
|
|
if utl is not None:
|
|
cpk_values.append((utl - mean) / (3 * sigma_within))
|
|
if ltl is not None:
|
|
cpk_values.append((mean - ltl) / (3 * sigma_within))
|
|
if cpk_values:
|
|
cpk = round(min(cpk_values), 4)
|
|
|
|
if std_dev_sample > 0:
|
|
if utl is not None and ltl is not None:
|
|
pp = round((utl - ltl) / (6 * std_dev_sample), 4)
|
|
|
|
ppk_values = []
|
|
if utl is not None:
|
|
ppk_values.append((utl - mean) / (3 * std_dev_sample))
|
|
if ltl is not None:
|
|
ppk_values.append((mean - ltl) / (3 * std_dev_sample))
|
|
if ppk_values:
|
|
ppk = round(min(ppk_values), 4)
|
|
|
|
return CapabilityData(
|
|
cp=cp, cpk=cpk, pp=pp, ppk=ppk,
|
|
mean=round(mean, 6),
|
|
std_dev=round(std_dev_sample, 6),
|
|
std_dev_within=round(sigma_within, 6) if sigma_within is not None else None,
|
|
n=n,
|
|
utl=utl, ltl=ltl, nominal=nominal,
|
|
)
|
|
|
|
|
|
def compute_control_chart(
|
|
values: list[float],
|
|
timestamps: list[datetime],
|
|
utl: float | None,
|
|
uwl: float | None,
|
|
lwl: float | None,
|
|
ltl: float | None,
|
|
nominal: float | None,
|
|
) -> ControlChartData:
|
|
"""Compute control chart data (individuals chart) with UCL/LCL and
|
|
out-of-control detection.
|
|
|
|
UCL = mean + 3*sigma_within = mean + 2.66*MRbar
|
|
LCL = mean - 3*sigma_within = mean - 2.66*MRbar
|
|
with sigma_within = MRbar/1.128 (see `_sigma_within`).
|
|
Out-of-control: points outside UCL/LCL.
|
|
|
|
Args:
|
|
values: Measured values in chronological order.
|
|
timestamps: Corresponding timestamps.
|
|
utl/uwl/lwl/ltl: Tolerance/warning limits.
|
|
nominal: Nominal value.
|
|
|
|
Returns:
|
|
ControlChartData with values, limits, and OOC indices.
|
|
"""
|
|
n = len(values)
|
|
if n == 0:
|
|
return ControlChartData(
|
|
values=[], timestamps=[], mean=0.0, ucl=0.0, lcl=0.0,
|
|
utl=utl, uwl=uwl, lwl=lwl, ltl=ltl, nominal=nominal,
|
|
out_of_control=[],
|
|
)
|
|
|
|
mean = stats.mean(values)
|
|
|
|
sigma_within = _sigma_within(values)
|
|
sigma = sigma_within if sigma_within is not None else 0.0
|
|
|
|
ucl = mean + 3 * sigma
|
|
lcl = mean - 3 * sigma
|
|
|
|
# Detect out-of-control points (outside UCL/LCL)
|
|
out_of_control = []
|
|
for i, v in enumerate(values):
|
|
if v > ucl or v < lcl:
|
|
out_of_control.append(i)
|
|
|
|
return ControlChartData(
|
|
values=values,
|
|
timestamps=timestamps,
|
|
mean=round(mean, 6),
|
|
ucl=round(ucl, 6),
|
|
lcl=round(lcl, 6),
|
|
utl=utl,
|
|
uwl=uwl,
|
|
lwl=lwl,
|
|
ltl=ltl,
|
|
nominal=nominal,
|
|
out_of_control=out_of_control,
|
|
)
|
|
|
|
|
|
def compute_histogram(
|
|
values: list[float],
|
|
n_bins: int = 20,
|
|
) -> HistogramData:
|
|
"""Compute histogram bin data and normal curve overlay.
|
|
|
|
Args:
|
|
values: Measured values.
|
|
n_bins: Number of histogram bins (default 20).
|
|
|
|
Returns:
|
|
HistogramData with bins, counts, and normal curve points.
|
|
"""
|
|
n = len(values)
|
|
if n == 0:
|
|
return HistogramData(
|
|
bins=[], counts=[], normal_x=[], normal_y=[],
|
|
mean=0.0, std_dev=0.0, n=0,
|
|
)
|
|
|
|
mean = stats.mean(values)
|
|
std_dev = stats.pstdev(values) if n >= 2 else 0.0
|
|
|
|
min_val = min(values)
|
|
max_val = max(values)
|
|
|
|
# Avoid zero-width range
|
|
if max_val == min_val:
|
|
max_val = min_val + 1.0
|
|
|
|
bin_width = (max_val - min_val) / n_bins
|
|
bins = [round(min_val + i * bin_width, 6) for i in range(n_bins + 1)]
|
|
|
|
# Count values per bin
|
|
counts = [0] * n_bins
|
|
for v in values:
|
|
idx = int((v - min_val) / bin_width)
|
|
if idx >= n_bins:
|
|
idx = n_bins - 1
|
|
counts[idx] += 1
|
|
|
|
# Normal curve overlay (100 points)
|
|
normal_x: list[float] = []
|
|
normal_y: list[float] = []
|
|
if std_dev > 0:
|
|
n_curve_points = 100
|
|
x_min = mean - 4 * std_dev
|
|
x_max = mean + 4 * std_dev
|
|
x_step = (x_max - x_min) / (n_curve_points - 1)
|
|
|
|
for i in range(n_curve_points):
|
|
x = x_min + i * x_step
|
|
# Normal PDF: (1 / (sigma * sqrt(2*pi))) * exp(-0.5 * ((x-mu)/sigma)^2)
|
|
y = (1.0 / (std_dev * math.sqrt(2 * math.pi))) * math.exp(
|
|
-0.5 * ((x - mean) / std_dev) ** 2
|
|
)
|
|
# Scale to match histogram: y * n * bin_width
|
|
y_scaled = y * n * bin_width
|
|
normal_x.append(round(x, 6))
|
|
normal_y.append(round(y_scaled, 4))
|
|
|
|
return HistogramData(
|
|
bins=bins,
|
|
counts=counts,
|
|
normal_x=normal_x,
|
|
normal_y=normal_y,
|
|
mean=round(mean, 6),
|
|
std_dev=round(std_dev, 6),
|
|
n=n,
|
|
)
|