2. 加速度记录去趋势与基线校正的模态预处理

为什么海洋平台加速度记录必须先做去趋势与基线校正

导管架式海洋平台的加速度响应中，结构模态信号通常只占很小的幅值，而测量链路上叠加的低频成分往往更强：传感器直流偏置、积分器温漂、长周期波浪引起的准静态刚体运动、以及采集系统零点漂移，都会在时程上表现为缓慢的基线漂移或低频趋势。这些成分在频谱上集中在 0 Hz 附近，却会通过泄漏和窗函数旁瓣污染低阶模态频段（海洋平台一阶模态常在 0.1–0.5 Hz 量级），使峰值拾取偏移、半功率带宽估计失真。

本例取一段教学示例数据：采样频率 fs = 50 Hz，时长 600 s，共 30000 点。信号由三部分合成——目标模态分量（0.32 Hz 与 0.78 Hz 两个衰减正弦）、叠加的线性趋势与二次趋势（模拟基线漂移）、以及零均值高斯噪声。下面用同一套流程完成去趋势、基线校正与零均值化，并对比处理前后的频谱与模态参数估计。

处理流程与各环节的作用位置

关键点在于：去趋势与基线校正必须放在滤波和频谱估计之前，否则低频趋势会被后续带通滤波器的过渡带"搬移"到模态频段附近，形成虚假峰。同时，去趋势只能去除确定性趋势，不能替代高通滤波；对积分得到的位移/速度量，还需做基线校正以避免积分漂移累积。

可运行的预处理实现

下面用 Python 完成合成数据、去趋势、基线校正、零均值化，并对比处理前后的频谱峰值与半功率带宽阻尼比估计。所有参数均为教学示例值，结果属于由给定输入直接推出的预期现象，不是实测数据。

import numpy as np

fs = 50.0
T = 600.0
N = int(fs * T)
t = np.arange(N) / fs
rng = np.random.default_rng(2024)

# 目标模态分量：两个衰减正弦（模拟结构低阶模态响应）
f1, z1 = 0.32, 0.02
f2, z2 = 0.78, 0.03
modal = (1.0 * np.exp(-2*np.pi*z1*f1*t) * np.sin(2*np.pi*f1*t)
         + 0.6 * np.exp(-2*np.pi*z2*f2*t) * np.sin(2*np.pi*f2*t))

# 基线漂移：线性 + 二次趋势 + 慢变正弦（模拟长周期波浪准静态分量）
trend = 0.004 * t + 1.2e-6 * t**2 + 0.8 * np.sin(2*np.pi*0.008*t)

# 噪声
noise = 0.05 * rng.standard_normal(N)

a_raw = modal + trend + noise

def detrend_poly(x, order=2):
    """多项式最小二乘去趋势，返回残差"""
    tt = np.arange(len(x)) / fs
    coef = np.polyfit(tt, x, order)
    return x - np.polyval(coef, tt)

def baseline_correct(x, seg_sec=60.0):
    """分段均值基线校正：每段减去段内均值，抑制残余慢漂移"""
    seg = int(seg_sec * fs)
    nseg = len(x) // seg
    y = x.copy()
    for k in range(nseg):
        s = slice(k*seg, (k+1)*seg)
        y[s] = y[s] - np.mean(y[s])
    return y

def zero_mean(x):
    return x - np.mean(x)

# 处理链：去趋势 -> 基线校正 -> 零均值化
a_dt = detrend_poly(a_raw, order=2)
a_bc = baseline_correct(a_dt, seg_sec=60.0)
a_pre = zero_mean(a_bc)

def psd_peak(x, f_lo=0.1, f_hi=1.2):
    """Welch PSD，返回目标频段内峰值频率与半功率带宽阻尼比"""
    from scipy.signal import welch
    f, Pxx = welch(x, fs=fs, nperseg=4096, noverlap=2048)
    band = (f >= f_lo) & (f <= f_hi)
    fb, Pb = f[band], Pxx[band]
    i = np.argmax(Pb)
    fp = fb[i]
    half = Pb[i] / 2.0
    # 半功率点：峰值两侧首次低于半功率处
    left = np.where(Pb[:i] <= half)[0]
    right = np.where(Pb[i:] <= half)[0]
    if len(left) == 0 or len(right) == 0:
        return fp, np.nan
    fl = fb[left[-1]]
    fr = fb[i + right[0]]
    zeta = (fr - fl) / (2 * fp)
    return fp, zeta

for name, sig in [("原始", a_raw), ("去趋势", a_dt), ("+基线校正", a_bc), ("+零均值", a_pre)]:
    fp, zeta = psd_peak(sig)
    print(f"{name:10s} 峰值频率={fp:.3f} Hz  阻尼比估计={zeta:.4f}")

预期输出（示例数据下的典型结果，实际数值随随机种子略有浮动）：

如何判断预处理是否到位

• 时域检查：处理后的时程应围绕零轴波动，首尾段均值差应远小于模态分量幅值（本例中应小于 0.05）。
• 频域检查：0–0.1 Hz 频段的 PSD 应比处理前下降 1–2 个数量级，且目标模态峰（0.32 Hz、0.78 Hz）位置不移动。
• 参数一致性：去趋势阶数从 1 增到 2 时，峰值频率变化应小于频率分辨率（本例 4096 点、fs=50 Hz 时约 0.012 Hz）；若变化明显，说明趋势阶数选择不当或存在更强非线性漂移。
• 不要过度处理：基线校正分段过短（如小于目标模态周期的 5 倍）会削掉模态低频成分，导致阻尼比被高估。