Coverage for scanpath_studio/alignment.py: 97%

380 statements  

« prev     ^ index     » next       coverage.py v7.16.2, created at 2026-10-07 21:10 +0000

1"""Vertical drift-correction (line-assignment) algorithms. 

2 

3Native port of the ten algorithms surveyed by 

4 

5 Carr, J. W., Pescuma, V. N., Furlan, M., Ktori, M., & Crepaldi, D. (2021). 

6 Algorithms for the automated correction of vertical drift in eye-tracking 

7 data. *Behavior Research Methods*, 54, 287–310. 

8 https://doi.org/10.3758/s13428-021-01554-0 

9 

10The reference implementations live in the companion repository 

11``https://github.com/jwcarr/drift`` (``algorithms.py``), released under 

12**CC BY 4.0**. They are adapted/modified here (returning *line assignments* 

13rather than mutating coordinates, plus DataFrame plumbing for this app); see 

14:data:`NOTICE`. We reimplement natively rather than depend on ``eyekit``, which 

15is GPL-3.0 and incompatible with this MIT-licensed, PyPI-distributed project. 

16 

17Each fixation in a trial is assigned the index of the text line it most likely 

18belongs to. The module ships the ten methods surveyed by Carr et al. (2021), a 

19native implementation of the later run-based ``slice`` method, and a consensus 

20vote: 

21 

22 attach, chain, cluster, compare, merge, regress, segment, split, stretch, 

23 warp, slice, consensus 

24 

25Pure functions, no Streamlit. 

26""" 

27 

28from __future__ import annotations 

29 

30import numpy as np 

31import pandas as pd 

32 

33from .measures import cluster_word_lines, word_box_bounds 

34from .model_scanpaths import _ordered_word_rows 

35 

36NOTICE = ( 

37 "Vertical drift-correction algorithms (alignment.py) adapted from " 

38 "Carr, Pescuma, Furlan, Ktori & Crepaldi (2021), 'Algorithms for the " 

39 "automated correction of vertical drift in eye-tracking data', Behavior " 

40 "Research Methods. Reference code: https://github.com/jwcarr/drift " 

41 "(CC BY 4.0). Adapted/modified." 

42) 

43 

44METHOD_CITATIONS = { 

45 "attach": "Carr et al. (2021), nearest-line baseline", 

46 "chain": "Carr et al. (2021), chain", 

47 "cluster": "Carr et al. (2021), cluster", 

48 "compare": "Carr et al. (2021), compare", 

49 "merge": "Carr et al. (2021), merge", 

50 "regress": "Carr et al. (2021), regress", 

51 "segment": "Carr et al. (2021), segment", 

52 "split": "Carr et al. (2021), split", 

53 "stretch": "Carr et al. (2021), stretch", 

54 "warp": "Carr et al. (2021), warp", 

55 "slice": "Schroeder (2022), run-based Slice line assignment", 

56 "consensus": "Carr et al. (2021), wisdom-of-the-crowd vote", 

57} 

58 

59ALGORITHMS = ( 

60 "attach", 

61 "chain", 

62 "cluster", 

63 "compare", 

64 "merge", 

65 "regress", 

66 "segment", 

67 "split", 

68 "stretch", 

69 "warp", 

70 "slice", 

71 "consensus", 

72) 

73 

74# Deterministic seed for the k-means based algorithms (cluster / split) so the 

75# same trial always yields the same assignment across reruns. 

76_KMEANS_SEED = 0 

77 

78 

79# --------------------------------------------------------------------------- 

80# Shared helpers 

81# --------------------------------------------------------------------------- 

82 

83 

84def _dynamic_time_warping( 

85 sequence1: np.ndarray, sequence2: np.ndarray 

86) -> tuple[list[list[int]], np.ndarray]: 

87 """DTW between two point sequences. 

88 

89 Returns ``(path, cost)`` where ``path[i]`` lists the indices of 

90 ``sequence2`` aligned to ``sequence1[i]`` and ``cost`` is the accumulated 

91 cost matrix (total cost is ``cost[-1, -1]``). Ported from ``jwcarr/drift``. 

92 """ 

93 n1 = len(sequence1) 

94 n2 = len(sequence2) 

95 dtw_cost = np.zeros((n1 + 1, n2 + 1)) 

96 dtw_cost[0, :] = np.inf 

