feat(plot): nice-number ticks, axes, grid, traces, cursors, measurement queries

src/plot/axis.cpp

@@ -0,0 +1,18 @@
+#include <plot/plot/axis.hpp>
+
+namespace plot {
+
+void Axis::refresh(f32 px, f32 min_px, bool include_minors) {
+  generate_ticks(range, px, min_px, ticks, include_minors);
+  if (format != TickFormat::Time) return;
+  // generate_ticks labels with format_si; a time axis relabels its majors here,
+  // because the tick generator has no business knowing about a clock.
+  char buf[64];
+  for (Tick& tick : ticks) {
+    if (!tick.major) continue;
+    format_time(tick.value, buf, sizeof buf);
+    tick.label = buf;
+  }
+}
+
+}  // namespace plot

src/plot/grid.cpp

@@ -0,0 +1,45 @@
+#include <plot/plot/grid.hpp>
+
+#include <algorithm>
+#include <cmath>
+
+namespace plot {
+namespace {
+
+/// True when `v` lies in `r`, allowing the rounding slack a generated tick can
+/// carry at the range's edge: a tick that is within 1e-9 relative of the edge is
+/// the edge, and dropping it would leave the frame un-drawn.
+bool inside(f64 v, Range r) noexcept {
+  const f64 tolerance = std::max(std::fabs(r.min), std::fabs(r.max)) * 1e-9 + 1e-12;
+  return v >= r.min - tolerance && v <= r.max + tolerance;
+}
+
+}  // namespace
+
+void Grid::build(const std::vector<Tick>& x_ticks, Range x_range,
+                 const std::vector<Tick>& y_ticks, Range y_range) {
+  lines_.clear();
+  if (!std::isfinite(x_range.min) || !std::isfinite(x_range.max) ||
+      !std::isfinite(y_range.min) || !std::isfinite(y_range.max)) {
+    return;
+  }
+
+  // Segments stay in data space: the renderer applies the transform, so a pan or
+  // zoom needs no grid rebuild other than for a changed tick set.
+  if (vertical) {
+    for (const Tick& tick : x_ticks) {
+      if (!tick.major && !show_minor) continue;
+      if (!std::isfinite(tick.value) || !inside(tick.value, x_range)) continue;
+      lines_.push_back(GridLine{tick.value, y_range.min, tick.value, y_range.max, tick.major});
+    }
+  }
+  if (horizontal) {
+    for (const Tick& tick : y_ticks) {
+      if (!tick.major && !show_minor) continue;
+      if (!std::isfinite(tick.value) || !inside(tick.value, y_range)) continue;
+      lines_.push_back(GridLine{x_range.min, tick.value, x_range.max, tick.value, tick.major});
+    }
+  }
+}
+
+}  // namespace plot

src/plot/measure.cpp

@@ -0,0 +1,43 @@
+#include <plot/plot/measure.hpp>
+
+#include <cmath>
+
+namespace plot {
+
+Measurement Measure::evaluate(const data::Series& s, u64 begin, u64 end) {
+  Measurement m;
+  // Both guards run before any query, so an empty, inverted or out-of-range
+  // interval returns an invalid result without touching the series.
+  if (end <= begin || end > s.size()) return m;
+
+  // Exact extremes and sums come from Series -- never from the pyramid's
+  // *_buckets forms, which answer in whole buckets only and would drop up to two
+  // partial edge buckets, silently returning wrong extrema and sums for any
+  // interval that is not bucket-aligned. Each query here is O(kBucket) worst
+  // case, so a measurement never scales with series length.
+  const Range extremes = s.reduce(begin, end);
+  const f64 total = s.sum(begin, end);
+  const f64 squared = s.sumsq(begin, end);
+  const f64 n = static_cast<f64>(end - begin);
+  const f64 step = s.dx();
+
+  m.ymin = extremes.min;
+  m.ymax = extremes.max;
+  m.mean = total / n;
+  m.rms = std::sqrt(squared / n);
+  m.integral = total * step;  // left-Riemann sum: see the note in measure.hpp
+
+  if (end - begin >= 2) {
+    // First and last sample of the interval: a single-sample query reduces to
+    // that sample's value, and costs the same O(kBucket) as any other.
+    const f64 first = s.reduce(begin, begin + 1).min;
+    const f64 last = s.reduce(end - 1, end).min;
+    m.dy = last - first;
+    m.dx = static_cast<f64>(end - 1 - begin) * step;
+    m.slope = m.dx != 0.0 ? m.dy / m.dx : 0.0;
+  }
+  m.valid = true;
+  return m;
+}
+
+}  // namespace plot

