1. Import and Standardize Injury Time-History Data

Why standardization must precede overlay

Injury time-history data rarely arrives in a comparable form. A head acceleration channel exported from LS-DYNA is typically a *DATABASE_ASCII or nodout-style text file with a time column in milliseconds and acceleration in g; a PAM-CRASH run may deliver TH curves through its time-history output with time in seconds and acceleration in mm/ms²; a Radioss run commonly writes T01/TH ASCII with time in milliseconds but acceleration in m/s². Chest deflection adds another axis of variation: some solvers report absolute deflection in mm, others report a normalized deflection ratio, and sign conventions differ depending on whether compression is positive or negative. Overlaying these curves without normalization produces plots that look plausible but are physically meaningless — a 1000 g head pulse and a 1000 m/s² head pulse differ by a factor of roughly 9.81, and a deflection curve in ratio units cannot be compared against one in millimeters.

This case establishes the ingestion layer for the whole 30-case series: a single, explicit, auditable pipeline that reads heterogeneous solver exports, maps each channel to a canonical name and unit, resamples everything onto one common time base, and verifies alignment before any curve is drawn. The deliverable is a small, self-contained Python module plus a verification report, not a plot.

Canonical channel and unit contract

Fix the target representation first, because every downstream decision depends on it. The contract used here is deliberately simple and solver-agnostic:

Two rules make the contract enforceable. First, the time base is always milliseconds and always starts at zero; any solver offset is removed during import, not during plotting. Second, the common sampling rate is a single scalar, FS_HZ = 10000.0 (0.1 ms spacing), chosen because it is an integer multiple of the typical 0.1 ms and 0.05 ms solver output steps and is fast enough to preserve the sharp peaks that dominate HIC. Resampling to a rate lower than the native step is forbidden by an explicit guard, since decimation silently attenuates peak acceleration.

Import adapters and unit normalization

Each solver gets a thin adapter that returns a two-column array and a unit descriptor; all unit conversion happens in one place so it can be audited. The example below is complete and runnable with only NumPy and pandas, and it uses synthetic sample data so the logic can be exercised without a solver installation.

import numpy as np
import pandas as pd

G = 9.80665                      # m/s^2 per g
FS_HZ = 10000.0                  # common sampling rate, Hz
DT_MS = 1000.0 / FS_HZ           # 0.1 ms

# Unit conversion factors to canonical units (time -> ms, accel -> g, defl -> mm)
ACCEL_TO_G = {"g": 1.0, "m/s^2": 1.0 / G, "mm/ms^2": 1.0 / G * 1.0}
DEFL_TO_MM = {"mm": 1.0, "m": 1000.0, "ratio": None}   # ratio needs a reference depth

def load_two_column(path, time_unit, value_unit, value_kind):
    """Read a whitespace/comma separated 2-column export and normalize units."""
    df = pd.read_csv(path, sep=r"[,\s]+", engine="python",
                     comment="#", header=None, names=["t", "v"])
    df = df.apply(pd.to_numeric, errors="coerce").dropna()
    t = df["t"].to_numpy()
    v = df["v"].to_numpy()
    if time_unit == "s":
        t = t * 1000.0
    elif time_unit != "ms":
        raise ValueError(f"unsupported time unit: {time_unit}")
    if value_kind == "accel":
        v = v * ACCEL_TO_G[value_unit]
    elif value_kind == "defl":
        factor = DEFL_TO_MM[value_unit]
        if factor is None:
            raise ValueError("ratio deflection requires a reference depth")
        v = v * factor
    else:
        raise ValueError(f"unknown value_kind: {value_kind}")
    return t, v

def to_common_base(t, v, fs_hz=FS_HZ):
    """Shift to t=0, enforce monotonic time, resample onto the common grid."""
    order = np.argsort(t)
    t, v = t[order], v[order]
    t = t - t[0]
    if np.any(np.diff(t) <= 0):
        raise ValueError("duplicate or non-increasing time stamps")
    native_dt = np.median(np.diff(t))
    if native_dt > DT_MS * 1.001:
        raise ValueError(
            f"native step {native_dt:.4f} ms is coarser than target {DT_MS:.4f} ms; "
            "refusing to upsample-interpolate a coarse signal")
    n = int(np.floor(t[-1] / DT_MS)) + 1
    grid = np.arange(n) * DT_MS
    return grid, np.interp(grid, t, v)

The guard on native_dt is the single most important line in the module. Interpolating a 0.5 ms signal onto a 0.1 ms grid does not add information and can create the false impression of a well-resolved peak; rejecting it forces the analyst to either accept a coarser common rate or re-export at higher resolution. The ratio deflection case is intentionally left as an explicit error rather than a silent guess, because converting a normalized deflection to millimeters requires the chest depth used by that specific model, which must be supplied by the analyst.

Alignment and unit verification before overlay

After normalization, three checks decide whether the curves may be overlaid. They are cheap, deterministic, and catch the failure modes that actually occur in practice.

1. Time-base identity. All channels must share the same grid length and the same DT_MS; a mismatch means one adapter used a different rate.
2. Event alignment. Locate the time of peak acceleration, t_peak, on each head channel. For simulations of the same impact configuration, t_peak should agree within one or two native steps. A systematic offset of tens of milliseconds usually indicates a solver output offset or a different trigger time, not a physical difference.
3. Unit plausibility. Compare peak magnitudes against expected ranges for the load case. A frontal head pulse in the 50–150 g range and a chest deflection in the 20–60 mm range are typical for occupant sled tests; a head peak near 10 or near 1000 is a strong signal of a missed or double-applied unit conversion.

Expected output and pass criteria

Running the module over the three sample exports should produce a single tidy table with columns time, head_accel_resultant, and chest_deflection, one row per common-grid sample, plus a verification summary. For the synthetic data used in this example, the expected observable outcomes are: identical grid length across all three sources, t_peak agreement within 0.2 ms, and head peaks that fall inside the 50–150 g plausibility band after conversion. If any of these fail, the summary should name the failing gate and the offending source file rather than proceeding to overlay.

The practical rule to carry forward: a curve is only overlay-ready when its time base, unit, sign convention, and event timing have each been explicitly checked and recorded. Anything less turns a comparison plot into a source of silent error.