97 dtw_cost[:, 0] = np.inf 

98 dtw_cost[0, 0] = 0 

99 for i in range(n1): 

100 for j in range(n2): 

101 this_cost = np.sqrt(np.sum((sequence1[i] - sequence2[j]) ** 2)) 

102 dtw_cost[i + 1, j + 1] = this_cost + min( 

103 dtw_cost[i, j + 1], dtw_cost[i + 1, j], dtw_cost[i, j] 

104 ) 

105 dtw_cost = dtw_cost[1:, 1:] 

106 dtw_path = [[] for _ in range(n1)] 

107 i = n1 - 1 

108 j = n2 - 1 

109 while i > 0 or j > 0: 

110 dtw_path[i].append(j) 

111 if i == 0: 

112 j -= 1 

113 elif j == 0: 

114 i -= 1 

115 else: 

116 neighbours = ( 

117 dtw_cost[i - 1, j - 1], 

118 dtw_cost[i - 1, j], 

119 dtw_cost[i, j - 1], 

120 ) 

121 best = min(neighbours) 

122 if dtw_cost[i - 1, j] == best: 

123 i -= 1 

124 elif dtw_cost[i, j - 1] == best: 

125 j -= 1 

126 else: 

127 i -= 1 

128 j -= 1 

129 dtw_path[0].append(0) 

130 return dtw_path, dtw_cost 

131 

132 

133def _mode(values: np.ndarray) -> float: 

134 """Most common value (lowest wins ties) — used by ``warp``.""" 

135 uniq, counts = np.unique(values, return_counts=True) 

136 return float(uniq[np.argmax(counts)]) 

137 

138 

139def _nearest_line_indices(fixation_Y: np.ndarray, line_Y: np.ndarray) -> np.ndarray: 

140 """Index of the nearest line center for each fixation y.""" 

141 return np.abs(fixation_Y[:, None] - line_Y[None, :]).argmin(axis=1) 

142 

143 

144# --------------------------------------------------------------------------- 

145# The ten algorithms — each returns an int array of line indices (0..m-1), 

146# one per fixation, aligned to the rows of ``fixation_XY``. 

147# --------------------------------------------------------------------------- 

148 

149 

150def _attach(fixation_XY: np.ndarray, line_Y: np.ndarray, word_XY) -> np.ndarray: 

151 return _nearest_line_indices(fixation_XY[:, 1], line_Y) 

152 

153 

154def _chain( 

155 fixation_XY: np.ndarray, 

156 line_Y: np.ndarray, 

157 word_XY, 

158 *, 

159 x_thresh: float = 192, 

160 y_thresh: float = 32, 

161) -> np.ndarray: 

162 n = len(fixation_XY) 

163 assignment = np.zeros(n, dtype=int) 

164 dist_X = np.abs(np.diff(fixation_XY[:, 0])) 

165 dist_Y = np.abs(np.diff(fixation_XY[:, 1])) 

166 end_chain_indices = list(np.where((dist_X > x_thresh) | (dist_Y > y_thresh))[0] + 1) 

167 end_chain_indices.append(n) 

168 start_of_chain = 0 

169 for end_of_chain in end_chain_indices: 

170 mean_y = np.mean(fixation_XY[start_of_chain:end_of_chain, 1]) 

171 line_i = int(np.argmin(np.abs(line_Y - mean_y))) 

172 assignment[start_of_chain:end_of_chain] = line_i 

173 start_of_chain = end_of_chain 

174 return assignment 

175 

176 

177def _cluster(fixation_XY: np.ndarray, line_Y: np.ndarray, word_XY) -> np.ndarray: 

178 # ENG-91: scipy is imported where it is used, never at module import — 

179 # `controls` imports this module on every app run, and Windows Smart App 

180 # Control can block scipy's compiled extensions in the desktop bundle. 

181 from scipy.cluster.vq import kmeans2 

182 

183 m = len(line_Y) 

184 n = len(fixation_XY) 

185 fixation_Y = fixation_XY[:, 1].reshape(-1, 1).astype(float) 

186 try: 

187 centers, labels = kmeans2( 

188 fixation_Y, m, iter=100, minit="++", missing="raise", seed=_KMEANS_SEED 

189 ) 

190 except Exception: 

191 # Degenerate clustering (e.g. an empty cluster) → naive fallback. 

