Project existing simulation field outputs onto an fsaverage template.
Complementary to :mod:tit.source.forward: where the forward pipeline lets a
researcher reconstruct EEG sources onto fsaverage, this pipeline puts the
stimulation fields on the same grid, so source activity and TI exposure can be
compared vertexwise.
For an already-completed TI simulation it reads the run's own central surface
(<TI|mTI>/mesh/surfaces/<sim>_<TI|mTI>_central.msh) and channel meshes and
morphs the requested scalar fields to fsaverage:
TI_max -- orientation-maximized TI envelope |E| (central surface)TI_normal -- directional TI envelope along the cortical normal (central surface)hf_peak -- peak carrier field, max over channel sign combinations (Cassarà 2025, safety)hf_rms -- RMS carrier field over time, sqrt(tit.calc.hf_power) (V/m)
This path covers both standard TI and mTI. The two differ only in what the
finished run is called on disk, and each resolver here takes either:
:func:_ti_max_overlay reads TI_max from TI/mesh/surfaces or the same
quantity spelled mTI_max under mTI/mesh/surfaces, :func:_ti_normal_surface
reads TI_normal from the same surface, and :func:_carrier_volume_meshes collects every
carrier volume mesh under high_Frequency/ -- TDCS_1/TDCS_2 for TI, and
TDCS_A..TDCS_Z for mTI's N channels. hf_peak/hf_rms are then computed
over all N channel fields, not a hardcoded two.
TI_max / TI_normal are read from the run's central surface;
hf_peak / hf_rms are interpolated from the carrier volume meshes
using the same formulas as the volume output (:func:tit.calc.hf_peak,
:func:tit.calc.hf_power).
Unlike SimNIBS's native map_to_fsavg (which only runs at simulation time and
only emits TI_max), this works post-hoc on any finished simulation and on
the derived TI_normal / hf_peak / hf_rms quantities.
Output, in <sim>/fsaverage/: sub-<id>_sim-<sim>_space-fsaverage<N>_fields.msh
-- the fsaverage central surfaces (lh nodes, then rh) with each field as a node
field, openable in Gmsh or the Viewer -- plus its .msh.opt Gmsh view and a
same-stem .json provenance sidecar (subject_id, simulation,
carrier_only, fields).
Runs under simnibs_python (reads SimNIBS meshes)::
simnibs_python -m tit.source fsavg_config.json
See Also
tit.source.config.FsavgMapConfig : Parameter object.
tit.source.forward.prepare_forward : The companion forward pipeline.
write_fsaverage_msh
Write [lh; rh] fsaverage values as a viewable SimNIBS surface mesh.
Follows SimNIBS's own fsavg_overlays/<name>_fsavg.msh (join_and_write in
:func:simnibs.utils.transformations.middle_gm_interpolation): the template
central surfaces joined lh then rh, one node field per quantity, and the Gmsh
.opt view beside it. The template shares its vertex numbering with the
sphere :func:cross_subject_map resamples onto, so value i lands on node
i (tests/numerical/test_fsaverage_msh.py). .msh carries no free
metadata, so provenance goes to a same-stem .json sidecar (BIDS style).
Source code in tit/source/fsaverage.py
| def write_fsaverage_msh(
path: Path, maps: dict[str, np.ndarray], spacing: int, provenance: dict
) -> None:
"""Write ``[lh; rh]`` fsaverage values as a viewable SimNIBS surface mesh.
Follows SimNIBS's own ``fsavg_overlays/<name>_fsavg.msh`` (``join_and_write`` in
:func:`simnibs.utils.transformations.middle_gm_interpolation`): the template
central surfaces joined lh then rh, one node field per quantity, and the Gmsh
``.opt`` view beside it. The template shares its vertex numbering with the
sphere :func:`cross_subject_map` resamples onto, so value ``i`` lands on node
``i`` (``tests/numerical/test_fsaverage_msh.py``). ``.msh`` carries no free
metadata, so *provenance* goes to a same-stem ``.json`` sidecar (BIDS style).
"""
from simnibs.mesh_tools import mesh_io
surfs = mesh_io.load_fsaverage_template("central", spacing)
mesh = surfs["lh"].join_mesh(surfs["rh"])
for name, values in maps.items():
mesh.add_node_field(np.asarray(values, dtype=float), name)
path.parent.mkdir(parents=True, exist_ok=True)
mesh_io.write_msh(mesh, str(path))
mesh.view(visible_fields=next(iter(maps))).write_opt(str(path))
path.with_suffix(".json").write_text(
json.dumps({**provenance, "fields": list(maps)}, indent=2) + "\n"
)
|
read_fsaverage_fields
Every node field of a mesh written by :func:write_fsaverage_msh.
Uses SimNIBS's reader: every caller (tit.stats, tit.plotting,
scripts) runs under simnibs_python, like the writer.
Source code in tit/source/fsaverage.py
| def read_fsaverage_fields(path: Path) -> dict[str, np.ndarray]:
"""Every node field of a mesh written by :func:`write_fsaverage_msh`.
Uses SimNIBS's reader: every caller (``tit.stats``, ``tit.plotting``,
scripts) runs under ``simnibs_python``, like the writer.
"""
from simnibs.mesh_tools import mesh_io
mesh = mesh_io.read_msh(str(path))
return {
name: np.asarray(data.value, dtype=np.float64).reshape(-1)
for name, data in mesh.field.items()
}
|
read_fsaverage_field
One node field of a mesh written by :func:write_fsaverage_msh.
Source code in tit/source/fsaverage.py
| def read_fsaverage_field(path: Path, field: str) -> np.ndarray:
"""One node field of a mesh written by :func:`write_fsaverage_msh`."""
fields = read_fsaverage_fields(path)
if field not in fields:
raise KeyError(
f"{path.name} has no field {field!r}; available: {sorted(fields)}"
)
return fields[field]
|
project_subject
Project one (subject, simulation) and write it as an fsaverage .msh.
Returns (subject_id, status, message) where status is one of
{"ok", "partial", "cached", "failed"} so batch runs can record and
continue; "partial" means some requested fields could not be projected.
Source code in tit/source/fsaverage.py
| def project_subject(
subject_id: str,
sim: str,
cfg: FsavgMapConfig,
) -> tuple[str, str, str]:
"""Project one (subject, simulation) and write it as an fsaverage ``.msh``.
Returns ``(subject_id, status, message)`` where status is one of
``{"ok", "partial", "cached", "failed"}`` so batch runs can record and
continue; "partial" means some requested fields could not be projected.
"""
from dataclasses import replace
from tit.source.config import VALID_FSAVG_FIELDS
from tit.tools.msh_fields import msh_field_names
pm = get_path_manager()
out_path = Path(pm.sim_fsaverage_fields(subject_id, sim, cfg.fsaverage_spacing))
try:
carrier_only = _recorded_carrier_only(pm, subject_id, sim)
if carrier_only:
# A carrier-only control has no TI envelope: TI_max / TI_normal cannot
# exist. The default (all-fields) request means "what this run has";
# an explicit selection naming a TI field is refused.
ti = {"TI_max", "TI_normal"} & set(cfg.fields)
if ti and tuple(cfg.fields) != VALID_FSAVG_FIELDS:
raise ValueError(
f"{sorted(ti)} cannot be projected: {sim} is a carrier-only run "
"with no TI envelope; request hf_peak / hf_rms."
)
cfg = replace(cfg, fields=tuple(f for f in cfg.fields if f not in ti))
if out_path.exists() and not cfg.overwrite:
names = msh_field_names(out_path) # None: could not tell, keep the file
if names is None or set(cfg.fields) <= names:
return subject_id, "cached", out_path.name
maps = _compute_fields(pm, subject_id, sim, cfg, carrier_only)
except Exception as exc: # noqa: BLE001 - record per-subject and continue
return subject_id, "failed", repr(exc)
write_fsaverage_msh(
out_path,
maps,
cfg.fsaverage_spacing,
{"subject_id": subject_id, "simulation": sim, "carrier_only": carrier_only},
)
medians = ", ".join(f"{k} med={np.median(v):.3f}" for k, v in maps.items())
missing = [f for f in cfg.fields if f not in maps]
if missing:
return subject_id, "partial", f"{medians} -> {out_path.name}; missing {missing}"
return subject_id, "ok", f"{medians} -> {out_path.name}"
|
project_fields_to_fsaverage
Project field outputs to fsaverage for many (subject, simulation) pairs.
Parameters
subjects : list of (str, str)
(subject_id, simulation_name) pairs to project.
cfg : FsavgMapConfig
Field selection, fsaverage spacing, worker count, overwrite flag.
Returns
list of (str, str, str)
Per-pair (subject_id, status, message) results.
Source code in tit/source/fsaverage.py
| def project_fields_to_fsaverage(
subjects: list[tuple[str, str]],
cfg: FsavgMapConfig,
) -> list[tuple[str, str, str]]:
"""Project field outputs to fsaverage for many (subject, simulation) pairs.
Parameters
----------
subjects : list of (str, str)
``(subject_id, simulation_name)`` pairs to project.
cfg : FsavgMapConfig
Field selection, fsaverage spacing, worker count, overwrite flag.
Returns
-------
list of (str, str, str)
Per-pair ``(subject_id, status, message)`` results.
"""
workers = max(1, min(cfg.workers, len(subjects))) if subjects else 1
logger.info("Projecting %d simulation(s) with %d worker(s)", len(subjects), workers)
results: list[tuple[str, str, str]] = []
if workers == 1:
for subject_id, sim in subjects:
res = project_subject(subject_id, sim, cfg)
_log_result(res)
results.append(res)
else:
with ProcessPoolExecutor(max_workers=workers) as pool:
futures = {
pool.submit(project_subject, sid, sim, cfg): (sid, sim)
for sid, sim in subjects
}
for future in as_completed(futures):
res = future.result()
_log_result(res)
results.append(res)
return results
|