3. Zero-Phase Low-Pass Filtering for Noise Reduction

Why Zero-Phase Matters for Crash Pulses

A crash acceleration pulse is judged by its peak value, its time of occurrence, and the shape of the leading edge. Any causal (real-time) IIR filter delays the signal, and the delay is frequency-dependent: the group delay of a low-pass Butterworth filter is largest near the cutoff and smaller in the passband. When you filter a sharp impact transient with a causal low-pass, the peak is not only attenuated but also shifted in time, and the pulse can be visibly skewed. For injury-criterion work (HIC, 3 ms clip, chest deflection corridors) that time shift and shape distortion are unacceptable, because the metrics are defined on the filtered signal's timing.

The standard remedy is forward–backward filtering: run the same digital filter over the data once in the forward direction, then again over the reversed result. The magnitude response is squared (so the effective cutoff moves), but the phase response cancels exactly, giving zero phase distortion and no net time shift. In Python this is scipy.signal.filtfilt. The price is a non-causal result: the output at sample n depends on samples after n, so the first and last few samples are contaminated by edge transients and must be handled with padding.

Designing the Filter for a Crash Pulse

Two parameters drive the design: the cutoff frequency fc and the order N. Crash accelerometers are typically sampled at 10 kHz or higher, and the mechanical content of interest (occupant-relevant acceleration) sits well below a few hundred hertz. A common teaching choice is a 4th-order Butterworth with fc = 300 Hz at fs = 10 kHz. Because filtfilt applies the filter twice, the effective magnitude response is |H(f)|², i.e. the −3 dB point of the combined filter is lower than the design cutoff. To land the combined −3 dB point at the intended frequency, pre-warp the design cutoff:

fc_eff = fc / sqrt(2**(1/N) - 1)   # for filtfilt, N = order of the single-pass filter

For N = 4 and fc = 300 Hz, fc_eff ≈ 300 / sqrt(2**0.25 − 1) ≈ 300 / 0.435 ≈ 690 Hz. Designing the single-pass filter at ~690 Hz and applying it forward–backward yields a combined −3 dB point near 300 Hz. If you skip this step the combined response is over-attenuated and the peak drops more than intended — a classic source of "my filtered peak is too low" complaints.

Reference Implementation

The snippet below is self-contained: it synthesizes a representative crash pulse (a half-sine impact plus a decaying oscillation and additive high-frequency noise), designs the pre-warped Butterworth filter, applies zero-phase filtering, and checks peak preservation and timing. It runs with only NumPy and SciPy.

import numpy as np
from scipy import signal

# --- 1. Synthetic crash pulse (example data, not a measured record) ---
fs = 10_000.0                 # Hz
t = np.arange(0, 0.12, 1/fs)  # 120 ms window
impact = 80.0 * np.sin(np.pi * t / 0.04) * (t < 0.04)          # 80 g half-sine, 40 ms
ring   = 12.0 * np.exp(-t/0.02) * np.sin(2*np.pi*450*t)         # 450 Hz structural ring
noise  = 1.5 * np.random.default_rng(0).standard_normal(t.size) # broadband sensor noise
raw    = impact + ring + noise

# --- 2. Zero-phase Butterworth design with cutoff pre-warp ---
fc, N = 300.0, 4
fc_design = fc / np.sqrt(2**(1/N) - 1)          # ~690 Hz for N=4, fc=300
b, a = signal.butter(N, fc_design, btype='low', fs=fs)

# --- 3. Forward-backward filtering (zero phase) ---
# padlen must exceed the filter's transient length; default is fine here,
# but set it explicitly for reproducibility.
padlen = 3 * max(len(a), len(b))
filt = signal.filtfilt(b, a, raw, padlen=padlen)

# --- 4. Verification: peak value and peak time ---
i_raw, i_filt = np.argmax(raw), np.argmax(filt)
print(f"raw  peak = {raw[i_raw]:.2f} g at t = {t[i_raw]*1e3:.2f} ms")
print(f"filt peak = {filt[i_filt]:.2f} g at t = {t[i_filt]*1e3:.2f} ms")
print(f"peak shift = {(t[i_filt]-t[i_raw])*1e3:+.3f} ms")
print(f"peak retention = {filt[i_filt]/raw[i_raw]*100:.1f} %")

With the seed fixed, the raw peak is dominated by the 450 Hz ring riding on the half-sine, so it sits slightly above 80 g and slightly off the half-sine apex. After filtering, the ring and noise are attenuated and the peak should land close to the half-sine apex near 20 ms, with a peak shift of a few tenths of a millisecond at most (residual shift comes from the ring's contribution to the raw peak, not from filter delay). The peak retention figure is the number to watch: for a 300 Hz combined cutoff and a 40 ms half-sine, retention should be high (roughly 95–99 %), because the pulse's dominant energy is far below cutoff. If retention drops below ~90 %, either the cutoff is too low for the pulse width or the pre-warp was omitted.

Reading the Result: What the Filter Should and Should Not Change

The figure below contrasts the raw and zero-phase filtered pulse and shows the two failure modes to check for: peak attenuation (cutoff too low) and edge ringing (insufficient padding).

Verification Checklist and Failure Symptoms

Run these checks on every zero-phase filtering job before the filtered pulse is used for metrics:

• Peak retention: filtered peak / raw peak. For a 40 ms half-sine at 300 Hz combined cutoff, expect roughly 95–99 %. A large drop means the cutoff is too low or the pre-warp was skipped.
• Peak timing: the filtered peak should sit at the half-sine apex (≈20 ms here). A consistent shift of several milliseconds indicates a causal filter was used by mistake, or the wrong axis/sign convention.
• Edge behavior: inspect the first and last padlen samples. Visible overshoot or a non-zero start means padding is too short; increase padlen or use method='gust' in filtfilt.
• Combined response: plot |H(f)|² and confirm the −3 dB point matches the intended fc. This is the single most common documentation error.

Once these four checks pass, the filtered pulse is safe to feed into downstream metric computation: the peak value is preserved within the expected band, the timing is unshifted, and the high-frequency ring and sensor noise have been removed without introducing artificial transients at the record boundaries.