192 return _attach(fixation_XY, line_Y, word_XY) 

193 order = np.argsort(centers[:, 0]) 

194 line_for_cluster = np.empty(m, dtype=int) 

195 line_for_cluster[order] = np.arange(m) 

196 return line_for_cluster[labels.astype(int)][:n] 

197 

198 

199def _segment(fixation_XY: np.ndarray, line_Y: np.ndarray, word_XY) -> np.ndarray: 

200 n = len(fixation_XY) 

201 m = len(line_Y) 

202 assignment = np.zeros(n, dtype=int) 

203 diff_X = np.diff(fixation_XY[:, 0]) 

204 # The m-1 largest return sweeps (most negative Δx) mark the line changes. 

205 line_change_indices = {int(i) for i in np.argsort(diff_X)[: m - 1]} 

206 current_line_i = 0 

207 for fixation_i in range(n): 

208 assignment[fixation_i] = min(current_line_i, m - 1) 

209 if fixation_i in line_change_indices: 

210 current_line_i += 1 

211 return assignment 

212 

213 

214def _split(fixation_XY: np.ndarray, line_Y: np.ndarray, word_XY) -> np.ndarray: 

215 from scipy.cluster.vq import kmeans2 

216 

217 n = len(fixation_XY) 

218 assignment = np.zeros(n, dtype=int) 

219 diff_X = np.diff(fixation_XY[:, 0]).reshape(-1, 1).astype(float) 

220 if len(diff_X) < 2: 

221 return _attach(fixation_XY, line_Y, word_XY) 

222 try: 

223 centers, labels = kmeans2( 

224 diff_X, 2, iter=100, minit="++", missing="raise", seed=_KMEANS_SEED 

225 ) 

226 except Exception: 

227 return _attach(fixation_XY, line_Y, word_XY) 

228 # The cluster with the most negative mean Δx holds the return sweeps. 

229 sweep_marker = int(np.argmin(centers[:, 0])) 

230 end_line_indices = list(np.where(labels == sweep_marker)[0] + 1) 

231 end_line_indices.append(n) 

232 start_of_line = 0 

233 for end_of_line in end_line_indices: 

234 mean_y = np.mean(fixation_XY[start_of_line:end_of_line, 1]) 

235 line_i = int(np.argmin(np.abs(line_Y - mean_y))) 

236 assignment[start_of_line:end_of_line] = line_i 

237 start_of_line = end_of_line 

238 return assignment 

239 

240 

241def _merge( 

242 fixation_XY: np.ndarray, 

243 line_Y: np.ndarray, 

244 word_XY, 

245 *, 

246 y_thresh: float = 32, 

247 gradient_thresh: float = 0.1, 

248 error_thresh: float = 20, 

249) -> np.ndarray: 

250 n = len(fixation_XY) 

251 m = len(line_Y) 

252 diff_X = np.diff(fixation_XY[:, 0]) 

253 dist_Y = np.abs(np.diff(fixation_XY[:, 1])) 

254 boundaries = list(np.where((diff_X < 0) | (dist_Y > y_thresh))[0] + 1) 

255 starts = [0] + boundaries 

256 ends = boundaries + [n] 

257 sequences = [list(range(s, e)) for s, e in zip(starts, ends) if e > s] 

258 # Phases progressively relax the minimum sequence lengths, ending with an 

259 # unconstrained pass that always merges down to m lines (paper defaults). 

260 phases = [ 

261 (3, 3, False), 

262 (1, 3, False), 

263 (1, 1, False), 

264 (1, 1, True), 

265 ] 

266 for min_i, min_j, no_constraints in phases: 

267 while len(sequences) > m: 

268 best_merger = None 

269 best_error = np.inf 

270 for i in range(len(sequences) - 1): 

271 if len(sequences[i]) < min_i: 

272 continue 

273 for j in range(i + 1, len(sequences)): 

274 if len(sequences[j]) < min_j: 

275 continue 

276 candidate = sequences[i] + sequences[j] 

277 xy = fixation_XY[candidate] 

278 if len(np.unique(xy[:, 0])) < 2: 

279 gradient = 0.0 

280 intercept = float(np.mean(xy[:, 1])) 

281 else: 

282 gradient, intercept = np.polyfit(xy[:, 0], xy[:, 1], 1) 

