"""Plot daily C3/C4 condition averages for the TCH recordings."""
from __future__ import annotations
import csv
import re
import shutil
from datetime import datetime
from pathlib import Path
import matplotlib.pyplot as plt
import numpy as np
ROOT = Path(__file__).resolve().parent.parent
CACHE = Path(__file__).resolve().parent / "data_cache"
OUTPUT = ROOT / "輸出呈現用" / "20260721_中榮data_每日C3C4平均波形"
CHANNELS = ("Fp1", "Fp2", "Fz", "C3", "C4", "Pz", "O1", "O2")
COLORS = {0: "#2c7a7b", 1: "#d46a6a"}
NAMES = {0: "Left (-1)", 1: "Right (1)"}
def session_dates() -> dict[str, str]:
result = {}
for path in (ROOT / "中榮data").glob("*.npz"):
matched = re.match(r"(s\d_\d)_(\d{4})_data", path.stem)
if matched:
result[matched.group(1)] = datetime.strptime(matched.group(2), "%m%d").replace(year=2026).strftime("%Y-%m-%d")
return result
def plot_session(session: str, day: str, values: np.ndarray, labels: np.ndarray, output: Path) -> dict[str, object]:
time = np.arange(values.shape[-1]) / 200 - 1
# Remove each trial's pre-event level; do not alter within-trial dynamics.
corrected = values - values[:, :, 100:200].mean(axis=-1, keepdims=True)
figure, axes = plt.subplots(2, 1, figsize=(10.5, 6.3), sharex=True, constrained_layout=True)
counts = {}
for axis, (channel, index) in zip(axes, (("C3", CHANNELS.index("C3")), ("C4", CHANNELS.index("C4")))):
for label in (0, 1):
trials = corrected[labels == label, index]
counts[label] = len(trials)
mean = trials.mean(axis=0)
sem = trials.std(axis=0, ddof=1) / np.sqrt(len(trials))
axis.plot(time, mean, color=COLORS[label], lw=1.8, label=f"{NAMES[label]} (n={len(trials)})")
axis.fill_between(time, mean - sem, mean + sem, color=COLORS[label], alpha=.18, linewidth=0)
axis.axvspan(-1, 0, color="#e8edf2", alpha=.65, label="Fixation" if channel == "C3" else None)
axis.axvspan(0, 2, color="#edf6f0", alpha=.6, label="ME" if channel == "C3" else None)
axis.axvline(0, color="#526273", lw=1, ls="--")
axis.set_ylabel(f"{channel}\nbaseline-corrected nV")
axis.grid(axis="y", alpha=.22)
axis.spines[["top", "right"]].set_visible(False)
axes[0].legend(ncol=4, frameon=False, loc="upper right")
axes[1].set_xlabel("Time relative to event (s)")
figure.suptitle(f"{session} | {day} | Daily valid-trial average, 1-40 Hz IIR", fontsize=14)
filename = f"{session}_{day.replace('-', '')}_C3C4.png"
figure.savefig(output / "img" / filename, dpi=160, bbox_inches="tight")
plt.close(figure)
return {"session_id": session, "date": day, "left_trials": counts[0], "right_trials": counts[1], "image": f"img/{filename}"}
def main() -> None:
(OUTPUT / "img").mkdir(parents=True, exist_ok=True)
(OUTPUT / "code").mkdir(exist_ok=True)
dates = session_dates()
rows = []
for path in sorted((CACHE / "sessions").glob("*_x.npy")):
session = path.stem[:-2]
if session not in dates:
continue
values = np.load(path, allow_pickle=False)
labels = np.load(path.with_name(f"{session}_y.npy"), allow_pickle=False)
rows.append(plot_session(session, dates[session], values, labels, OUTPUT))
rows.sort(key=lambda row: (row["date"], row["session_id"]))
with (OUTPUT / "daily_summary.csv").open("w", encoding="utf-8-sig", newline="") as handle:
writer = csv.DictWriter(handle, fieldnames=rows[0].keys()); writer.writeheader(); writer.writerows(rows)
shutil.copy2(Path(__file__), OUTPUT / "code" / Path(__file__).name)
cards = "".join(f"{r['session_id']} | {r['date']}
有效 trials：左腳 {r['left_trials']}、右腳 {r['right_trials']}
![\"{r['session_id']}](\"{r['image']}\")
" for r in rows)
page = f"""中榮資料每日 C3/C4 平均波形中榮資料：每日 C3/C4 平均波形
每張圖是單一受試者單一日期的有效 trial 平均，不跨受試者混合。原始資料先連續做 1–40 Hz fourth-order zero-phase SOS IIR，降到 200 Hz；每個 trial 再以 -0.5–0 s baseline correction。陰影為 trial-level SEM。
注意：s1_4 與 s2_1 同為 2026-06-24、s1_5 與 s2_2 同為 2026-07-01，但分屬不同受試者，因此分圖呈現而不直接平均。
{cards}"""
(OUTPUT / "index.html").write_text(page, encoding="utf-8")
print(f"Wrote {OUTPUT}")
if __name__ == "__main__":
main()