src/plot/plot.cpp

@@ -0,0 +1,371 @@
+#include <plot/plot/plot.hpp>
+
+#include <algorithm>
+#include <cmath>
+
+namespace plot {
+namespace {
+
+/// Layout passes. The left strip's width depends on the tick labels, which depend
+/// on the plot width the strip leaves; the loop stops as soon as the strip stops
+/// changing, so a normal plot converges on the second pass and never reaches this.
+constexpr u32 kMaxLayoutPasses = 4;
+
+const Axis& default_axis() noexcept {
+  static const Axis axis{};
+  return axis;
+}
+
+const Measurement& invalid_measurement() noexcept {
+  static const Measurement measurement{};
+  return measurement;
+}
+
+Axis axis_from(const AxisDesc& desc) {
+  Axis axis;
+  axis.label = desc.label;
+  axis.range = desc.range;
+  axis.auto_scale = desc.auto_scale;
+  axis.place = desc.place;
+  axis.log = desc.log;
+  axis.format = desc.format;
+  return axis;
+}
+
+/// Adds the margin an auto-fitted range needs. A zero-span union (a single-sample
+/// series) is widened instead of padded, because a degenerate range has no span to
+/// take 2% of and would leave the transform's inverse ambiguous.
+Range fit_range(Range r, bool log) {
+  if (!std::isfinite(r.min) || !std::isfinite(r.max)) return Range{0.0, 1.0};
+  f64 lo = std::min(r.min, r.max);
+  f64 hi = std::max(r.min, r.max);
+
+  if (log && lo > 0.0) {
+    if (hi > lo) {
+      const f64 factor = std::pow(hi / lo, 0.02);  // 2% of the span, in log space
+      return Range{lo / factor, hi * factor};
+    }
+    return Range{lo / 10.0, hi * 10.0};            // one decade either side
+  }
+
+  if (hi - lo > 0.0) {
+    const f64 pad = (hi - lo) * 0.02;
+    return Range{lo - pad, hi + pad};
+  }
+  const f64 pad = std::fabs(lo) * 1e-3 + 0.5;
+  return Range{lo - pad, hi + pad};
+}
+
+}  // namespace
+
+u32 Plot::add_x_axis(const AxisDesc& desc) {
+  x_axes_.push_back(axis_from(desc));
+  return static_cast<u32>(x_axes_.size() - 1);
+}
+
+u32 Plot::add_y_axis(const AxisDesc& desc) {
+  y_axes_.push_back(axis_from(desc));
+  return static_cast<u32>(y_axes_.size() - 1);
+}
+
+u32 Plot::add_trace(const TraceDesc& desc) {
+  Trace trace;
+  trace.series = desc.series;
+  trace.color = desc.color;
+  trace.width = desc.width;
+  trace.axis_x = desc.axis_x;
+  trace.axis_y = desc.axis_y;
+  traces_.push_back(trace);
+  mark_measures_dirty();
+  return static_cast<u32>(traces_.size() - 1);
+}
+
+u32 Plot::add_cursor(f64 x, Color color) {
+  Cursor cursor;
+  cursor.x = x;
+  cursor.color = color;
+  cursor.x = snapped_x(cursor);
+  cursors_.push_back(cursor);
+  mark_measures_dirty();
+  return static_cast<u32>(cursors_.size() - 1);
+}
+
+u32 Plot::add_measure(u32 cursor_a, u32 cursor_b, u32 trace) {
+  MeasureSlot slot;
+  slot.cursor_a = cursor_a;
+  slot.cursor_b = cursor_b;
+  slot.trace = trace;
+  measures_.push_back(slot);
+  return static_cast<u32>(measures_.size() - 1);
+}
+
+u32 Plot::add_annotation(const Annotation& annotation) {
+  annotations_.push_back(annotation);
+  return static_cast<u32>(annotations_.size() - 1);
+}
+
+const Axis& Plot::x_axis(u32 axis) const {
+  return axis < x_axes_.size() ? x_axes_[axis] : default_axis();
+}
+
+const Axis& Plot::y_axis(u32 axis) const {
+  return axis < y_axes_.size() ? y_axes_[axis] : default_axis();
+}
+
+const Trace& Plot::trace(u32 id) const {
+  static const Trace fallback{};
+  return id < traces_.size() ? traces_[id] : fallback;
+}
+
+const Cursor& Plot::cursor(u32 id) const {
+  static const Cursor fallback{};
+  return id < cursors_.size() ? cursors_[id] : fallback;
+}
+
+const Annotation& Plot::annotation(u32 id) const {
+  static const Annotation fallback{};
+  return id < annotations_.size() ? annotations_[id] : fallback;
+}
+
+Range Plot::x_range(u32 axis) const { return x_axis(axis).range; }
+
+Range Plot::y_range(u32 axis) const { return y_axis(axis).range; }
+
+void Plot::set_x_range(u32 axis, Range r) {
+  if (axis >= x_axes_.size()) return;
+  x_axes_[axis].range = r;
+  x_axes_[axis].auto_scale = false;
+}
+
+void Plot::set_y_range(u32 axis, Range r) {
+  if (axis >= y_axes_.size()) return;
+  y_axes_[axis].range = r;
+  y_axes_[axis].auto_scale = false;
+}
+
+void Plot::set_x_auto_scale(u32 axis, bool on) {
+  if (axis < x_axes_.size()) x_axes_[axis].auto_scale = on;
+}
+
+void Plot::set_y_auto_scale(u32 axis, bool on) {
+  if (axis < y_axes_.size()) y_axes_[axis].auto_scale = on;
+}
+
+Transform Plot::x_transform(u32 axis) const {
+  const Axis& a = x_axis(axis);
+  const Rect r = layout_.axes;
+  // Only the x half of the transform is meaningful; the dummy y range keeps the
+  // returned Transform's y mapping an identity rather than a second, conflicting
+  // axis mapping. A collapsed rect or range is widened inside from_range.
+  return Transform::from_range(a.range, Range{}, static_cast<f64>(r.x),
+                               static_cast<f64>(r.right()), 0.0, 1.0, a.log, false);
+}
+
+Transform Plot::y_transform(u32 axis) const {
+  const Axis& a = y_axis(axis);
+  const Rect r = layout_.axes;
+  // Screen y grows downwards, so the range's minimum maps to the rect's bottom
+  // and its maximum to the top: py0/py1 are passed bottom-first on purpose.
+  return Transform::from_range(Range{}, a.range, 0.0, 1.0, static_cast<f64>(r.bottom()),
+                               static_cast<f64>(r.y), false, a.log);
+}
+
+bool Plot::fit_x_axis(u32 axis) {
+  Range union_range{};
+  bool any = false;
+  for (const Trace& t : traces_) {
+    if (t.axis_x != axis || !t.visible()) continue;  // null or empty series: skip
+    const Range bounds = t.series->x_bounds();
+    if (!std::isfinite(bounds.min) || !std::isfinite(bounds.max)) continue;
+    if (!any) {
+      union_range = bounds;
+      any = true;
+    } else {
+      union_range.min = std::min(union_range.min, bounds.min);
+      union_range.max = std::max(union_range.max, bounds.max);
+    }
+  }
+  if (!any) return false;  // nothing bound to this axis: leave the caller's range
+
+  Axis& a = x_axes_[axis];
+  Range fitted = fit_range(union_range, a.log);
+  if (a.log && fitted.min <= 0.0) fitted.min = fitted.max > 0.0 ? fitted.max * 1e-3 : 1.0;
+  a.range = fitted;
+  return true;
+}
+
+bool Plot::fit_y_axis(u32 axis) {
+  Range union_range{};
+  bool any = false;
+  for (const Trace& t : traces_) {
+    if (t.axis_y != axis || !t.visible()) continue;
+    const Range bounds = t.series->y_bounds();
+    if (!std::isfinite(bounds.min) || !std::isfinite(bounds.max)) continue;
+    if (!any) {
+      union_range = bounds;
+      any = true;
+    } else {
+      union_range.min = std::min(union_range.min, bounds.min);
+      union_range.max = std::max(union_range.max, bounds.max);
+    }
+  }
+  if (!any) return false;
+
+  Axis& a = y_axes_[axis];
+  Range fitted = fit_range(union_range, a.log);
+  if (a.log && fitted.min <= 0.0) fitted.min = fitted.max > 0.0 ? fitted.max * 1e-3 : 1.0;
+  a.range = fitted;
+  return true;
+}
+
+void Plot::auto_fit_x() {
+  for (u32 i = 0; i < x_axes_.size(); ++i) {
+    // Fitting re-arms auto scaling, so the next layout keeps the axis framed. An
+    // axis with nothing bound to it has nothing to fit and nothing to re-arm.
+    if (fit_x_axis(i)) x_axes_[i].auto_scale = true;
+  }
+}
+
+void Plot::auto_fit_y(u32 axis) {
+  if (axis >= y_axes_.size()) return;
+  if (fit_y_axis(axis)) y_axes_[axis].auto_scale = true;
+}
+
+void Plot::autoscale_axes() {
+  for (u32 i = 0; i < x_axes_.size(); ++i) {
+    if (x_axes_[i].auto_scale) fit_x_axis(i);
+  }
+  for (u32 i = 0; i < y_axes_.size(); ++i) {
+    if (y_axes_[i].auto_scale) fit_y_axis(i);
+  }
+}
+
+void Plot::set_cursor_x(u32 id, f64 x) {
+  if (id >= cursors_.size()) return;
+  cursors_[id].x = x;
+  cursors_[id].x = snapped_x(cursors_[id]);
+  mark_measures_dirty();
+}
+
+void Plot::set_cursor_snap(u32 id, SnapMode mode) {
+  if (id >= cursors_.size()) return;
+  cursors_[id].snap = mode;
+  cursors_[id].x = snapped_x(cursors_[id]);
+  mark_measures_dirty();
+}
+
+const Measurement& Plot::measurement(u32 id) const {
+  if (id >= measures_.size()) return invalid_measurement();
+  const MeasureSlot& slot = measures_[id];
+  const data::Series* series = nullptr;
+  if (slot.trace < traces_.size()) series = traces_[slot.trace].series;
+  const u64 samples = series != nullptr ? series->size() : 0;
+
+  // Recompute only when someone moved a cursor, replaced the trace's series, or
+  // grew it. Everything else is a cache hit, which is what keeps a steady frame
+  // from re-evaluating every measure.
+  if (slot.dirty || series != slot.last_series || samples != slot.last_size) {
+    slot.result = Measurement{};
+    if (series != nullptr && samples > 0 && slot.cursor_a < cursors_.size() &&
+        slot.cursor_b < cursors_.size()) {
+      const u64 a = series->pyramid().snap(cursors_[slot.cursor_a].x, series->x0(), series->dx());
+      const u64 b = series->pyramid().snap(cursors_[slot.cursor_b].x, series->x0(), series->dx());
+      // snap() clamps to [0, size-1], so [a, b+1) is inside the series whenever
+      // the cursor pair is not inverted. An inverted pair is not an interval.
+      if (a <= b) slot.result = Measure::evaluate(*series, a, b + 1);
+    }
+    slot.last_series = series;
+    slot.last_size = samples;
+    slot.dirty = false;
+  }
+  return slot.result;
+}
+
+void Plot::invalidate_measurements() noexcept { mark_measures_dirty(); }
+
+const PlotLayout& Plot::layout() {
+  // Axes still on auto scale refit here, once per layout, so a growing stream
+  // stays framed. A caller gesture clears auto_scale (set_x_range/set_y_range),
+  // which is precisely what stops this from fighting the range the user chose.
+  autoscale_axes();
+
+  const f32 pad = style_.grid_label_pad;
+  // Seeds: PlotLayout's defaults when nobody can measure text, otherwise the
+  // metrics the text backend does know (its line height, and just the padding
+  // until the first tick list says how wide the labels are).
+  f32 y_width = PlotLayout{}.y_axis_width;
+  f32 x_height = PlotLayout{}.x_axis_height;
+  if (text_ != nullptr) {
+    y_width = y_axes_.empty() ? 0.0f : 2.0f * pad;
+    x_height = text_->height() + 2.0f * pad;
+  }
+
+  for (u32 pass = 0;; ++pass) {
+    // Clamped to the viewport so a huge label set can never push the axes rect
+    // outside it (or make its width negative).
+    y_width = std::clamp(y_width, 0.0f, viewport_.w);
+    x_height = std::clamp(x_height, 0.0f, viewport_.h);
+    layout_.y_axis_width = y_width;
+    layout_.x_axis_height = x_height;
+    layout_.axes = Rect{viewport_.x + y_width, viewport_.y, viewport_.w - y_width,
+                        viewport_.h - x_height};
+    layout_.y_label_area = Rect{viewport_.x, layout_.axes.y, y_width, layout_.axes.h};
+    layout_.x_label_area = Rect{layout_.axes.x, layout_.axes.bottom(), layout_.axes.w, x_height};
+
+    refresh_ticks();
+
+    if (text_ == nullptr || pass + 1 >= kMaxLayoutPasses) break;
+    const f32 need_y = std::clamp(y_label_width(), 0.0f, viewport_.w);
+    const f32 need_x = std::clamp(text_->height() + 2.0f * pad, 0.0f, viewport_.h);
+    if (need_y == y_width && need_x == x_height) break;  // strips are a fixed point
+    y_width = need_y;
+    x_height = need_x;
+  }
+
+  build_grid();
+  return layout_;
+}
+
+void Plot::refresh_ticks() {
+  // Ticks are generated for the final axes size: an x axis' density depends on
+  // the plot's width and a y axis' on its height.
+  for (Axis& a : x_axes_) a.refresh(layout_.axes.w, style_.tick_min_px, true);
+  for (Axis& a : y_axes_) a.refresh(layout_.axes.h, style_.tick_min_px, true);
+}
+
+void Plot::build_grid() {
+  if (x_axes_.empty() || y_axes_.empty()) {
+    grid_.clear();
+    return;
+  }
+  grid_.build(x_axes_[0].ticks, x_axes_[0].range, y_axes_[0].ticks, y_axes_[0].range);
+}
+
+void Plot::mark_measures_dirty() noexcept {
+  for (MeasureSlot& slot : measures_) slot.dirty = true;
+}
+
+f64 Plot::snapped_x(const Cursor& c) const {
+  if (c.snap != SnapMode::NearestSample) return c.x;
+  for (const Trace& t : traces_) {
+    if (t.axis_x != c.axis_x || !t.visible()) continue;
+    const data::Series& s = *t.series;
+    const u64 i = s.pyramid().snap(c.x, s.x0(), s.dx());
+    return s.x0() + static_cast<f64>(i) * s.dx();
+  }
+  return c.x;
+}
+
+f32 Plot::y_label_width() const {
+  if (text_ == nullptr || y_axes_.empty()) return 0.0f;
+  f32 widest = 0.0f;
+  for (const Axis& a : y_axes_) {
+    for (const Tick& tick : a.ticks) {
+      if (!tick.major || tick.label.empty()) continue;
+      widest = std::max(widest, text_->measure(tick.label.c_str()));
+    }
+  }
+  return widest + 2.0f * style_.grid_label_pad;
+}
+
+}  // namespace plot