283 residuals = xy[:, 1] - (gradient * xy[:, 0] + intercept) 

284 error = float(np.sqrt(np.mean(residuals**2))) 

285 if not no_constraints and ( 

286 abs(gradient) > gradient_thresh or error > error_thresh 

287 ): 

288 continue 

289 if error < best_error: 

290 best_error = error 

291 best_merger = (i, j) 

292 if best_merger is None: 

293 break 

294 i, j = best_merger 

295 sequences[i] = sequences[i] + sequences[j] 

296 del sequences[j] 

297 assignment = np.zeros(n, dtype=int) 

298 order = sorted( 

299 range(len(sequences)), key=lambda s: np.mean(fixation_XY[sequences[s], 1]) 

300 ) 

301 for line_i, seq_i in enumerate(order): 

302 assignment[sequences[seq_i]] = min(line_i, m - 1) 

303 return assignment 

304 

305 

306def _regress( 

307 fixation_XY: np.ndarray, 

308 line_Y: np.ndarray, 

309 word_XY, 

310 *, 

311 slope_bounds: tuple[float, float] = (-0.1, 0.1), 

312 offset_bounds: tuple[float, float] = (-50, 50), 

313 std_bounds: tuple[float, float] = (1, 20), 

314) -> np.ndarray: 

315 from scipy.optimize import minimize 

316 from scipy.stats import norm 

317 

318 n = len(fixation_XY) 

319 m = len(line_Y) 

320 fixation_X = fixation_XY[:, 0] 

321 fixation_Y = fixation_XY[:, 1] 

322 density = np.zeros((n, m)) 

323 

324 def fit_lines(params, return_assignment=False): 

325 k = slope_bounds[0] + (slope_bounds[1] - slope_bounds[0]) * norm.cdf(params[0]) 

326 o = offset_bounds[0] + (offset_bounds[1] - offset_bounds[0]) * norm.cdf( 

327 params[1] 

328 ) 

329 s = std_bounds[0] + (std_bounds[1] - std_bounds[0]) * norm.cdf(params[2]) 

330 predicted = np.array([fixation_X * k + (line_y + o) for line_y in line_Y]) 

331 for line_i in range(m): 

332 density[:, line_i] = norm.logpdf(fixation_Y, predicted[line_i], s) 

333 if return_assignment: 

334 return density.argmax(axis=1) 

335 return -np.sum(density.max(axis=1)) 

336 

337 best_fit = minimize(fit_lines, [0, 0, 0], method="powell") 

338 return np.asarray(fit_lines(best_fit.x, return_assignment=True), dtype=int) 

339 

340 

341def _stretch( 

342 fixation_XY: np.ndarray, 

343 line_Y: np.ndarray, 

344 word_XY, 

345 *, 

346 scale_bounds: tuple[float, float] = (0.9, 1.1), 

347 offset_bounds: tuple[float, float] = (-50, 50), 

348) -> np.ndarray: 

349 from scipy.optimize import minimize 

350 

351 fixation_Y = fixation_XY[:, 1] 

352 

353 def fit_lines(params, return_assignment=False): 

354 candidate_Y = fixation_Y * params[0] + params[1] 

355 nearest = _nearest_line_indices(candidate_Y, line_Y) 

356 if return_assignment: 

357 return nearest 

358 return float(np.sum(np.abs(candidate_Y - line_Y[nearest]))) 

359 

360 best_fit = minimize( 

361 fit_lines, 

362 [1.0, 0.0], 

363 method="nelder-mead", 

364 bounds=[scale_bounds, offset_bounds], 

365 ) 

366 return np.asarray(fit_lines(best_fit.x, return_assignment=True), dtype=int) 

367 

368 

369def _warp(fixation_XY: np.ndarray, line_Y: np.ndarray, word_XY) -> np.ndarray: 

370 if word_XY is None or len(word_XY) == 0: 

371 return _attach(fixation_XY, line_Y, word_XY) 

372 n = len(fixation_XY) 

373 dtw_path, _ = _dynamic_time_warping(fixation_XY, word_XY) 

374 assignment = np.zeros(n, dtype=int) 

375 for fixation_i, mapped_words in enumerate(dtw_path): 

376 candidate_Y = word_XY[mapped_words, 1] 

377 line_y = _mode(candidate_Y) 

