Skills › Scientific Agent Skills › neurokit2 › plan_epochs.py
plan_epochs.py
290 lines · skills/neurokit2/scripts/plan_epochs.py at 5f2bc03
1#!/usr/bin/env python3
2"""Plan sample-exact event epochs without loading NeuroKit2."""
4from __future__ import annotations
6import argparse
7import math
8from typing import Any
10from _common import (
11 MAX_CSV_BYTES,
12 MAX_ROWS,
13 CliError,
14 checked_input_file,
15 emit_json,
16 finite_float,
17 read_numeric_columns,
18 require_deidentified,
19 run_cli,
20)
22MAX_EVENTS = 100_000
25def _parse_events(value: str) -> list[float]:
26 pieces = [piece.strip() for piece in value.split(",")]
27 if not pieces or any(not piece for piece in pieces):
28 raise CliError("--events must be a comma-separated numeric list")
29 if len(pieces) > MAX_EVENTS:
30 raise CliError(f"event count exceeds {MAX_EVENTS}")
31 events: list[float] = []
32 for index, piece in enumerate(pieces, start=1):
33 try:
34 event = float(piece)
35 except ValueError as exc:
36 raise CliError(f"event {index} is not numeric") from exc
37 if not math.isfinite(event) or event < 0:
38 raise CliError(f"event {index} must be finite and nonnegative")
39 events.append(event)
40 return events
43def _sample_offset(seconds: float, sampling_rate: float, *, name: str) -> int:
44 exact = seconds * sampling_rate
45 if not math.isfinite(exact):
46 raise CliError(f"{name} must map to a finite sample offset")
47 rounded = round(exact)
48 if not math.isclose(exact, rounded, rel_tol=0.0, abs_tol=1e-8):
49 raise CliError(
50 f"{name}={seconds} s is not sample-aligned at {sampling_rate} Hz "
51 f"({exact} samples)"
52 )
53 return int(rounded)
56def _event_samples(
57 values: list[float],
58 *,
59 unit: str,
60 sampling_rate: float,
61) -> list[int]:
62 samples: list[int] = []
63 for index, value in enumerate(values, start=1):
64 if not math.isfinite(value) or value < 0:
65 raise CliError(f"event {index} must be finite and nonnegative")
66 if unit == "samples":
67 rounded = round(value)
68 if not math.isclose(value, rounded, rel_tol=0.0, abs_tol=1e-9):
69 raise CliError(f"event {index} is not an integer sample index: {value}")
70 sample = int(rounded)
71 else:
72 sample = _sample_offset(value, sampling_rate, name=f"event {index} onset")
73 samples.append(sample)
74 return samples
77def plan(
78 events: list[int],
79 *,
80 recording_samples: int,
81 sampling_rate: float,
82 epoch_start_s: float,
83 epoch_end_s: float,
84 baseline_start_s: float | None,
85 baseline_end_s: float | None,
86 boundary_policy: str,
87) -> dict[str, Any]:
88 start_offset = _sample_offset(epoch_start_s, sampling_rate, name="--epoch-start")
89 end_offset = _sample_offset(epoch_end_s, sampling_rate, name="--epoch-end")
90 if end_offset <= start_offset:
91 raise CliError("--epoch-end must be after --epoch-start")
92 if recording_samples < 1 or recording_samples > MAX_ROWS:
93 raise CliError(f"--recording-samples must be between 1 and {MAX_ROWS}")
95 baseline: dict[str, Any] | None = None
96 if (baseline_start_s is None) != (baseline_end_s is None):
97 raise CliError("provide both --baseline-start and --baseline-end")
98 if baseline_start_s is not None and baseline_end_s is not None:
99 baseline_start = _sample_offset(
100 baseline_start_s, sampling_rate, name="--baseline-start"
101 )
102 baseline_end = _sample_offset(
103 baseline_end_s, sampling_rate, name="--baseline-end"
104 )
105 if not start_offset <= baseline_start < baseline_end <= end_offset:
106 raise CliError("baseline must be a nonempty subset of the epoch")
107 if baseline_end > 0:
108 raise CliError("baseline must not extend after event onset")
109 baseline = {
110 "end_offset_samples_exclusive": baseline_end,
111 "end_s_exclusive": baseline_end_s,
112 "sample_count": baseline_end - baseline_start,
113 "start_offset_samples": baseline_start,
114 "start_s": baseline_start_s,
115 }
117 rows: list[dict[str, Any]] = []
118 invalid = 0
119 for ordinal, onset in enumerate(events, start=1):
120 if isinstance(onset, bool) or not isinstance(onset, int):
121 raise CliError("event onsets must be integer sample indices")
122 if not 0 <= onset < recording_samples:
123 raise CliError(f"event {ordinal} onset is outside the recording")
124 start = onset + start_offset
125 end = onset + end_offset
126 window_samples = end_offset - start_offset
127 before = min(window_samples, max(0, -start))
128 after = min(window_samples, max(0, end - recording_samples))
129 complete = before == 0 and after == 0
130 if not complete:
131 invalid += 1
132 if boundary_policy == "drop" and not complete:
133 continue
134 rows.append(
135 {
136 "complete": complete,
137 "end_sample_exclusive": end,
138 "event_ordinal": ordinal,
139 "onset_sample": onset,
140 "pad_after_samples": after,
141 "pad_before_samples": before,
142 "start_sample": start,
143 }
144 )
145 if boundary_policy == "error" and invalid:
146 raise CliError(f"{invalid} event(s) would cross recording boundaries")
148 warnings: list[str] = []
149 if events != sorted(events):
150 warnings.append("event onsets are not sorted")
151 if len(events) != len(set(events)):
152 warnings.append("duplicate event onsets are present")
153 if invalid:
154 warnings.append(
155 f"{invalid} event(s) cross a boundary; NeuroKit2 0.2.13 pads "
156 "non-integer signal columns with NaN"
157 )
159 return {
160 "baseline": baseline,
161 "boundary_policy": boundary_policy,
162 "event_count_input": len(events),
163 "event_count_output": len(rows),
164 "events": rows,
165 "indexing": {
166 "event_onsets": "zero-based sample indices",
167 "window": "[start_sample, end_sample_exclusive)",
168 },
169 "neurokit2_epoch_note": (
170 "epochs_create() slices the end sample exclusively but labels its "
171 "floating time index with an inclusive epochs_end endpoint. Its "
172 "baseline_correction=True subtracts the mean from epoch start through "
173 "t=0; use manual correction for a narrower baseline."
174 ),
175 "recording_samples": recording_samples,
176 "sampling_rate_hz": sampling_rate,
177 "window": {
178 "end_offset_samples_exclusive": end_offset,
179 "end_s_exclusive": epoch_end_s,
180 "sample_count": end_offset - start_offset,
181 "start_offset_samples": start_offset,
182 "start_s": epoch_start_s,
183 },
184 "warnings": warnings,
185 }
188def build_parser() -> argparse.ArgumentParser:
189 parser = argparse.ArgumentParser(
190 description=(
191 "Plan zero-based, sample-exact event epochs and boundary handling. "
192 "No signal values or NeuroKit2 imports are required."
193 )
194 )
195 source = parser.add_mutually_exclusive_group(required=True)
196 source.add_argument("--events", help="comma-separated onsets")
197 source.add_argument("--events-csv", help="bounded local CSV containing onsets")
198 parser.add_argument("--onset-column", default="onset")
199 parser.add_argument(
200 "--event-unit", choices=("samples", "seconds"), default="samples"
201 )
202 parser.add_argument("--sampling-rate", type=float, required=True, help="Hz")
203 parser.add_argument("--recording-samples", type=int, required=True)
204 parser.add_argument("--epoch-start", type=float, required=True, help="seconds")
205 parser.add_argument("--epoch-end", type=float, required=True, help="seconds")
206 parser.add_argument("--baseline-start", type=float, help="seconds")
207 parser.add_argument("--baseline-end", type=float, help="seconds")
208 parser.add_argument(
209 "--boundary-policy",
210 choices=("report", "drop", "error"),
211 default="report",
212 )
213 parser.add_argument("--root", default=".")
214 parser.add_argument(
215 "--deidentified",
216 action="store_true",
217 help="required with --events-csv",
218 )
219 parser.add_argument("--max-events", type=int, default=MAX_EVENTS)
220 parser.add_argument("--output", help="optional local .json report")
221 parser.add_argument("--force", action="store_true")
222 return parser
225def main() -> None:
226 args = build_parser().parse_args()
227 sampling_rate = finite_float(
228 args.sampling_rate,
229 name="--sampling-rate",
230 minimum=0.01,
231 maximum=1_000_000.0,
232 )
233 epoch_start = finite_float(args.epoch_start, name="--epoch-start")
234 epoch_end = finite_float(args.epoch_end, name="--epoch-end")
235 baseline_start = (
236 finite_float(args.baseline_start, name="--baseline-start")
237 if args.baseline_start is not None
238 else None
239 )
240 baseline_end = (
241 finite_float(args.baseline_end, name="--baseline-end")
242 if args.baseline_end is not None
243 else None
244 )
245 if not 1 <= args.max_events <= MAX_EVENTS:
246 raise CliError(f"--max-events must be between 1 and {MAX_EVENTS}")
248 if args.events_csv:
249 require_deidentified(args.deidentified)
250 path = checked_input_file(
251 args.events_csv,
252 root=args.root,
253 suffixes={".csv"},
254 max_bytes=MAX_CSV_BYTES,
255 )
256 data, _, row_count = read_numeric_columns(
257 path,
258 [args.onset_column],
259 max_rows=args.max_events,
260 )
261 if row_count > args.max_events:
262 raise CliError(f"event count exceeds {args.max_events}")
263 values = [float(value) for value in data[args.onset_column]]
264 else:
265 values = _parse_events(args.events)
266 if len(values) > args.max_events:
267 raise CliError(f"event count exceeds {args.max_events}")
269 event_samples = _event_samples(
270 values,
271 unit=args.event_unit,
272 sampling_rate=sampling_rate,
273 )
274 report = plan(
275 event_samples,
276 recording_samples=args.recording_samples,
277 sampling_rate=sampling_rate,
278 epoch_start_s=epoch_start,
279 epoch_end_s=epoch_end,
280 baseline_start_s=baseline_start,
281 baseline_end_s=baseline_end,
282 boundary_policy=args.boundary_policy,
283 )
284 report["input_event_unit"] = args.event_unit
285 report["path_redacted"] = True
286 emit_json(report, output=args.output, root=args.root, force=args.force)
289if __name__ == "__main__":
290 raise SystemExit(run_cli(main))