src/plot/ticks.cpp

@@ -0,0 +1,181 @@
+#include <plot/plot/ticks.hpp>
+
+#include <cmath>
+#include <cstdio>
+#include <cstring>
+
+namespace plot {
+namespace {
+
+/// Longest label a tick keeps before falling back to `%.3g`.
+constexpr usize kMaxLabelChars = 12;
+/// Scratch for one label; longest possible output is a `%.3g` fallback plus a
+/// prefix plus a unit, far under this.
+constexpr usize kLabelBuf = 64;
+
+/// SI prefixes indexed by engineering exponent / 3, from 10^-24 to 10^24. The
+/// micro sign is UTF-8 U+00B5 (0xC2 0xB5), not the Greek mu U+03BC: the bitmap
+/// font has a glyph for U+00B5 and would draw a missing-glyph box for U+03BC.
+struct SiPrefix {
+  int index;
+  const char* text;
+};
+constexpr SiPrefix kPrefixes[] = {
+    {-8, "y"}, {-7, "z"}, {-6, "a"},  {-5, "f"}, {-4, "p"},  {-3, "n"}, {-2, "\xC2\xB5"},
+    {-1, "m"}, {0, ""},   {1, "k"},   {2, "M"},  {3, "G"},   {4, "T"},  {5, "P"},
+    {6, "E"},  {7, "Z"},  {8, "Y"},
+};
+
+const char* prefix_for(int index) noexcept {
+  for (const SiPrefix& p : kPrefixes) {
+    if (p.index == index) return p.text;
+  }
+  return "";
+}
+
+/// `%.3g` with a NUL-terminated, never-overflowing write. Used for every value
+/// the SI path refuses.
+void fallback_number(f64 v, const char* unit, char* buf, usize n) noexcept {
+  std::snprintf(buf, n, "%.3g%s", v, unit);
+}
+
+}  // namespace
+
+void format_si(f64 v, char* buf, usize n, const char* unit) {
+  if (buf == nullptr || n == 0) return;
+  buf[0] = '\0';
+  const char* const u = unit != nullptr ? unit : "";
+
+  if (!std::isfinite(v)) {
+    fallback_number(v, u, buf, n);
+    return;
+  }
+  if (v == 0.0) {
+    std::snprintf(buf, n, "0%s", u);
+    return;
+  }
+
+  // Engineering exponent: |v| == scaled * 10^(3*e3) with scaled in [1, 1000).
+  const f64 magnitude = std::fabs(v);
+  const int e3 = static_cast<int>(std::floor(std::log10(magnitude) / 3.0));
+  if (e3 < kPrefixes[0].index || e3 > kPrefixes[sizeof(kPrefixes) / sizeof(kPrefixes[0]) - 1].index) {
+    fallback_number(v, u, buf, n);
+    return;
+  }
+
+  f64 scaled = v / std::pow(10.0, 3.0 * static_cast<f64>(e3));
+  int index = e3;
+  // Rounding to three significant digits must not spill the scale: 999.9 goes
+  // to 1.00 of the next prefix rather than to `1e+03` of this one.
+  if (std::fabs(scaled) >= 999.5) {
+    index = e3 + 1;
+    if (index > kPrefixes[sizeof(kPrefixes) / sizeof(kPrefixes[0]) - 1].index) {
+      fallback_number(v, u, buf, n);
+      return;
+    }
+    scaled = v / std::pow(10.0, 3.0 * static_cast<f64>(index));
+  }
+
+  std::snprintf(buf, n, "%.3g%s%s", scaled, prefix_for(index), u);
+  if (std::strlen(buf) > kMaxLabelChars) {
+    // A 12-character budget keeps axis labels readable; anything longer is a
+    // scientific-notation value the caller is better served by verbatim.
+    fallback_number(v, u, buf, n);
+  }
+}
+
+void format_time(f64 seconds, char* buf, usize n) {
+  if (buf == nullptr || n == 0) return;
+  buf[0] = '\0';
+  if (!std::isfinite(seconds)) {
+    std::snprintf(buf, n, "%.3g", seconds);
+    return;
+  }
+
+  const char* const sign = seconds < 0.0 ? "-" : "";
+  const f64 t = std::fabs(seconds);
+  if (t >= 86400.0) {
+    std::snprintf(buf, n, "%s%.3gd", sign, t / 86400.0);
+  } else if (t >= 3600.0) {
+    const auto hours = static_cast<unsigned long long>(t / 3600.0);
+    const auto mins = static_cast<unsigned long long>((t - 3600.0 * static_cast<f64>(hours)) / 60.0);
+    const auto secs =
+        static_cast<unsigned long long>(std::llround(t - 3600.0 * static_cast<f64>(hours) -
+                                                     60.0 * static_cast<f64>(mins)));
+    std::snprintf(buf, n, "%s%llu:%02llu:%02llu", sign, hours, mins, secs);
+  } else if (t >= 60.0) {
+    const auto mins = static_cast<unsigned long long>(t / 60.0);
+    const auto secs = static_cast<unsigned long long>(std::llround(t - 60.0 * static_cast<f64>(mins)));
+    std::snprintf(buf, n, "%s%llu:%02llu", sign, mins, secs);
+  } else if (t >= 1.0) {
+    std::snprintf(buf, n, "%s%.3gs", sign, t);
+  } else {
+    std::snprintf(buf, n, "%s%.3gms", sign, t * 1000.0);
+  }
+}
+
+void generate_ticks(Range range, f32 px, f32 min_px, std::vector<Tick>& out, bool include_minors) {
+  out.clear();
+
+  const f64 lo = range.min;
+  const f64 hi = range.max;
+  const f64 span = hi - lo;
+  // Non-finite endpoints, an inverted range and a zero span (a collapsed range,
+  // e.g. the x bounds of a single-sample series) all have nothing to tick.
+  if (!std::isfinite(lo) || !std::isfinite(hi) || !(span > 0.0)) return;
+  if (!(px > 0.0f) || !(min_px > 0.0f)) return;
+
+  const f64 target = std::max(1.0, static_cast<f64>(px) / static_cast<f64>(min_px));
+  const f64 raw_step = span / target;
+  if (!std::isfinite(raw_step) || !(raw_step > 0.0)) return;
+
+  const f64 magnitude = std::pow(10.0, std::floor(std::log10(raw_step)));
+  if (!std::isfinite(magnitude) || !(magnitude > 0.0)) return;
+  const f64 normalised = raw_step / magnitude;  // in [1, 10)
+  // Round the step UP to the next 1/2/5 x 10^n, so the major spacing really does
+  // meet min_px instead of landing just under it.
+  const f64 step = magnitude * (normalised <= 1.0 ? 1.0
+                                  : normalised <= 2.0 ? 2.0
+                                  : normalised <= 5.0 ? 5.0
+                                                      : 10.0);
+  if (!std::isfinite(step) || !(step > 0.0)) return;  // span so large the step overflowed
+
+  const f64 first_major = std::ceil(lo / step) * step;
+  if (!std::isfinite(first_major)) return;
+
+  // Minor ticks subdivide the major step: by 4 when the step's leading digit is
+  // 2, by 5 for 1 and 5 (and for a step of 10, i.e. a leading 1 one decade up).
+  const f64 lead = step / magnitude;
+  const f64 minor_ratio = (lead > 1.5 && lead < 3.0) ? 4.0 : 5.0;
+  const f64 minor_step = step / minor_ratio;
+  if (!std::isfinite(minor_step) || !(minor_step > 0.0)) return;
+
+  const f64 stride = include_minors ? minor_step : step;
+  const f64 first = include_minors ? std::ceil(lo / minor_step) * minor_step : first_major;
+  if (!std::isfinite(first)) return;
+
+  char label[kLabelBuf];
+  for (f64 v = first; v <= hi + stride * 1e-9;) {
+    if (out.size() >= kMaxTicks) break;
+    const bool major = !include_minors ||
+                       std::fabs(v - std::round(v / step) * step) <= step * 1e-9;
+    Tick tick;
+    tick.value = v;
+    tick.major = major;
+    if (major) {
+      // format_si owns both label guards: it never writes past the buffer, and it
+      // falls back to %.3g for a non-finite value or a label over 12 characters.
+      format_si(v, label, sizeof label);
+      tick.label = label;
+    }
+    out.push_back(std::move(tick));
+
+    const f64 next = v + stride;
+    // No progress would mean an infinite loop on a denormal stride far from the
+    // origin; stop instead of spinning.
+    if (!(next > v)) break;
+    v = next;
+  }
+}
+
+}  // namespace plot