Add harp.data.synchronization for aligning non-Harp devices to the Harp clock

[bruno-f-cruz]

Add harp.data.synchronization for aligning non-Harp devices to the Harp clock

Closes #14. Supersedes harp-tech/harp-python#38, and is the analysis-side counterpart to the spec discussion in harp-tech/protocol#128.

Important

Blocked by harp-tech/protocol#128. The downsampled clock transmission is not specified yet, so this decoder is written against what current emitters happen to do. The open question below (which edge of the packet the whole second belongs to) has to be settled in the spec before this merges — it decides the default and, on today's hardware, a fixed 40 ms of alignment error either way.

Why

Devices outside the Harp bus timestamp their data on their own clock. To make them alignable post-hoc, some emitters (the Behavior board today, White Rabbit's aux UART) mirror the Synchronization Clock on a digital output at a much lower baud rate — typically 1 kbps instead of 100 kbps — so that acquisition systems sampling at 5–30 kHz can record it on a spare digital or analog input. Reading that recording back is currently done with ad-hoc scripts, one per lab.

This PR adds the decoding half of that story to harp.data, in its own synchronization namespace.

What

from harp.data.synchronization import decode_clock_from_samples, decode_clock_from_transitions

decode_clock_from_samples(samples, sample_rate, *, threshold=None,
                          baud_rate=1000.0, anchor="last_bit", max_drift=0.1)
decode_clock_from_transitions(timestamps, states, *,
                              baud_rate=1000.0, anchor="last_bit", max_drift=0.1)

Two entry points, one per shape of recording, each keyed on the axis that system already timestamps its own data on:

input	index	column
uniformly sampled waveform (digital, or analog with a threshold)	Sample (int64)	Time
line transitions (time, 0/1), as event-based systems report	LocalTime (seconds)	Time

Time holds the decoded whole Harp second — the same name the readers give the Harp time axis. The rows are anchors: placing the device's own timestamps on the Harp axis is left to the caller, e.g. np.interp(spike_samples, clock.index, clock["Time"]), since the choice between interpolating neighbouring anchors (absorbs drift) and a global fit (rejects noise) belongs to the analysis, not the decoder.

The wire format decoded is plain RS-232 without parity — idle high, one low start bit, eight data bits LSB-first, one high stop bit — carrying the Harp time as a little-endian uint32. Packets are accepted with or without the 0xAA 0xAF header of the full protocol packet; the framing is detected per packet, so a future emitter that downsamples the complete 6-byte packet needs no flag.

Open question: which edge carries the second (needs resolving)

The 100 kbps protocol is explicit that the last byte of the packet is the synchronization event: its transmission starts 672 µs before the second lapses, and receivers align the next whole second to that. Both existing downsampled emitters do the opposite — per the Behavior board notes in protocol#128 and White Rabbit's aux UART docs, "the first falling edge occurs when a new second has elapsed", i.e. transmission starts on the boundary and the value describes the second that just began.

At 100 kbps the distinction is nearly free: the whole packet is 600 µs, comparable to the 672 µs offset the spec already names. Downsampled it is not — at 1 kbps a 4-byte packet spans 40 ms, so picking the wrong end puts every anchor 40 ms off. It is a constant offset, so recovered drift and interpolation are unaffected and any downstream analysis stays self-consistent; but absolute Harp time is wrong by ~40 ms, which is two orders of magnitude above the sub-millisecond accuracy the clock protocol exists to deliver, and it silently varies with the configured baud rate.

The decoder therefore exposes both, and does not pretend to know which is right: (tho this will change if we managed to standardize in the protocol)

• anchor="last_bit" (current default) — anchors on the end of the last transmitted bit, mirroring the protocol's synchronization event.
• anchor="first_edge" — anchors on the first falling edge, matching what the Behavior board and White Rabbit emit today. This is what reproduces harp-python#38's expected values, and what anyone using existing hardware should pass right now.

What needs deciding in protocol#128: whether the downsampled transmission starts or ends on the second boundary, and whether the residual offset (the equivalent of the 672 µs) is specified or left to the receiver. Once settled, this PR needs at most a default flip and a docstring change — no API change. Related: if the answer is "ends on the boundary", emitters change, not this decoder, so existing recordings will need anchor="first_edge" indefinitely and that should be documented as a legacy mode.

Other decoding decisions worth reviewing

• Per-frame re-synchronization. Bits are read at bit-period centers, and each frame re-syncs on its own start bit like a real UART, rather than predicting all 40 bit positions from the packet start. Without it a ~3% baud error walks off the bits by the fourth byte; with it that error stays a fraction of a bit.
• Two independent rejection layers. Packets whose start/stop bits do not check out are dropped silently — a glitched or half-recorded packet costs one anchor, not the packets around it. max_drift then drops packets whose seconds do not advance consistently with elapsed local time, by growing the longest mutually consistent chain of anchors, so a single corrupt value (including the first) does not take out everything after it. That one warns; max_drift=None disables it.
• Sub-Nyquist input raises rather than silently returning garbage: fewer than two samples per bit is a ValueError naming both rates (five per bit is the practical recommendation).

Tests

35 tests in tests/data/test_clock.py, driven by a synthetic emitter that renders packets to either transitions or a sampled waveform:

• framing: both anchor modes, headered and bare packets, byte patterns whose runs of all-zeros/all-ones leave long stretches without a transition, 5 kHz and 30 kHz sampling, analog input with a threshold, empty and transition-free recordings;
• robustness: idle-line glitches, packets truncated at either end of the recording, ±0.2-bit jitter on every edge, 3% baud error, corrupt packets (first and middle), dropped packets, wrong baud rate, invalid arguments;
• a minute-long 30 kHz recording combining 0.1% clock drift, jitter and one lost packet;
• regression against real data: two Open Ephys traces of a Behavior board emitting at 1 kbps, lifted from harp-python#38, decoding to the Harp times that PR expects (3806874-5 and 2600957-60) — including its corrupt fifth packet being rejected.

Docs

harp.data.synchronization entries on the Data API page, a "Align a non-Harp device to the Harp clock" section in the harp-data README, and a new example page, Aligning Local Timestamps to the Harp Clock, wired into the nav.

ruff format/ruff check, pyright and codespell are clean; the full suite passes (392 tests).

[bruno-f-cruz]

@jsiegle Leaving this as a placeholder and describing what the current blockers are that must be resolved at the specification level before this can be merged.