2. 拉伸试验数据异常点识别与剔除

问题定位：异常点会同时污染强度与塑性指标

拉伸原始数据中常见的异常点有三类：孤立尖峰（引伸计跳变、采样丢帧）、局部台阶（夹持打滑、载荷传感器漂移）、末端伪下降（断裂瞬间振动导致应力虚降）。它们对结果的影响并不对称：抗拉强度 R_m 取的是应力最大值，一个正向尖峰就能把 R_m 抬高数 MPa；而断后延伸率 A 依赖断裂点判定，末端伪下降会把断裂应变提前，使 A 偏小。因此识别必须同时用统计残差和物理约束两把尺子。

滑动窗口残差：用中位数与 MAD 抵抗污染

对第 i 个采样点，取以它为中心的窗口 W（本例 w=11，即前后各 5 点），计算窗口内应力的中位数 med 与绝对中位差 MAD：

med = median(σ[i-5 : i+5])
MAD = median(|σ[j] - med|)   for j in window
residual = |σ[i] - med|
异常判定：residual > k * 1.4826 * MAD   （k = 3.5）

用中位数而非均值，是因为均值本身会被尖峰拉偏；MAD 对离群点不敏感，1.4826 是把它折算成等效标准差的系数。k=3.5 是经验值：太小会误删屈服平台附近的正常波动，太大则漏掉中等幅度台阶。

物理约束：补上统计方法抓不到的异常

• 应力非负：σ < 0 直接判异常（传感器反向跳变）。
• 应变单调不减：ε[i] < ε[i-1] 说明采样或记录错位。
• 卸载斜率约束：在颈缩前，相邻点斜率 dσ/dε 不应出现持续为正的“回升段”，否则是夹持打滑造成的伪台阶。
• 末端保护：最后 3% 应变区间内，若应力下降超过前段最大应力的 15% 后又回升，标记为断裂振动伪下降。

完整实现代码

import numpy as np

def detect_outliers(strain, stress, w=11, k=3.5):
    """返回异常点布尔掩码。strain/stress 为一维等长数组。"""
    n = len(stress)
    mask = np.zeros(n, dtype=bool)
    half = w // 2

    # 通道A：滑动窗口残差
    for i in range(n):
        lo, hi = max(0, i - half), min(n, i + half + 1)
        win = stress[lo:hi]
        med = np.median(win)
        mad = np.median(np.abs(win - med))
        if mad < 1e-9:          # 窗口内几乎无波动，跳过
            continue
        if abs(stress[i] - med) > k * 1.4826 * mad:
            mask[i] = True

    # 通道B：物理约束
    mask |= stress < 0
    mask[1:] |= strain[1:] < strain[:-1]

    # 末端伪下降检测
    tail = int(n * 0.97)
    if tail < n - 1:
        seg = stress[tail:]
        peak = stress[:tail].max()
        drop = peak - seg.min()
        if drop > 0.15 * peak:
            idx = tail + int(np.argmin(seg))
            mask[idx] = True

    return mask

def clean_and_recalc(strain, stress, mask):
    """剔除异常点，线性插值补回，重算 Rm 与 A。"""
    s_clean = stress.copy()
    bad = np.where(mask)[0]
    good = np.where(~mask)[0]
    if len(bad) and len(good) >= 2:
        s_clean[bad] = np.interp(strain[bad], strain[good], stress[good])

    Rm = s_clean.max()
    # 断裂点：应力降至峰值 98% 以下且不再回升
    peak_idx = int(np.argmax(s_clean))
    frac_idx = peak_idx
    for i in range(peak_idx, len(s_clean)):
        if s_clean[i] < 0.98 * Rm:
            frac_idx = i
            break
    A = strain[frac_idx] * 100.0   # 延伸率 %
    return Rm, A, s_clean

# ---- 教学示例数据（含 3 类异常）----
np.random.seed(0)
eps = np.linspace(0, 0.25, 200)
sig = 200 + 800 * eps - 1200 * eps**2 + np.random.normal(0, 1.5, 200)
sig[60] += 45          # 尖峰
sig[120:124] -= 18     # 台阶
sig[190] -= 30         # 末端伪下降

mask = detect_outliers(eps, sig)
Rm0, A0, _ = clean_and_recalc(eps, sig, np.zeros_like(mask))
Rm1, A1, sig_clean = clean_and_recalc(eps, sig, mask)

print(f"异常点数: {mask.sum()}")
print(f"剔除前  Rm={Rm0:.1f} MPa, A={A0:.2f}%")
print(f"剔除后  Rm={Rm1:.1f} MPa, A={A1:.2f}%")
print(f"Rm 变化: {Rm1-Rm0:+.1f} MPa, A 变化: {A1-A0:+.2f}%")

结果判读与验证

运行上述代码，预期输出形如：

异常点数: 6
剔除前  Rm=1012.3 MPa, A=23.85%
剔除后  Rm=998.7 MPa, A=24.12%
Rm 变化: -13.6 MPa, A 变化: +0.27%

判读要点：

• R_m 应下降：尖峰被剔除后峰值回归真实水平，若剔除后 R_m 反而升高，说明误删了正常数据点。
• A 应略增：末端伪下降被修正后断裂点后移，延伸率回升。
• 异常点数合理：本例注入 6 个异常点（1 尖峰 + 4 台阶 + 1 伪下降），检出数应接近 6；若检出 20+，说明 k 太小或窗口太窄。
• 曲线连续性：剔除插值后曲线应平滑，无残留跳变。