3. Reading CAE result files with Python

In a Python-based injury extraction workflow, the reader layer is the boundary between solver-native binary data and the tabular structures that downstream logic can reason about. Every subsequent stage — injury criterion calculation, threshold comparison, and Excel reporting — depends on this layer producing consistent, well-typed arrays. The central engineering challenge is that each CAE format was designed for a different purpose: Abaqus .odb is a hierarchical output database, LS-DYNA .d3plot is a compact binary state dump, Radioss/OpenRadioss .h3d is a time-history container, and CSV is a human-readable export with no schema guarantees. A robust reader must normalize all four into a common in-memory contract.

Defining the internal data contract first

Before writing any parser, decide what a "result set" means in your pipeline. A practical contract for injury work is a dictionary keyed by part or node-set name, where each entry holds a time vector and one or more response arrays aligned to it:

{
  "HEAD": {
    "time": np.ndarray,        # shape (n_frames,), seconds
    "accel": np.ndarray,       # shape (n_frames, 3), m/s^2, global XYZ
    "units": {"time": "s", "accel": "m/s^2"}
  },
  "CHEST_BELT": { ... }
}

Two conventions must be fixed here and never changed later: the time origin (solver start time, typically 0.0) and the acceleration sign convention (positive along the global axis, gravity not subtracted). Injury criteria such as HIC and 3 ms clip are sign-sensitive in their integration, so an undocumented flip between readers will silently corrupt results. Record the convention in the result-set metadata and assert it in the reader.

Reading Abaqus .odb files

The .odb format is accessed through the Abaqus Scripting Interface, which requires a licensed Abaqus installation and must run under abaqus python or abaqus cae noGUI=script.py. Pure CPython cannot open it. The typical extraction pattern opens the database, selects a step and frame, and reads a field output:

from odbAccess import openOdb

odb = openOdb(path="run_job.odb", readOnly=True)
step = odb.steps["Impact"]
frame = step.frames[-1]
accel_field = frame.fieldOutputs["A"]
region = odb.rootAssembly.nodeSets["HEAD_SET"]

values = accel_field.getSubset(region=region).values
# values[i].data -> tuple of components; values[i].nodeLabel -> node id

For time histories, iterate step.frames and collect frame.frameValue as the time stamp. Two practical pitfalls dominate: node ordering in getSubset follows the set definition, not numeric order, so build a label-to-index map if you need stable alignment; and field output must actually be requested in the input deck (*OUTPUT, FIELD with the relevant variables), otherwise the key lookup raises KeyError. Because the Abaqus Python interpreter is often Python 2.7, keep the ODB extraction script isolated and serialize results to .npz or JSON for the main pipeline to consume.

Reading LS-DYNA .d3plot files

The .d3plot family is a sequence of binary state files (d3plot, d3plot01, d3plot02, …). The open-source lasso-python library reads them without a solver license:

import lasso
db = lasso.Dyna.D3plot("d3plot")
# db.arrays is keyed by lasso.Dyna.ArrayType
node_ids = db.arrays[lasso.Dyna.ArrayType.NODE_IDS]
disp = db.arrays[lasso.Dyna.ArrayType.NODE_DISPLACEMENT]
times = db.get_timesteps()

Acceleration is usually not stored directly; it is obtained by differentiating the displacement history twice, or by reading NODE_ACCELERATION if the deck requested it. Numerical differentiation amplifies noise, so apply a consistent filter (e.g., SAE J211 CFC 1000 for head impact) before computing injury metrics, and document the filter class alongside the data. Part and node-set membership comes from db.get_part_list() and the keyword file; cross-reference node IDs rather than positional indices, because LS-DYNA renumbers internally.

Reading Radioss/OpenRadioss .h3d files

.h3d is a time-history container read with pyradi or the OpenRadioss Python API. The access pattern is request-based: you declare which variables and entity sets you want, then read per time step.

import pyradi.rdreader as rd
rd.Reader("run.h3d")
rd.QueryCurves("HEAD_SET", "ACCX")
times, values = rd.ReadCurve("HEAD_SET", "ACCX")

The main risk with .h3d is unit ambiguity: Radioss decks frequently use millimetres, milliseconds, and tonnes, so acceleration may arrive in mm/ms^2. Convert explicitly at the reader boundary and store the source units in metadata rather than assuming SI.

CSV exports and schema validation

CSV is the least structured input and therefore needs the most defensive handling. A typical export has a header row, a time column, and one column per channel, sometimes with a units row and comment lines prefixed by #. Read with explicit dtype control and validate the schema before use:

import pandas as pd

df = pd.read_csv("channels.csv", comment="#", skip_blank_lines=True)
required = {"time", "HEAD_ACC_X", "HEAD_ACC_Y", "HEAD_ACC_Z"}
missing = required - set(df.columns)
if missing:
    raise ValueError(f"CSV missing channels: {missing}")
df = df.sort_values("time").drop_duplicates("time")
assert df["time"].is_monotonic_increasing

Common failure symptoms and their causes are worth tabulating, because they recur across projects:

Normalizing into one result set

Each reader should end by emitting the same contract dictionary, with units converted to SI and the sign convention asserted. A thin dispatcher selects the reader by file extension and returns the normalized object, so the injury-calculation stage never branches on format. The diagram below shows this flow: four heterogeneous sources converge on a single normalization boundary, and only the normalized result set crosses into the analysis layer.

Verification before handing off

A reader is only trustworthy once it has been checked against a known reference. For each format, run a small sanity suite: confirm the time vector starts at the solver start time and ends at the termination time; confirm the number of frames matches the requested output interval; and compare a single channel against the solver's own plot or a CSV export of the same run. A practical acceptance check is that the peak resultant acceleration from the .d3plot path and the CSV path agree within 1–2%, with any larger gap traced to filtering or unit conversion rather than accepted. Once these checks pass, the normalized result set can be cached to .npz so that repeated injury-calculation runs do not re-parse large binaries.