Stampe (1993) gaze-mapping polynomial

Per-eye Stampe (1993) gaze-mapping polynomial in Python, generalised to arbitrary degree.

This package fits a polynomial without cross terms (1, x, y, x², y², ..., xᵈ, yᵈ) plus an optional per-quadrant corner correction that supplies the missing xy cross-term coefficient on outer points.

Corner correction is independent of d: pass exactly 4 outer points (one per quadrant relative to the centroid) and the model adds one xy coefficient per quadrant. With no outer points, the polynomial alone is used.

The staged inner/outer fit matches what the EyeLink 1000 Plus does at calibration and has been verified against EyeLink HREF P-CR data.

Reference: Stampe, D. M. (1993). "Heuristic filtering and reliable calibration methods for video-based pupil-tracking systems." Behavior Research Methods, 25(2), 137-142.

Stampe's published mapping function (paper eq., p. 141, HV9 case):

x₁ = a + b·x + c·y + d·x² + e·y²
y₁ = f + g·x + h·y + i·x² + j·y²
X  = x₁ + m[q]·x₁·y₁
Y  = y₁ + n[q]·x₁·y₁

where (x, y) is the tracker feature, (X, Y) is the screen coordinate, and q is the quadrant into which (x₁, y₁) falls.

Installation

From PyPI:

pip install stampe1993-gaze-mapping

With uv:

uv add stampe1993-gaze-mapping

For a local checkout, install editable:

pip install -e .
# or
uv add --editable .

Quick start

import numpy as np
from stampe1993_gaze_mapping import StampeModel, angular_error

# HV9: 5 inner + 4 outer corner targets
inner_pcr = np.array([...])  # (5, 2) P-CR features at centre + 4 edges
inner_targets = np.array([...])  # (5, 2) screen-pixel targets
outer_pcr = np.array([...])  # (4, 2) P-CR features at the 4 corners
outer_targets = np.array([...])  # (4, 2) screen targets at the 4 corners

model = StampeModel(degree=2)
model.fit(inner_pcr, inner_targets, outer_pcr, outer_targets)

# Predict screen coordinates from a new (or training) P-CR vector
predicted_xy = model.predict(inner_pcr)

# Pixel error → degrees of visual angle (per-axis signed + magnitude)
err_px = predicted_xy - inner_targets
deg_x, deg_y, mag_deg = angular_error(
    err_px,
    pitch_x_mm=screen_width_mm / screen_res_x,
    pitch_y_mm=screen_height_mm / screen_res_y,
    viewing_dist_mm=eye_to_screen_mm,
)
print(f"calibration residual: mean {mag_deg.mean():.3f}° max {mag_deg.max():.3f}°")

For HV5 (no corners) skip the outer arguments — outer_* defaults to None:

model = StampeModel(degree=2)
model.fit(inner5_pcr, inner5_targets)

For HV3 (bilinear):

model = StampeModel(degree=1)
model.fit(inner3_pcr, inner3_targets)

For HV13 (bicubic), either use the model as a pure polynomial fitter on all 13 points, or split 9 inner + 4 outer for Stampe-style correction. The caller decides the split:

# pure polynomial: all 13 as inner, no corners
model = StampeModel(degree=3)
model.fit(all13_pcr, all13_targets)

# Stampe-style: caller picks which 9 are inner and which 4 are outer
model.fit(inner9_pcr, inner9_targets, corners4_pcr, corners4_targets)

Cross terms above xy (x²y, xy², etc.) are not modelled in any configuration. The Stampe-style polynomial drops them and corner correction only puts back xy per quadrant. If you need a full bicubic with all cross terms, use a different fitter.

The fitted model is fully described by these attributes after fit:

• model.degree — polynomial degree.
• model.coef_x — shape (1 + 2*degree,), polynomial coefficients for the X axis.
• model.coef_y — shape (1 + 2*degree,), same for Y.
• model.centroid — shape (2,), the screen-coord centroid of the calibration targets (used as the quadrant origin).
• model.corner — shape (4, 2), the per-quadrant (cx, cy) corner-correction coefficients. Quadrants are indexed (top-left, top-right, bottom-left, bottom-right) in image-y-down coordinates. All zeros if no outer points were supplied.

API

fit raises ValueError on: shape mismatch, fewer than 1 + 2*degree inner points, asymmetric outer (one None, one not), outer_* length ≠ 4, outer points not spanning all 4 quadrants around the centroid, or an outer-predicted point on a cardinal axis through the centroid.

Internals (also exported): design_row(pcr, degree=2) returns the [1, x, y, x², y², ..., xᵈ, yᵈ] feature row for one P-CR vector. quadrant_code(dx, dy) returns 0..3 from the sign of (dx, dy), image-y-down.

Acknowledgments

This project has received funding from the European Union's Horizon Europe research and innovation funding program under grant agreement No 101072410, Eyes4ICU project.