378 assignment[fixation_i] = int(np.argmin(np.abs(line_Y - line_y))) 

379 return assignment 

380 

381 

382def _compare( 

383 fixation_XY: np.ndarray, 

384 line_Y: np.ndarray, 

385 word_XY, 

386 *, 

387 x_thresh: float = 512, 

388 n_nearest_lines: int = 3, 

389) -> np.ndarray: 

390 if word_XY is None or len(word_XY) == 0: 

391 return _attach(fixation_XY, line_Y, word_XY) 

392 n = len(fixation_XY) 

393 assignment = np.zeros(n, dtype=int) 

394 diff_X = np.diff(fixation_XY[:, 0]) 

395 end_line_indices = list(np.where(diff_X < -x_thresh)[0] + 1) 

396 end_line_indices.append(n) 

397 start_of_line = 0 

398 for end_of_line in end_line_indices: 

399 gaze_line = fixation_XY[start_of_line:end_of_line] 

400 mean_y = np.mean(gaze_line[:, 1]) 

401 nearest = np.argsort(np.abs(line_Y - mean_y))[:n_nearest_lines] 

402 costs = np.zeros(len(nearest)) 

403 for k, candidate_line_i in enumerate(nearest): 

404 text_line = word_XY[word_XY[:, 1] == line_Y[candidate_line_i]] 

405 if len(text_line) == 0: 

406 costs[k] = np.inf 

407 continue 

408 _, cost = _dynamic_time_warping(gaze_line[:, 0:1], text_line[:, 0:1]) 

409 costs[k] = cost[-1, -1] 

410 line_i = int(nearest[int(np.argmin(costs))]) 

411 assignment[start_of_line:end_of_line] = line_i 

412 start_of_line = end_of_line 

413 return assignment 

414 

415 

416def _slice( 

417 fixation_XY: np.ndarray, 

418 line_Y: np.ndarray, 

419 word_XY, 

420 *, 

421 run_y_factor: float = 0.65, 

422 run_x_factor: float = 0.60, 

423 same_line_factor: float = 0.45, 

424 adjacent_factor: float = 1.50, 

425) -> np.ndarray: 

426 """Assign locally coherent reading runs to lines. 

427 

428 This is an independent implementation of Slice's published algorithmic 

429 idea: split the scanpath at large horizontal/vertical jumps, seed the run 

430 with the broadest horizontal coverage, then grow line labels by comparing 

431 each run with already labelled fixations at nearby x positions. The final 

432 relative labels are aligned to the absolute stimulus line centres. 

433 

434 Keeping a run together is the important distinction from ``attach``: a 

435 locally drifted fixation does not jump to another line merely because it is 

436 a few pixels closer to that line's centre. 

437 """ 

438 n = len(fixation_XY) 

439 if n == 0: 

440 return np.zeros(0, dtype=int) 

441 

442 sorted_lines = np.sort(np.asarray(line_Y, dtype=float)) 

443 line_gaps = np.diff(sorted_lines) 

444 line_gaps = line_gaps[np.isfinite(line_gaps) & (line_gaps > 0)] 

445 if not len(line_gaps): 

446 return np.zeros(n, dtype=int) 

447 line_height = float(np.median(line_gaps)) 

448 

449 if word_XY is not None and len(word_XY): 

450 text_x = np.asarray(word_XY, dtype=float)[:, 0] 

451 else: 

452 text_x = fixation_XY[:, 0] 

453 text_x = text_x[np.isfinite(text_x)] 

454 text_span = float(np.ptp(text_x)) if len(text_x) > 1 else 0.0 

455 # A line-height-derived floor keeps narrow stimuli from fragmenting at 

456 # ordinary within-line saccades. 

457 run_x = max(text_span * run_x_factor, line_height * 3.0) 

458 run_y = line_height * run_y_factor 

459 

460 dx = np.abs(np.diff(fixation_XY[:, 0])) 

461 dy = np.abs(np.diff(fixation_XY[:, 1])) 

462 boundaries = np.flatnonzero((dx >= run_x) | (dy >= run_y)) + 1 

463 runs = [run for run in np.split(np.arange(n), boundaries) if len(run)] 

464 

465 def horizontal_span(run: np.ndarray) -> float: 

466 return float(np.ptp(fixation_XY[run, 0])) if len(run) > 1 else 0.0 

