2. Python environment for crash data analysis

Crash signal processing sits at the intersection of mechanical engineering, signal theory, and data engineering. Before a single acceleration trace can be filtered or a HIC value computed, the analyst needs a Python environment that is stable, versioned, and capable of reading the file formats produced by crash-test data acquisition systems. This chapter establishes that foundation: which libraries matter and why, how to install them reproducibly, how to load the three formats you will encounter most often (ISO MME, DIAdem, and CSV), and how to structure a workflow so that a result computed today can be recomputed identically six months from now.

Why the Environment Is Part of the Analysis

In crash analysis, a filtered acceleration signal feeds directly into injury criteria such as HIC, NIC, or chest deflection. A change in the filter implementation, the sampling-rate handling, or the unit conversion can shift a HIC value by several percent — enough to change a pass/fail judgment. Reproducibility is therefore not administrative hygiene; it is a technical requirement. The practical consequences are concrete: pin your library versions, record the Python interpreter version, and keep raw data immutable. Everything else in this chapter follows from those three rules.

The Core Stack and What Each Library Does

The four libraries form a natural division of labour, and understanding that division prevents the common mistake of reaching for pandas when NumPy is the right tool, or vice versa.

• NumPy — the numerical substrate. Crash signals are dense, uniformly sampled arrays; NumPy's ndarray and its vectorized operations (cumulative integration via np.cumsum, differentiation via np.gradient, resampling via np.interp) are what make channel processing fast and predictable. NumPy also defines the dtype discipline that keeps float32 sensor data from silently degrading to float64 or vice versa.
• SciPy — the signal-processing and interpolation layer. scipy.signal.butter plus sosfiltfilt implements the zero-phase Butterworth filters required by SAE J211 / ISO 6487 channel classes; scipy.signal.resample_poly handles rational resampling; scipy.integrate supports the integration steps behind velocity and displacement derivation.
• pandas — the tabular and metadata layer. Channel headers, test metadata, event markers, and multi-channel time series align naturally with DataFrame and DatetimeIndex. pandas is where you join a channel table to a test log, not where you run a 100 kHz filter.
• matplotlib — the verification layer. Every processing step should be visually checked: raw versus filtered overlay, residual plots, and annotated event markers. A filter that looks correct numerically but rings at the boundary is a classic failure that only a plot reveals.

Installing a Reproducible Environment

Use a dedicated environment per project. With conda (recommended when you also need HDF5, netCDF, or vendor SDKs), create and pin:

conda create -n crash-py311 python=3.11 numpy=1.26 scipy=1.11 pandas=2.1 matplotlib=3.8
conda activate crash-py311
conda env export --no-builds > environment.yml

With plain pip, the equivalent is a virtual environment plus a fully pinned requirements file:

python -m venv .venv
source .venv/bin/activate        # Windows: .venv\Scripts\activate
pip install numpy==1.26.4 scipy==1.11.4 pandas==2.1.4 matplotlib==3.8.4
pip freeze > requirements.txt

Commit environment.yml or requirements.txt alongside your analysis code. When a colleague reproduces your HIC calculation, they install from that file and get bit-identical library behaviour. Record the interpreter version in the output report as well — NumPy's random number handling and pandas' default dtypes have changed across minor versions.

Loading the Three Common Crash Data Formats

The formats differ in structure, and each needs a small adapter that returns the same internal representation: a time vector plus a dictionary of named channels, all in SI units with an explicit sampling rate.

ISO MME (ISO 13499, the "Measurement Data Management" exchange format) is a binary container with a structured header describing channels, units, and scaling. It is the standard interchange format between crash laboratories. Python support is typically via a vendor library or a small parser; the key discipline is to read the header first and never assume channel order. A typical adapter returns channel names, units, and the sample interval from the header, then reads the data block.

DIAdem files (.dat, .tdm, .tdms) come from National Instruments tooling widely used in test cells. TDMS is hierarchical: groups contain channels, and each channel carries properties such as unit and description. The nptdms package reads them cleanly:

from nptdms import TdmsFile
tdms = TdmsFile.read("test_0421.tdms")
group = tdms["Channels"]
time = group["Time"].time_track()
acc_x = group["Accel_X"].data   # check .properties["unit"] before use

CSV is the least standardized and therefore the most dangerous. Delimiters, decimal separators (comma versus point), header rows, and unit annotations vary by exporter. Always load with explicit parameters and validate the result:

import pandas as pd
df = pd.read_csv("run_07.csv", sep=";", decimal=",", skiprows=3,
                 names=["time_s", "accel_x_g", "accel_y_g"], header=None)
assert df["time_s"].is_monotonic_increasing, "time axis not sorted"
dt = df["time_s"].diff().median()
print(f"median dt = {dt:.6f} s  ->  {1/dt:.1f} Hz")

The assert and the printed sampling rate are not decoration. A non-monotonic time axis usually means the exporter wrote event rows into the data block; a sampling rate that does not match the test log means a unit or scaling error upstream. Catching both at load time is far cheaper than debugging a HIC value later.

Normalizing to a Common Internal Representation

Whatever the source format, convert immediately to one internal structure. A minimal, effective convention is a dataclass holding the time array, a dict of channel arrays, a dict of units, and the nominal sampling rate. Convert g to m/s² at this boundary (multiply by 9.80665) and record the conversion in the metadata, so that downstream code never has to guess whether a channel is in g or m/s². This single decision eliminates a large share of unit-mismatch bugs in injury calculation.

Reproducible Workflows in Practice

Three habits carry most of the reproducibility weight. First, separate raw from derived data: keep the original MME/TDMS/CSV files read-only, and write every processed artefact to a new path with a hash of the input recorded in its metadata. Second, make the pipeline a script, not a notebook session. Notebooks are excellent for exploration, but the final calculation should live in a module with a main() that takes a config file and writes a result file. Third, log the parameters that affect the numbers: filter class, cutoff frequency, resampling method, and the unit-conversion constant. A short JSON sidecar next to each output is sufficient.

import json, hashlib, pathlib
def file_sha256(path):
    h = hashlib.sha256()
    with open(path, "rb") as f:
        for chunk in iter(lambda: f.read(1 << 20), b""):
            h.update(chunk)
    return h.hexdigest()

meta = {
    "input": str(pathlib.Path("test_0421.tdms")),
    "sha256": file_sha256("test_0421.tdms"),
    "python": "3.11",
    "numpy": "1.26.4", "scipy": "1.11.4",
    "filter": {"type": "butterworth", "order": 4, "fc_hz": 1650, "zero_phase": True},
    "unit_conversion": {"g_to_m_s2": 9.80665}
}
pathlib.Path("test_0421.meta.json").write_text(json.dumps(meta, indent=2))

With this sidecar in place, a reviewer can confirm that the HIC value in the report came from a specific input file, processed with a specific filter, under a specific library version. That is the level of traceability crash analysis demands, and it costs only a few lines of code.

Verifying the Environment Before Real Work

Run a short smoke test after installation: generate a synthetic 10 kHz sine wave, apply the J211-class filter you intend to use, and confirm that the filtered amplitude at the passband is within a fraction of a percent of the input. If that test passes, the numerical stack is behaving. If it fails, the problem is in the environment, not in the crash data — and you have isolated it before spending hours on a real test file. This smoke test belongs in the repository as tests/test_environment.py and should run in CI if the project has one.

The environment described here is deliberately minimal. It gives you array math, filtering, tabular metadata, and plotting — exactly the four capabilities the rest of this course builds on. Later chapters add domain-specific tooling, but every one of them assumes the pinned, validated, traceable foundation established in this chapter.