4. Apply CFC Filtering to Injury Curves

Why CFC Filtering Changes the Numbers You Report

SAE J211 defines channel frequency classes (CFC) as a way to make injury-relevant time histories comparable across labs, sleds, and simulation solvers. A CFC class is not a "smoothing preference" — it is a specified low-pass characteristic with a defined magnitude response and phase behavior. When you apply CFC 180 to a chest acceleration channel and CFC 600 to a head acceleration channel, you are not just cleaning up noise; you are changing the peak amplitude, the time of the peak, and therefore every downstream injury metric that depends on those two quantities.

The practical consequence is that an unfiltered peak and a CFC-filtered peak are different physical claims. Reporting a HIC15 or a 3 ms clip value without stating the CFC class is an incomplete result. This case builds a small, self-contained filtering routine, applies CFC 180 and CFC 600 to representative acceleration and force histories, and quantifies the peak shift so the effect on injury metrics is documented rather than assumed.

Filter Definition and Discretization

SAE J211 specifies a second-order (two-pole) Butterworth low-pass characteristic per CFC class. The class number is the cutoff frequency in hertz at which the magnitude response is −3 dB. For CFC 180, f_c = 180 Hz; for CFC 600, f_c = 600 Hz. The analog prototype is:

H(s) = ωc² / (s² + √2·ωc·s + ωc²),   ωc = 2π·fc

For sampled data the standard approach is a bilinear-transform (or matched) second-order IIR section. A robust, dependency-free implementation is a direct-form biquad with coefficients derived from the bilinear transform and pre-warped cutoff:

import numpy as np

def cfc_coeffs(fc, fs):
    """Second-order Butterworth low-pass, bilinear transform, pre-warped."""
    wc = 2.0 * np.pi * fc
    k  = wc / (2.0 * fs)                 # pre-warped analog-to-digital factor
    # Bilinear transform of the 2-pole Butterworth prototype
    K  = np.tan(np.pi * fc / fs)
    norm = 1.0 / (1.0 + np.sqrt(2.0) * K + K * K)
    b0 = K * K * norm
    b1 = 2.0 * b0
    b2 = b0
    a1 = 2.0 * (K * K - 1.0) * norm
    a2 = (1.0 - np.sqrt(2.0) * K + K * K) * norm
    return np.array([b0, b1, b2]), np.array([1.0, a1, a2])

def cfc_filter(x, fc, fs):
    """Zero-phase forward-backward biquad (filtfilt-equivalent, no scipy)."""
    b, a = cfc_coeffs(fc, fs)
    y = _biquad(x, b, a)          # forward pass
    y = _biquad(y[::-1], b, a)[::-1]  # backward pass
    return y

def _biquad(x, b, a):
    y = np.zeros_like(x, dtype=float)
    x1 = x2 = y1 = y2 = 0.0
    for n in range(len(x)):
        xn = x[n]
        yn = b[0]*xn + b[1]*x1 + b[2]*x2 - a[1]*y1 - a[2]*y2
        x2, x1 = x1, xn
        y2, y1 = y1, yn
        y[n] = yn
    return y

Two implementation details matter for reproducibility. First, the pre-warped tangent form keeps the −3 dB point at the specified class frequency even when f_c/f_s is not small; a naive bilinear transform without pre-warping will place the cutoff slightly low. Second, the forward-backward (zero-phase) pass is what most injury-metric workflows expect, because a single causal pass introduces a phase lag that shifts the peak in time and biases window-based metrics such as HIC15. If your reference procedure specifies a causal single pass, document that choice explicitly — the peak values will differ.

Applying the Filter to Acceleration and Force Channels

Load the time history as two columns (time in seconds, signal in g or kN) and resample to a uniform sample rate before filtering. Non-uniform sampling is a common source of silent error: the biquad assumes constant Δt. A minimal, runnable driver:

import numpy as np

# --- example input: uniform 10 kHz head resultant acceleration (g) ---
fs = 10000.0
t  = np.arange(0, 0.060, 1.0/fs)
# synthetic pulse: 120 g half-sine, 8 ms, plus 40 Hz noise
pulse = 120.0 * np.sin(np.pi * t / 0.008) * (t < 0.008)
noise = 6.0 * np.sin(2*np.pi*40*t) + 3.0 * np.sin(2*np.pi*900*t)
a_raw = pulse + noise

a_cfc180 = cfc_filter(a_raw, 180.0, fs)
a_cfc600 = cfc_filter(a_raw, 600.0, fs)

for name, sig in [("raw", a_raw), ("CFC180", a_cfc180), ("CFC600", a_cfc600)]:
    i = int(np.argmax(sig))
    print(f"{name:7s} peak={sig[i]:7.2f} g  t_peak={t[i]*1000:6.2f} ms")

Run this and you will see the raw peak inflated by the 900 Hz component, CFC 600 retaining most of the high-frequency content, and CFC 180 attenuating it. The exact numbers depend on the synthetic pulse, so treat them as expected behavior, not measured data. The same routine applies unchanged to a force channel (e.g., a seat-belt load cell in kN) — only the units and the CFC class change. Belt force channels are typically filtered at CFC 600 or per the test standard's channel table; always confirm the class against the governing procedure rather than defaulting to the acceleration class.

Quantifying the Effect on Injury Metrics

Peak amplitude alone is not the deliverable — the point is how the filter propagates into the metric. Compute the metric on each filtered variant and tabulate the delta. For a head channel, HIC15 uses a 15 ms sliding window; for a chest channel, the 3 ms clip uses a 3 ms window. The window search is identical regardless of filter, so the only variable is the signal fed into it.

def hic15(a, dt, window=0.015):
    n = int(round(window / dt))
    best = 0.0
    for i in range(len(a) - n):
        seg = a[i:i+n]
        val = (seg.mean()) * (window) ** 2.5 * 1.0  # simplified; use integral form below
    # integral form (correct):
    best = 0.0
    for i in range(len(a) - n):
        seg = a[i:i+n]
        integ = np.trapz(seg, dx=dt)
        val = (integ / window) * window ** 2.5
        if val > best:
            best = val
    return best

dt = 1.0 / fs
results = {}
for name, sig in [("raw", a_raw), ("CFC180", a_cfc180), ("CFC600", a_cfc600)]:
    results[name] = {
        "peak_g": float(sig.max()),
        "t_peak_ms": float(t[int(np.argmax(sig))] * 1000),
        "hic15": hic15(sig, dt),
    }
for k, v in results.items():
    print(k, {kk: round(vv, 3) for kk, vv in v.items()})

Expected behavior for the synthetic pulse: the raw signal shows the highest peak and the highest HIC15 because the 900 Hz component adds energy inside the 15 ms window; CFC 600 sits between raw and CFC 180; CFC 180 gives the lowest peak and lowest HIC15. The peak time may shift by a fraction of a millisecond between variants — record it, because a shifted peak changes which window is selected and can move the reported HIC15 by more than the amplitude change alone suggests. Do not report these synthetic values as validation data; they are a sanity check that the filter is doing what the standard requires.

Documentation and Verification Checklist

The deliverable for this case is a short filter-effect record that a reviewer can audit. Populate it with the actual numbers from your run:

Once the table is filled, the effect of CFC filtering on the injury metric is documented with the same rigor as the metric itself: class, sample rate, peak, peak time, and delta from the unfiltered baseline. That record is what makes the filtered result defensible in a report or a correlation study.