467 

468 seed = max(range(len(runs)), key=lambda i: (horizontal_span(runs[i]), len(runs[i]))) 

469 relative: dict[int, int] = {seed: 0} 

470 

471 def residual(run: np.ndarray, labelled_indices: np.ndarray) -> float: 

472 """Mean y offset from the closest-in-x labelled fixation.""" 

473 labelled = fixation_XY[labelled_indices] 

474 offsets = [] 

475 for point in fixation_XY[run]: 

476 nearest = int(np.argmin(np.abs(labelled[:, 0] - point[0]))) 

477 offsets.append(point[1] - labelled[nearest, 1]) 

478 return float(np.mean(offsets)) 

479 

480 while len(relative) < len(runs): 

481 candidates: list[tuple[float, int, int]] = [] 

482 for run_i, run in enumerate(runs): 

483 if run_i in relative: 

484 continue 

485 for anchor_label in sorted(set(relative.values())): 

486 anchor_runs = [ 

487 runs[i] for i, label in relative.items() if label == anchor_label 

488 ] 

489 anchor_indices = np.concatenate(anchor_runs) 

490 offset = residual(run, anchor_indices) 

491 if abs(offset) < same_line_factor * line_height: 

492 candidates.append((abs(offset), run_i, anchor_label)) 

493 elif ( 

494 same_line_factor * line_height 

495 <= offset 

496 < adjacent_factor * line_height 

497 ): 

498 candidates.append( 

499 (abs(offset - line_height), run_i, anchor_label + 1) 

500 ) 

501 elif ( 

502 -adjacent_factor * line_height 

503 < offset 

504 <= -same_line_factor * line_height 

505 ): 

506 candidates.append( 

507 (abs(offset + line_height), run_i, anchor_label - 1) 

508 ) 

509 if candidates: 

510 _, run_i, label = min(candidates) 

511 relative[run_i] = label 

512 continue 

513 

514 # A skipped line or unusually large drift can leave no adjacent 

515 # candidate. Extend from the closest labelled run by an integer number 

516 # of typical line heights, ensuring deterministic forward progress. 

517 best: tuple[float, int, int] | None = None 

518 for run_i, run in enumerate(runs): 

519 if run_i in relative: 

520 continue 

521 run_mean = float(np.mean(fixation_XY[run, 1])) 

522 for anchor_i, anchor_label in relative.items(): 

523 anchor_mean = float(np.mean(fixation_XY[runs[anchor_i], 1])) 

524 delta = (run_mean - anchor_mean) / line_height 

525 step = int(np.rint(delta)) 

526 step = step if step else (1 if delta > 0 else -1) 

527 score = abs(delta - step) 

528 candidate = (score, run_i, anchor_label + step) 

529 if best is None or candidate < best: 

530 best = candidate 

531 assert best is not None 

532 _, run_i, label = best 

533 relative[run_i] = label 

534 

535 # Slice labels are relative to the seed. Find the vertical shift that best 

536 # aligns those labels with the available absolute stimulus lines while 

537 # preserving their ordering and spacing. 

538 labels = np.array([relative[i] for i in range(len(runs))], dtype=int) 

539 unique_labels = np.unique(labels) 

540 lowest = int(unique_labels.min()) 

541 label_range = int(unique_labels.max() - lowest) 

542 run_means = np.array([np.mean(fixation_XY[run, 1]) for run in runs]) 

543 if label_range < len(sorted_lines): 

544 fits = [] 

545 for shift in range(len(sorted_lines) - label_range): 

546 absolute = labels - lowest + shift 

547 error = float(np.sum((run_means - sorted_lines[absolute]) ** 2)) 

548 fits.append((error, shift)) 

549 shift = min(fits)[1] 

550 run_assignment = labels - lowest + shift 

551 else: 

552 # More inferred gaze lines than stimulus lines: merge excess relative 

553 # labels into the closest physical line instead of discarding runs. 

554 label_means = { 

555 label: float(np.mean(run_means[labels == label])) for label in unique_labels 

556 } 

557 mapping = { 

558 label: int(np.argmin(np.abs(sorted_lines - mean_y))) 

559 for label, mean_y in label_means.items() 

560 } 

561 run_assignment = np.array([mapping[label] for label in labels], dtype=int) 

562 

563 assignment = np.zeros(n, dtype=int) 

564 for run, line_i in zip(runs, run_assignment): 

565 assignment[run] = int(line_i) 

566 return assignment 

567 

568 

569_DISPATCH = { 

570 "attach": _attach, 

571 "chain": _chain, 

572 "cluster": _cluster, 

573 "compare": _compare, 

574 "merge": _merge, 

575 "regress": _regress, 

576 "segment": _segment, 

577 "split": _split, 

578 "stretch": _stretch, 

579 "warp": _warp, 

580 "slice": _slice, 

581} 

582 

583 

584# --------------------------------------------------------------------------- 

585# Public core + DataFrame wrapper 

586# --------------------------------------------------------------------------- 

587 

588 

589def assign_lines( 

590 fixation_XY: np.ndarray, 

591 line_Y: np.ndarray, 

592 *, 

593 word_XY: np.ndarray | None = None, 

594 method: str = "attach", 

595) -> np.ndarray: 

596 """Assign each fixation a line index in ``[0, len(line_Y))``. 

597 

598 ``fixation_XY`` is an ``(n, 2)`` float array, ``line_Y`` a sorted-ascending 

599 ``(m,)`` array of line centers, and ``word_XY`` an ``(k, 2)`` array of word 

600 centers in reading order (required by ``warp`` / ``compare``; ignored 

601 otherwise). Returns an int ``(n,)`` array of line indices. 

602 """ 

603 if method not in _DISPATCH and method != "consensus": 

604 raise ValueError( 

605 f"unknown alignment method {method!r}; choose from {ALGORITHMS}" 

606 ) 

607 fixation_XY = np.asarray(fixation_XY, dtype=float) 

608 line_Y = np.asarray(line_Y, dtype=float) 

609 n = len(fixation_XY) 

610 if n == 0: 

611 return np.zeros(0, dtype=int) 

612 if len(line_Y) < 2: 

613 return np.zeros(n, dtype=int) 

614 if word_XY is not None: 

615 word_XY = np.asarray(word_XY, dtype=float) 

616 if method == "consensus": 

617 votes = np.vstack( 

618 [ 

619 _DISPATCH[name](fixation_XY, line_Y, word_XY) 

620 for name in _DISPATCH 

621 if name not in {"slice"} 

622 ] 

623 ) 

624 assignment = np.apply_along_axis( 

625 lambda column: np.bincount( 

626 column.astype(int), minlength=len(line_Y) 

627 ).argmax(), 

628 0, 

629 votes, 

630 ) 

631 else: 

632 assignment = _DISPATCH[method](fixation_XY, line_Y, word_XY) 

633 return np.clip(np.asarray(assignment, dtype=int), 0, len(line_Y) - 1) 

634 

635 

636def _line_centers(words: pd.DataFrame) -> np.ndarray: 

637 """Sorted-ascending line center y-coordinates from a single-trial words frame.""" 

638 lines = cluster_word_lines(words) 

639 y_center = ( 

640 pd.to_numeric(words["y"], errors="coerce") 

641 + pd.to_numeric(words["height"], errors="coerce").fillna(0) / 2.0 

642 ) 

643 centers = y_center.groupby(lines).mean().sort_values() 

644 return centers.to_numpy(dtype=float) 

645 

646 

647def _word_centers_reading_order(words: pd.DataFrame) -> np.ndarray: 

648 """Word-box centers ``(x_center, y_center)`` in reading order, NaN-geometry dropped.""" 

649 ordered = _ordered_word_rows(words) 

650 # The interest-area centres, from the experiment's own boxes (BUG-83) — the 

651 # same rectangles the fixations are assigned against. 

652 x0, y0, x1, y1 = word_box_bounds(ordered) 

653 finite = np.isfinite(x0) & np.isfinite(y0) & np.isfinite(x1) & np.isfinite(y1) 

654 cx = ((x0 + x1) / 2.0)[finite] 

655 cy = ((y0 + y1) / 2.0)[finite] 

656 return np.column_stack([cx, cy]) if len(cx) else np.empty((0, 2)) 

657 

658 

659def correct( 

660 fixations: pd.DataFrame, 

661 words: pd.DataFrame, 

662 method: str, 

663 *, 

664 snap: bool = True, 

665) -> tuple[pd.DataFrame, pd.Series]: 

666 """Apply a drift-correction algorithm to a single trial's fixations. 

667 

668 Returns ``(corrected_fixations, assigned_line)``: 

669 

670 - ``corrected_fixations`` is a copy of ``fixations`` with ``y`` snapped to 

671 the assigned line center when ``snap=True`` (unchanged otherwise). 

672 - ``assigned_line`` is a float Series (0-based line index, NaN where 

673 unmappable) index-aligned to ``fixations`` — the same shape/semantics as 

674 :func:`measures.assign_fixation_lines`, so the plots' ``color_by_line`` 

675 path consumes it unchanged. 

676 

677 Passthrough (no change, all-NaN line) when fewer than two text lines, empty 

678 words, or empty fixations. Fixations with non-finite coordinates are dropped 

679 from the algorithm input and left unassigned (NaN, y unchanged). 

680 """ 

681 corrected = fixations.copy() 

682 assigned = pd.Series(np.nan, index=fixations.index, dtype="float64") 

683 if fixations.empty or words is None or words.empty: 

684 return corrected, assigned 

685 

686 line_Y = _line_centers(words) 

687 if len(line_Y) < 2: 

688 return corrected, assigned 

689 

690 fx = pd.to_numeric(fixations["x"], errors="coerce") 

691 fy = pd.to_numeric(fixations["y"], errors="coerce") 

692 finite = fx.notna() & fy.notna() 

693 if not finite.any(): 

694 return corrected, assigned 

695 

696 fixation_XY = np.column_stack( 

697 [fx[finite].to_numpy(dtype=float), fy[finite].to_numpy(dtype=float)] 

698 ) 

699 word_XY = None 

700 if method in ("warp", "compare", "slice", "consensus"): 

701 word_XY = _word_centers_reading_order(words) 

702 

703 line_idx = assign_lines(fixation_XY, line_Y, word_XY=word_XY, method=method) 

704 idx = fixations.index[finite.to_numpy()] 

705 assigned.loc[idx] = line_idx.astype(float) 

706 if snap: 

707 # Snap to (float) line centers — cast first so an int64 `y` column 

708 # doesn't raise pandas' lossy-upcast error. 

709 corrected["y"] = corrected["y"].astype(float) 

710 corrected.loc[idx, "y"] = line_Y[line_idx] 

711 corrected["y_original"] = pd.to_numeric(fixations["y"], errors="coerce") 

712 corrected["y_correction"] = ( 

713 pd.to_numeric(corrected["y"], errors="coerce") - corrected["y_original"] 

714 ) 

715 if method == "consensus": 

716 votes = np.vstack( 

717 [ 

718 assign_lines(fixation_XY, line_Y, word_XY=word_XY, method=name) 

719 for name in _DISPATCH 

720 if name != "slice" 

721 ] 

722 ) 

723 agreement = [ 

724 int(np.sum(votes[:, pos] == line_idx[pos])) for pos in range(len(line_idx)) 

725 ] 

726 corrected["alignment_agreement"] = np.nan 

727 corrected.loc[idx, "alignment_agreement"] = agreement 

728 return corrected, assigned 

729 

730 

731def correction_sensitivity( 

732 fixations: pd.DataFrame, 

733 words: pd.DataFrame, 

734 methods: tuple[str, ...] = ("attach", "slice", "consensus"), 

735) -> tuple[pd.DataFrame, pd.DataFrame]: 

736 """Carry several assignments together and report mean correction (PRE-18).""" 

737 combined = fixations.copy() 

738 report = [] 

739 for method in methods: 

740 corrected, assigned = correct(fixations, words, method) 

741 combined[f"y_{method}"] = corrected["y"] 

742 combined[f"y_{method}_correction"] = corrected["y_correction"] 

743 combined[f"line_{method}"] = assigned 

744 if "alignment_agreement" in corrected: 

745 combined[f"agreement_{method}"] = corrected["alignment_agreement"] 

746 report.append( 

747 { 

748 "algorithm": method, 

749 "average_y_correction": float(corrected["y_correction"].abs().mean()), 

750 "max_y_correction": float(corrected["y_correction"].abs().max()), 

751 } 

752 ) 

753 line_cols = [f"line_{method}" for method in methods] 

754 combined["assignment_disagreement"] = combined[line_cols].nunique(axis=1) > 1 

755 return combined, pd.DataFrame(report)