From 138fa478c3885fe1e54fe49c4a0b7e307a1c8d80 Mon Sep 17 00:00:00 2001 From: Martino Ferrari Date: Sun, 27 Sep 2026 21:29:59 +0200 Subject: [PATCH] feat(data): SPSC ring buffer, incremental min/max LOD pyramid, series types --- CMakeLists.txt | 4 + include/plot/data/array_series.hpp | 43 +++++ include/plot/data/minmax_pyramid.hpp | 82 +++++++++ include/plot/data/ring_buffer.hpp | 101 +++++++++++ include/plot/data/series.hpp | 74 +++++++++ include/plot/data/stream_series.hpp | 77 +++++++++ src/data/array_series.cpp | 68 ++++++++ src/data/minmax_pyramid.cpp | 240 +++++++++++++++++++++++++++ src/data/series.cpp | 153 +++++++++++++++++ src/data/stream_series.cpp | 115 +++++++++++++ tests/CMakeLists.txt | 3 + tests/harness.hpp | 27 ++- tests/test_minmax_pyramid.cpp | 106 ++++++++++++ tests/test_ring_buffer.cpp | 19 +++ tests/test_series.cpp | 78 +++++++++ 15 files changed, 1182 insertions(+), 8 deletions(-) create mode 100644 include/plot/data/array_series.hpp create mode 100644 include/plot/data/minmax_pyramid.hpp create mode 100644 include/plot/data/ring_buffer.hpp create mode 100644 include/plot/data/series.hpp create mode 100644 include/plot/data/stream_series.hpp create mode 100644 src/data/array_series.cpp create mode 100644 src/data/minmax_pyramid.cpp create mode 100644 src/data/series.cpp create mode 100644 src/data/stream_series.cpp create mode 100644 tests/test_minmax_pyramid.cpp create mode 100644 tests/test_ring_buffer.cpp create mode 100644 tests/test_series.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 3d4de1f..e848e45 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -7,6 +7,10 @@ include(plot_options) add_library(plot_core src/core/status.cpp src/math/transform.cpp + src/data/minmax_pyramid.cpp + src/data/series.cpp + src/data/array_series.cpp + src/data/stream_series.cpp src/gl/loader.cpp src/gl/handle.cpp ) diff --git a/include/plot/data/array_series.hpp b/include/plot/data/array_series.hpp new file mode 100644 index 0000000..9019d07 --- /dev/null +++ b/include/plot/data/array_series.hpp @@ -0,0 +1,43 @@ +#pragma once +#include + +#include +#include + +namespace plot::data { + +/// Series over a uniform-dx sample vector. +/// +/// The vector is taken by value and owned from then on, so the caller cannot +/// mutate the samples behind the series' back; the store changes only through +/// reset() (a whole new vector) or commit(), which rebuilds the pyramid and the +/// cached bounds for whatever y_ currently holds. +class ArraySeries final : public Series { + public: + ArraySeries(std::vector y, f64 x0, f64 dx); + /// Replaces samples and grid, then rebuilds derived state (== commit()). + void reset(std::vector y, f64 x0, f64 dx); + /// Rebuilds the pyramid and cached bounds from y_ as it stands now. + void commit(); + + u64 size() const override; + f64 x0() const override; + f64 dx() const override; + Range x_bounds() const override; + Range y_bounds() const override; + const MinMaxPyramid& pyramid() const override; + void gather(u32 level, u64 begin, u64 end, std::vector& out_xy) const override; + + protected: + void copy_raw(u64 begin, u64 end, f32* out) const override; + + private: + std::vector y_; + f64 x0_ = 0.0; + f64 dx_ = 1.0; + MinMaxPyramid pyramid_; + Range x_bounds_{}; + Range y_bounds_{}; +}; + +} // namespace plot::data diff --git a/include/plot/data/minmax_pyramid.hpp b/include/plot/data/minmax_pyramid.hpp new file mode 100644 index 0000000..d632b20 --- /dev/null +++ b/include/plot/data/minmax_pyramid.hpp @@ -0,0 +1,82 @@ +#pragma once +#include +#include + +#include + +namespace plot::data { + +/// Aggregate of one bucket: the extremes of the samples it covers. +struct MinMax { + f32 min = 0.0f; + f32 max = 0.0f; +}; + +/// Samples per level-0 bucket. Fixed so the LOD geometry is identical for every +/// series and the level maths is a shift. +inline constexpr u32 kBucket = 64; + +/// Incremental min/max LOD pyramid. +/// +/// Level 0 buckets kBucket consecutive samples and keeps one (min, max) per +/// bucket; level l buckets 2^l level-0 buckets. Memory is n/4 bytes -- about 6% +/// of 4-byte raw samples. append() merges into the last bucket at each level and +/// walks up the O(log n) chain above it, so building is amortised O(1) per sample +/// and nothing is ever recomputed from raw samples. +/// +/// A pyramid stores *aggregates*, so it can only answer exactly for whole +/// buckets; a range with partial edge buckets would silently pull in samples +/// outside the range. Every query here is therefore bucket-grained and named +/// *_buckets, and exact range queries live on Series, which owns the raw samples. +/// reduce() is the one non-*_buckets query: it is the *covering* aggregate over +/// the level-0 buckets that [begin, end) touches. Use it to ask which buckets a +/// range spans, never to measure a range. +class MinMaxPyramid { + public: + MinMaxPyramid() { add_level(); } + + /// Sizes every level (buckets and prefix sums) for `n` samples so a following + /// run of appends performs no reallocation. The pyramid itself stays empty. + void reserve(u64 n); + void append(f32 v); + void clear(); + + u64 size() const noexcept { return size_; } + /// At least 1, even when empty, so level_count() - 1 is always valid. + u32 level_count() const noexcept { return live_levels_; } + /// Level `l`'s buckets; empty span when `l` is beyond level_count(). + std::span level(u32 l) const; + + /// Min/max over the whole level-0 buckets fully contained in [begin, end). + /// Partial edge buckets are excluded, never approximated; an empty result is a + /// degenerate Range with span() == 0. + Range reduce_buckets(u64 begin, u64 end) const; + /// Sum of the samples in those same whole buckets, accumulated in double. + f64 sum_buckets(u64 begin, u64 end) const; + /// Sum of squares of the samples in those same whole buckets. + f64 sumsq_buckets(u64 begin, u64 end) const; + + /// Min/max over the covering level-0 buckets: every bucket that holds a sample + /// of [begin, end). Exact only when the range is bucket-aligned. + Range reduce(u64 begin, u64 end) const; + + /// Pyramid level whose buckets hold about `samples_per_pixel` samples, for a + /// window of `n` samples: 0 when a pixel covers at most one bucket. Clamped to + /// the levels that exist, to levels with more than one bucket, and to levels + /// that give the window at least two buckets. + u32 level_for(u64 n, f64 samples_per_pixel) const; + /// Index of the sample nearest x on the uniform grid (x0 + i*dx), clamped to + /// [0, size() - 1]. Returns 0 for an empty pyramid or a non-positive dx. + u64 snap(f64 x, f64 x0, f64 dx) const; + + private: + void add_level(); + + u64 size_ = 0; + u64 reserved_n_ = 0; + u32 live_levels_ = 1; + std::vector> levels_; // levels_[0] = kBucket granularity + std::vector> sums_, sumsqs_; // inclusive prefix sums per level +}; + +} // namespace plot::data diff --git a/include/plot/data/ring_buffer.hpp b/include/plot/data/ring_buffer.hpp new file mode 100644 index 0000000..9e0ecf1 --- /dev/null +++ b/include/plot/data/ring_buffer.hpp @@ -0,0 +1,101 @@ +#pragma once +#include +#include +#include +#include +#include + +#include + +namespace plot::data { + +/// Single-producer / single-consumer ring of at most `capacity` elements. +/// +/// The capacity is rounded up to a power of two so slot addressing is a mask +/// instead of a modulo. The producer writes a slot and then publishes the new +/// sample count with a release store; a consumer that observes that count with an +/// acquire load is guaranteed to see the slot contents written before it. +/// +/// A push into a full ring silently discards the oldest element: for a live +/// plotting stream the newest sample is the interesting one and the retained +/// window slides forward. +/// +/// Two indexing conventions, deliberately named apart: +/// - retained-window (operator[], index 0 = oldest element still held) is a +/// single-threaded view: it re-derives the window from the live count, so a +/// producer pushing concurrently would move it under the reader's feet. +/// - absolute (base(), copy_from()) fixes one count for a whole read, which is +/// what a consumer holding a snapshot must use. StreamSeries::acquire() takes +/// such a snapshot and then reads only absolute indices. +template +class RingBuffer { + public: + explicit RingBuffer(u64 capacity) + : cap_(std::bit_ceil(capacity < 1 ? u64{1} : capacity)), + mask_(cap_ - 1), + buf_(static_cast(cap_)) {} + + /// Producer thread only. + void push(T v) { + const u64 i = produced_.load(std::memory_order_relaxed); + buf_[static_cast(i & mask_)] = v; + produced_.store(i + 1, std::memory_order_release); + } + + /// Producer thread only. Equivalent to n calls to push(). + void push(const T* vs, u64 n) { + for (u64 i = 0; i < n; ++i) push(vs[i]); + } + + /// Consumer: retained element count, saturating at capacity(). + u64 size() const noexcept { + const u64 p = produced_.load(std::memory_order_acquire); + return p < cap_ ? p : cap_; + } + + /// Consumer: total number of elements ever pushed (monotonic). + u64 produced() const noexcept { return produced_.load(std::memory_order_acquire); } + + u64 capacity() const noexcept { return cap_; } + + /// Consumer: absolute index of the oldest element still retained. 0 until the + /// ring has wrapped. + u64 base() const noexcept { + const u64 p = produced_.load(std::memory_order_acquire); + return p > cap_ ? p - cap_ : 0; + } + + /// Single-threaded: retained-window index (0 = oldest retained). `i` must be + /// < size(). + const T& operator[](u64 i) const noexcept { + return buf_[static_cast((base() + i) & mask_)]; + } + + /// Consumer: copies the `n` elements at absolute indices [absolute, absolute+n) + /// into `out`, splitting the copy at the ring's wrap point. The caller must hold + /// a snapshot that retains those indices; this never consults the live count, so + /// a pushing producer cannot shift the read. + void copy_from(u64 absolute, u64 n, T* out) const { + if (n == 0) return; + const u64 s = absolute & mask_; + const u64 run = std::min(n, cap_ - s); + const auto first = buf_.begin() + static_cast(s); + std::copy(first, first + static_cast(run), out); + if (run < n) { + std::copy(buf_.begin(), buf_.begin() + static_cast(n - run), + out + static_cast(run)); + } + } + + /// Consumer thread only: discards every element (not safe while a producer is + /// pushing). + void clear() noexcept { produced_.store(0, std::memory_order_relaxed); } + + private: + u64 cap_; + u64 mask_; + std::vector buf_; + std::atomic produced_{0}; +}; + +} // namespace plot::data diff --git a/include/plot/data/series.hpp b/include/plot/data/series.hpp new file mode 100644 index 0000000..43fe9e8 --- /dev/null +++ b/include/plot/data/series.hpp @@ -0,0 +1,74 @@ +#pragma once +#include + +#include +#include + +namespace plot::data { + +/// A uniform-grid series of y samples: sample i lives at x = x0() + i*dx(). +/// +/// Every concrete series answers range queries exactly, and never by scanning the +/// whole range: the bucket-aligned interior comes from the pyramid (O(log n) +/// bucket aggregates) and at most two partial edge buckets are read through +/// copy_raw() and scanned. Cost is independent of series length and no query +/// allocates. +/// +/// Indexing is over the series' own window: index 0 is the oldest sample the +/// series currently holds (for a ring-backed series that is the oldest retained +/// one), so reduce(start, end) and gather() are in the same coordinates as +/// size(), x0() and pyramid(). +class Series { + public: + virtual ~Series() = default; + + virtual u64 size() const = 0; + /// x of sample 0, i.e. of the oldest sample currently held. + virtual f64 x0() const = 0; + /// Uniform x spacing between consecutive samples. + virtual f64 dx() const = 0; + virtual Range x_bounds() const = 0; + virtual Range y_bounds() const = 0; + virtual const MinMaxPyramid& pyramid() const = 0; + + /// Appends interleaved (x, y) pairs in data coordinates for bucket-aligned + /// level `level`; level 0 emits raw samples. begin/end are expected to be + /// bucket-aligned at that level; an unaligned end is widened to the bucket that + /// holds it (so passing end == size() includes the trailing partial bucket). + /// Pairs are appended to out_xy, which is never cleared. + virtual void gather(u32 level, u64 begin, u64 end, std::vector& out_xy) const = 0; + + /// Exact min/max, sum and sum-of-squares over [begin, end). + /// Preconditions: begin < end, end <= size(); a degenerate range yields a + /// zero-span Range and zero sums. Implemented once, non-virtually, over + /// copy_raw() in src/data/series.cpp. + Range reduce(u64 begin, u64 end) const; + f64 sum(u64 begin, u64 end) const; + f64 sumsq(u64 begin, u64 end) const; + + protected: + /// Copies raw y samples for [begin, end) into `out`, which the caller sized to + /// (end - begin). Handles ring wraparound. Only ever called with the small edge + /// spans of the queries above, never on a per-sample path. Out-of-range spans + /// are clamped and the unwritten tail is zero-filled. + virtual void copy_raw(u64 begin, u64 end, f32* out) const = 0; + + /// The shared body of the concrete gather() overrides: raw samples at level 0, + /// bucket min/max envelope above it. + void gather_impl(u32 level, u64 begin, u64 end, std::vector& out_xy) const; + + private: + /// [begin, end) split into the bucket-aligned interior the pyramid can answer + /// exactly, plus the at most two partial edge spans it cannot. + struct QueryParts { + u64 interior_begin = 0, interior_end = 0; // whole level-0 buckets, may be empty + u64 left_begin = 0, left_end = 0; // partial head bucket + u64 right_begin = 0, right_end = 0; // partial tail bucket + }; + static QueryParts split_query(u64 begin, u64 end); + /// Reads both edge spans through copy_raw() into `out` (which must hold at + /// least 2*kBucket samples) and returns how many were written. + u64 copy_edges(const QueryParts& parts, f32* out) const; +}; + +} // namespace plot::data diff --git a/include/plot/data/stream_series.hpp b/include/plot/data/stream_series.hpp new file mode 100644 index 0000000..3922df1 --- /dev/null +++ b/include/plot/data/stream_series.hpp @@ -0,0 +1,77 @@ +#pragma once +#include + +#include +#include +#include +#include + +namespace plot::data { + +/// Single-producer / single-consumer streaming series. +/// +/// The producer thread only ever calls push()/push(ys, n) (and produced()): those +/// touch the ring buffer and nothing else -- no lock, no pyramid, no GL, so the +/// acquisition thread can never stall the renderer. The consumer thread calls +/// acquire() to publish a coherent snapshot, then queries it. +/// +/// Indexing is the *retained window published by the last acquire()*: size(), +/// reduce(), sum(), sumsq(), gather() and pyramid() all use it, and index 0 is +/// the oldest sample still held. Once the producer has written more than +/// capacity() samples the oldest ones are gone, so acquire() notices that the +/// window's base moved and rebuilds the pyramid from the retained samples -- an +/// incremental append would splice a dropped prefix onto a retained suffix and +/// mis-align every bucket (bucket boundaries are absolute, not relative to the +/// window). Until the ring saturates, acquire() is O(new samples); once it wraps, +/// an acquire whose base moved costs O(capacity()). +/// +/// x0() is therefore the x of the oldest retained sample, x_bounds().min == x0(), +/// and y_bounds() covers exactly the published samples. +/// +/// THREADING NOTE. The ring overwrites the oldest samples -- that is what keeps a +/// live stream showing the newest data, and what the drop test pins down -- so a +/// producer racing past capacity() can write the very slot the consumer is reading +/// at the window's oldest edge. That read is not ordered against the write, so a +/// race detector flags it: the value is still a single f32 (no tearing on +/// supported targets) and the worst case is one stale or too-new sample at the +/// oldest edge, corrected by the next acquire() once it sees the base moved. The +/// pyramid itself never mixes windows, because acquire() derives the whole window +/// from one count. Ingest that is race-free even at the edge needs the producer to +/// stop before unread data (drop-on-full) or a claim protocol, i.e. a different +/// contract from the one this class specifies. +class StreamSeries final : public Series { + public: + explicit StreamSeries(u64 capacity, f64 x0, f64 dx); + + /// Producer thread only. + void push(f32 y); + /// Producer thread only. + void push(const f32* ys, u64 n); + /// Producer thread only: number of samples ever pushed (monotonic). + u64 produced() const; + /// Consumer thread only: snapshots the produced count and brings the pyramid + /// and bounds up to it. + void acquire(); + + u64 size() const override; + f64 x0() const override; + f64 dx() const override; + Range x_bounds() const override; + Range y_bounds() const override; + const MinMaxPyramid& pyramid() const override; + void gather(u32 level, u64 begin, u64 end, std::vector& out_xy) const override; + + protected: + void copy_raw(u64 begin, u64 end, f32* out) const override; + + private: + RingBuffer ring_; + f64 x0_ = 0.0; + f64 dx_ = 1.0; + u64 base_ = 0; // absolute index of retained index 0 in the published snapshot + u64 end_ = 0; // absolute index one past the published snapshot + MinMaxPyramid pyramid_; + Range y_bounds_{}; +}; + +} // namespace plot::data diff --git a/src/data/array_series.cpp b/src/data/array_series.cpp new file mode 100644 index 0000000..b21c3c1 --- /dev/null +++ b/src/data/array_series.cpp @@ -0,0 +1,68 @@ +#include + +#include +#include +#include + +namespace plot::data { + +ArraySeries::ArraySeries(std::vector y, f64 x0, f64 dx) { reset(std::move(y), x0, dx); } + +void ArraySeries::reset(std::vector y, f64 x0, f64 dx) { + y_ = std::move(y); + x0_ = x0; + dx_ = dx; + commit(); +} + +void ArraySeries::commit() { + const u64 n = static_cast(y_.size()); + pyramid_.clear(); + pyramid_.reserve(n); + + f32 lo = 0.0f, hi = 0.0f; + bool first = true; + for (const f32 v : y_) { + pyramid_.append(v); + if (first) { + lo = v; + hi = v; + first = false; + } else { + lo = std::min(lo, v); + hi = std::max(hi, v); + } + } + y_bounds_ = first ? Range{0.0, 0.0} : Range{f64(lo), f64(hi)}; + x_bounds_ = n == 0 ? Range{x0_, x0_} : Range{x0_, x0_ + f64(n - 1) * dx_}; +} + +u64 ArraySeries::size() const { return static_cast(y_.size()); } + +f64 ArraySeries::x0() const { return x0_; } + +f64 ArraySeries::dx() const { return dx_; } + +Range ArraySeries::x_bounds() const { return x_bounds_; } + +Range ArraySeries::y_bounds() const { return y_bounds_; } + +const MinMaxPyramid& ArraySeries::pyramid() const { return pyramid_; } + +void ArraySeries::gather(u32 level, u64 begin, u64 end, std::vector& out_xy) const { + gather_impl(level, begin, end, out_xy); +} + +void ArraySeries::copy_raw(u64 begin, u64 end, f32* out) const { + const u64 want = end > begin ? end - begin : 0; + if (want == 0) return; + const u64 n = static_cast(y_.size()); + const u64 b = std::min(begin, n); + const u64 e = std::min(std::max(end, b), n); + const auto first = y_.begin(); + std::copy(first + static_cast(b), first + static_cast(e), out); + std::fill(out + static_cast(e - b), + out + static_cast(want), 0.0f); +} + +} // namespace plot::data diff --git a/src/data/minmax_pyramid.cpp b/src/data/minmax_pyramid.cpp new file mode 100644 index 0000000..20f0e73 --- /dev/null +++ b/src/data/minmax_pyramid.cpp @@ -0,0 +1,240 @@ +#include + +#include +#include + +namespace plot::data { +namespace { + +/// Samples covered by one level-l bucket. +constexpr u64 level_span(u32 l) { return u64{kBucket} << l; } + +/// Level-l buckets needed to hold n samples (ceiling division). +constexpr u64 buckets_for(u64 n, u32 l) { + const u64 span = level_span(l); + return n / span + (n % span == 0 ? u64{0} : u64{1}); +} + +/// First level-0 bucket boundary at or after v. +constexpr u64 ceil_to_bucket(u64 v) { + return (v / kBucket + (v % kBucket == 0 ? u64{0} : u64{1})) * kBucket; +} + +/// Canonical, coarsest-first decomposition of the level-0 bucket span [b0, b1) +/// into pyramid entries, calling fn(level, index) for each. +/// +/// A level-l entry is emitted only when the whole level-l bucket exists (all +/// kBucket << l samples are present, so its aggregate is final). Callers +/// guarantee b1 <= ceil(size() / kBucket), so the only possible remainder is the +/// single trailing partial level-0 bucket, which is emitted at level 0. The +/// decomposition has O(log) entries, so every query below is O(log n). +template +void cover_buckets(const MinMaxPyramid& p, u64 b0, u64 b1, Fn&& fn) { + const u64 complete = p.size() / kBucket; + u64 b = b0; + for (u32 li = p.level_count(); li-- > 0 && b < b1;) { + const u64 lb = u64{1} << li; // level-0 buckets per level-li bucket + const u64 avail = p.level(li).size(); + while (b < b1 && (b % lb) == 0 && b + lb <= b1 && b + lb <= complete && + (b / lb) < avail) { + fn(li, b / lb); + b += lb; + } + } + if (b < b1) fn(0u, b); // trailing partial level-0 bucket +} + +/// Folds the decomposition of [b0, b1) into a min/max Range. Returns a +/// zero-span Range when there is nothing to fold. +Range fold_range(const MinMaxPyramid& p, u64 b0, u64 b1) { + Range r{}; + bool first = true; + cover_buckets(p, b0, b1, [&](u32 l, u64 i) { + const MinMax m = p.level(l)[static_cast(i)]; + if (first) { + r.min = f64(m.min); + r.max = f64(m.max); + first = false; + } else { + r.min = std::min(r.min, f64(m.min)); + r.max = std::max(r.max, f64(m.max)); + } + }); + return first ? Range{0.0, 0.0} : r; +} + +} // namespace + +void MinMaxPyramid::add_level() { + levels_.emplace_back(); + sums_.emplace_back(); + sumsqs_.emplace_back(); +} + +void MinMaxPyramid::reserve(u64 n) { + if (n > reserved_n_) reserved_n_ = n; + const u64 want = std::max(std::max(n, size_), u64{2}); + + // Materialise every level n samples can need -- buckets *and* prefix sums -- + // and size each one up front, so a run of appends allocates nothing at all. + u32 max_levels = 1; + for (u64 span = level_span(0); span <= want && max_levels < 63; span <<= 1) ++max_levels; + while (levels_.size() < max_levels) add_level(); + for (u32 l = 0; l < levels_.size(); ++l) { + const usize want_l = static_cast(buckets_for(want, l)); + levels_[l].reserve(want_l); + sums_[l].reserve(want_l); + sumsqs_[l].reserve(want_l); + } +} + +void MinMaxPyramid::append(f32 v) { + const u64 prev = size_; + const f64 vd = f64(v); + // sum_[l].back() is the sum over every sample appended so far, and the newest + // sample always lands in the last bucket of every level, so the last prefix + // entry of every level is exactly the running total. + const f64 total_before = sums_[0].empty() ? 0.0 : sums_[0].back(); + const f64 total = total_before + vd; + const f64 sq_before = sumsqs_[0].empty() ? 0.0 : sumsqs_[0].back(); + const f64 total_sq = sq_before + vd * vd; + + const u64 idx = prev / kBucket; + if (idx == levels_[0].size()) { + levels_[0].push_back(MinMax{v, v}); + sums_[0].push_back(total); + sumsqs_[0].push_back(total_sq); + } else { + MinMax& cur = levels_[0][static_cast(idx)]; + cur.min = std::min(cur.min, v); + cur.max = std::max(cur.max, v); + sums_[0][static_cast(idx)] = total; + sumsqs_[0][static_cast(idx)] = total_sq; + } + + for (u32 l = 1; l < 63; ++l) { + // A level only carries data once the level below has more than one bucket: a + // one-child parent would just duplicate it. + if (levels_[l - 1].size() < 2) break; + if (l == live_levels_) { + ++live_levels_; + if (l == levels_.size()) { // no reserve() sized this level up front + add_level(); + if (reserved_n_ != 0) { + const usize want = static_cast(buckets_for(reserved_n_, l)); + levels_[l].reserve(want); + sums_[l].reserve(want); + sumsqs_[l].reserve(want); + } + } + } + const u64 ci = prev / level_span(l - 1); // current bucket at level l-1 + const u64 parent = ci / 2; // its parent at level l + MinMax m; + if ((ci & 1u) != 0) { + const MinMax& a = levels_[l - 1][static_cast(ci - 1)]; + const MinMax& b = levels_[l - 1][static_cast(ci)]; + m.min = std::min(a.min, b.min); + m.max = std::max(a.max, b.max); + } else { + m = levels_[l - 1][static_cast(ci)]; + } + if (parent == levels_[l].size()) { + levels_[l].push_back(m); + sums_[l].push_back(total); + sumsqs_[l].push_back(total_sq); + } else { + levels_[l][static_cast(parent)] = m; + sums_[l][static_cast(parent)] = total; + sumsqs_[l][static_cast(parent)] = total_sq; + } + } + size_ = prev + 1; +} + +void MinMaxPyramid::clear() { + size_ = 0; + live_levels_ = 1; + for (auto& lv : levels_) lv.clear(); + for (auto& s : sums_) s.clear(); + for (auto& s : sumsqs_) s.clear(); +} + +std::span MinMaxPyramid::level(u32 l) const { + if (l >= level_count()) return {}; + return levels_[l]; +} + +Range MinMaxPyramid::reduce_buckets(u64 begin, u64 end) const { + const u64 e = std::min(end, size_); + const u64 b = std::min(begin, e); + const u64 b0 = ceil_to_bucket(b) / kBucket; // first whole bucket inside + const u64 b1 = e / kBucket; // one past the last whole bucket + if (b1 <= b0) return Range{0.0, 0.0}; + return fold_range(*this, b0, b1); +} + +f64 MinMaxPyramid::sum_buckets(u64 begin, u64 end) const { + const u64 e = std::min(end, size_); + const u64 b = std::min(begin, e); + const u64 b0 = ceil_to_bucket(b) / kBucket; + const u64 b1 = e / kBucket; + if (b1 <= b0) return 0.0; + f64 total = 0.0; + cover_buckets(*this, b0, b1, [&](u32 l, u64 i) { + const f64 prev = i == 0 ? 0.0 : sums_[l][static_cast(i - 1)]; + total += sums_[l][static_cast(i)] - prev; + }); + return total; +} + +f64 MinMaxPyramid::sumsq_buckets(u64 begin, u64 end) const { + const u64 e = std::min(end, size_); + const u64 b = std::min(begin, e); + const u64 b0 = ceil_to_bucket(b) / kBucket; + const u64 b1 = e / kBucket; + if (b1 <= b0) return 0.0; + f64 total = 0.0; + cover_buckets(*this, b0, b1, [&](u32 l, u64 i) { + const f64 prev = i == 0 ? 0.0 : sumsqs_[l][static_cast(i - 1)]; + total += sumsqs_[l][static_cast(i)] - prev; + }); + return total; +} + +Range MinMaxPyramid::reduce(u64 begin, u64 end) const { + const u64 e = std::min(end, size_); + const u64 b = std::min(begin, e); + if (e <= b) return Range{0.0, 0.0}; + const u64 b0 = b / kBucket; // first bucket holding a sample + const u64 b1 = (e - 1) / kBucket + 1; // one past the last one + return fold_range(*this, b0, b1); +} + +u32 MinMaxPyramid::level_for(u64 n, f64 samples_per_pixel) const { + // The ratio test also rejects NaN, and the clamp keeps a huge (or infinite) + // samples-per-pixel from making the floor-to-u32 conversion ill-formed. + const f64 ratio = samples_per_pixel / f64(kBucket); + if (n == 0 || !(ratio > 1.0)) return 0; + const f64 lg = std::floor(std::log2(ratio)); + u32 l = lg >= f64(62) ? 62u : static_cast(lg); + while (l > 0) { + const bool exists = l < level_count() && level(l).size() >= 2; + if (exists && buckets_for(n, l) >= 2) break; + --l; + } + return l; +} + +u64 MinMaxPyramid::snap(f64 x, f64 x0, f64 dx) const { + if (dx <= 0.0 || size_ == 0) return 0; + const f64 i = (x - x0) / dx; + // Bound the round before it: llround() is undefined for values a long long + // cannot hold, and this also maps NaN to index 0. + if (!(i > 0.0)) return 0; + if (i >= f64(size_)) return size_ - 1; + const u64 idx = static_cast(std::llround(i)); + return idx >= size_ ? size_ - 1 : idx; +} + +} // namespace plot::data diff --git a/src/data/series.cpp b/src/data/series.cpp new file mode 100644 index 0000000..cf4c893 --- /dev/null +++ b/src/data/series.cpp @@ -0,0 +1,153 @@ +#include + +#include +#include + +namespace plot::data { +namespace { + +/// Raw samples read per copy_raw() call while gathering level 0. Bounded so the +/// gather path never allocates, and large enough to amortise the call. +constexpr u64 kGatherChunk = 256; + +} // namespace + +Series::QueryParts Series::split_query(u64 begin, u64 end) { + QueryParts p; + // First level-0 bucket boundary at or after begin, and the last at or before + // end: [interior_begin, interior_end) is exactly the part the pyramid can + // answer with whole buckets. + const u64 up = (begin / kBucket + (begin % kBucket == 0 ? u64{0} : u64{1})) * kBucket; + const u64 down = (end / kBucket) * kBucket; + if (up < down) { + p.interior_begin = up; + p.interior_end = down; + } + // The two edge spans never overlap and never exceed one bucket each: when the + // interior is empty they tile [begin, end) together, otherwise they are the + // head and tail buckets half inside the range. + p.left_begin = begin; + p.left_end = std::min(end, up); + p.right_begin = std::max(down, p.left_end); + p.right_end = end; + return p; +} + +u64 Series::copy_edges(const QueryParts& p, f32* out) const { + u64 n = 0; + if (p.left_begin < p.left_end) { + copy_raw(p.left_begin, p.left_end, out); + n = p.left_end - p.left_begin; + } + if (p.right_begin < p.right_end) { + copy_raw(p.right_begin, p.right_end, out + static_cast(n)); + n += p.right_end - p.right_begin; + } + return n; +} + +Range Series::reduce(u64 begin, u64 end) const { + const u64 e = std::min(end, size()); + if (begin >= e) return Range{0.0, 0.0}; + const QueryParts p = split_query(begin, e); + + Range r{}; + bool first = true; + if (p.interior_begin < p.interior_end) { + r = pyramid().reduce_buckets(p.interior_begin, p.interior_end); + first = false; + } + + // At most 2*kBucket - 1 edge samples, so the buffer is enough without + // allocating: each edge span is under one bucket, and when the interior is + // empty the two spans tile a range shorter than two buckets. + f32 edge[2 * kBucket]; + const u64 n = copy_edges(p, edge); + for (u64 i = 0; i < n; ++i) { + const f64 v = f64(edge[i]); + if (first) { + r.min = v; + r.max = v; + first = false; + } else { + r.min = std::min(r.min, v); + r.max = std::max(r.max, v); + } + } + return first ? Range{0.0, 0.0} : r; +} + +f64 Series::sum(u64 begin, u64 end) const { + const u64 e = std::min(end, size()); + if (begin >= e) return 0.0; + const QueryParts p = split_query(begin, e); + + f64 total = 0.0; + if (p.interior_begin < p.interior_end) { + total = pyramid().sum_buckets(p.interior_begin, p.interior_end); + } + f32 edge[2 * kBucket]; + const u64 n = copy_edges(p, edge); + for (u64 i = 0; i < n; ++i) total += f64(edge[i]); + return total; +} + +f64 Series::sumsq(u64 begin, u64 end) const { + const u64 e = std::min(end, size()); + if (begin >= e) return 0.0; + const QueryParts p = split_query(begin, e); + + f64 total = 0.0; + if (p.interior_begin < p.interior_end) { + total = pyramid().sumsq_buckets(p.interior_begin, p.interior_end); + } + f32 edge[2 * kBucket]; + const u64 n = copy_edges(p, edge); + for (u64 i = 0; i < n; ++i) { + const f64 v = f64(edge[i]); + total += v * v; + } + return total; +} + +void Series::gather_impl(u32 level, u64 begin, u64 end, std::vector& out_xy) const { + const u64 n = size(); + const u64 b = std::min(begin, n); + const u64 e = std::min(std::max(end, b), n); + if (b >= e) return; + + const f64 x_0 = x0(); + const f64 step = dx(); + if (level == 0) { + f32 raw[kGatherChunk]; + for (u64 i = b; i < e;) { + const u64 take = std::min(kGatherChunk, e - i); + copy_raw(i, i + take, raw); + for (u64 k = 0; k < take; ++k) { + out_xy.push_back(f32(x_0 + f64(i + k) * step)); + out_xy.push_back(raw[k]); + } + i += take; + } + return; + } + + // Level >= 1: one envelope pair per bucket, x at the bucket's first sample + // (min) and at its last sample (max), which keeps x monotonic while covering + // the bucket's span. + const u64 span = u64{kBucket} << level; + const std::span buckets = pyramid().level(level); + const u64 first = b / span; + const u64 last = std::min((e - 1) / span + 1, static_cast(buckets.size())); + for (u64 bi = first; bi < last; ++bi) { + const MinMax m = buckets[static_cast(bi)]; + const u64 first_sample = bi * span; + const u64 last_sample = std::min(first_sample + span, n) - 1; + out_xy.push_back(f32(x_0 + f64(first_sample) * step)); + out_xy.push_back(m.min); + out_xy.push_back(f32(x_0 + f64(last_sample) * step)); + out_xy.push_back(m.max); + } +} + +} // namespace plot::data diff --git a/src/data/stream_series.cpp b/src/data/stream_series.cpp new file mode 100644 index 0000000..788dd55 --- /dev/null +++ b/src/data/stream_series.cpp @@ -0,0 +1,115 @@ +#include + +#include +#include + +namespace plot::data { +namespace { + +/// Samples copied per ring read on the ingest path, so rebuilding a saturated +/// window costs no allocation. +constexpr u64 kReadChunk = 256; + +} // namespace + +StreamSeries::StreamSeries(u64 capacity, f64 x0, f64 dx) + : ring_(capacity), x0_(x0), dx_(dx) {} + +void StreamSeries::push(f32 y) { ring_.push(y); } + +void StreamSeries::push(const f32* ys, u64 n) { ring_.push(ys, n); } + +u64 StreamSeries::produced() const { return ring_.produced(); } + +void StreamSeries::acquire() { + // One atomic read: every absolute index below is relative to this single count, + // so the producer may keep pushing without the snapshot shifting under us. + const u64 p = ring_.produced(); + const u64 cap = ring_.capacity(); + const u64 base = p > cap ? p - cap : 0; // oldest index the ring still holds + const u64 n = p - base; // retained, published sample count + + if (base != base_) { + // The window slid: samples the pyramid covered are gone, and because bucket + // boundaries are absolute, extending it incrementally would mis-align every + // bucket. Rebuild from the retained samples instead. + pyramid_.clear(); + pyramid_.reserve(n); + f32 lo = 0.0f, hi = 0.0f; + bool first = true; + f32 chunk[kReadChunk]; + for (u64 i = 0; i < n;) { + const u64 take = std::min(kReadChunk, n - i); + ring_.copy_from(base + i, take, chunk); + for (u64 k = 0; k < take; ++k) { + const f32 v = chunk[k]; + pyramid_.append(v); + if (first) { + lo = v; + hi = v; + first = false; + } else { + lo = std::min(lo, v); + hi = std::max(hi, v); + } + } + i += take; + } + y_bounds_ = first ? Range{0.0, 0.0} : Range{f64(lo), f64(hi)}; + base_ = base; + } else { + // Contiguous: only the samples produced since the last acquire() are new. + if (n == 0) y_bounds_ = Range{0.0, 0.0}; + f32 chunk[kReadChunk]; + for (u64 i = end_ - base_; i < n;) { + const u64 take = std::min(kReadChunk, n - i); + ring_.copy_from(base + i, take, chunk); + for (u64 k = 0; k < take; ++k) { + const f32 v = chunk[k]; + pyramid_.append(v); + if (i + k == 0) { + y_bounds_ = Range{f64(v), f64(v)}; + } else { + y_bounds_.min = std::min(y_bounds_.min, f64(v)); + y_bounds_.max = std::max(y_bounds_.max, f64(v)); + } + } + i += take; + } + } + end_ = p; +} + +u64 StreamSeries::size() const { return end_ - base_; } + +f64 StreamSeries::x0() const { return x0_ + f64(base_) * dx_; } + +f64 StreamSeries::dx() const { return dx_; } + +Range StreamSeries::x_bounds() const { + const f64 n = f64(size()); + return n == 0.0 ? Range{x0(), x0()} : Range{x0(), x0() + (n - 1.0) * dx_}; +} + +Range StreamSeries::y_bounds() const { return y_bounds_; } + +const MinMaxPyramid& StreamSeries::pyramid() const { return pyramid_; } + +void StreamSeries::gather(u32 level, u64 begin, u64 end, std::vector& out_xy) const { + gather_impl(level, begin, end, out_xy); +} + +void StreamSeries::copy_raw(u64 begin, u64 end, f32* out) const { + const u64 want = end > begin ? end - begin : 0; + if (want == 0) return; + const u64 n = size(); + const u64 b = std::min(begin, n); + const u64 e = std::min(std::max(end, b), n); + // Retained index b is absolute base_ + b: fixed for the whole copy, so a + // pushing producer cannot shift the window midway. + ring_.copy_from(base_ + b, e - b, out); + std::fill(out + static_cast(e - b), + out + static_cast(want), 0.0f); +} + +} // namespace plot::data diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 1ed6230..4a0309f 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -26,6 +26,9 @@ add_executable(plot_tests main.cpp test_transform.cpp test_gl_context.cpp + test_ring_buffer.cpp + test_minmax_pyramid.cpp + test_series.cpp egl_headless.cpp ) target_link_libraries(plot_tests PRIVATE plot::core ${plot_egl_library}) diff --git a/tests/harness.hpp b/tests/harness.hpp index 5a87ba0..4e1aa89 100644 --- a/tests/harness.hpp +++ b/tests/harness.hpp @@ -35,6 +35,23 @@ struct Failure { std::string msg; }; +/// Compares two doubles within `eps`. `PLOT_CHECK_NEAR(a, b)` leaves it to the +/// default below; pass an explicit tolerance for computed values. +inline constexpr double kNearEps = 1e-9; + +inline void check_near(double a, double b, double eps, const char* file, int line, + const char* lhs, const char* rhs) { + if (!(std::fabs(a - b) <= eps)) { + throw Failure{std::string(file) + ":" + std::to_string(line) + " CHECK_NEAR(" + lhs + ", " + + rhs + ") got " + std::to_string(a) + " vs " + std::to_string(b)}; + } +} + +inline void check_near(double a, double b, const char* file, int line, const char* lhs, + const char* rhs) { + check_near(a, b, kNearEps, file, line, lhs, rhs); +} + inline int run_all() { int failed = 0; for (auto& c : cases()) { @@ -64,11 +81,5 @@ inline int run_all() { " CHECK(" #cond ")"}; \ } while (0) -#define PLOT_CHECK_NEAR(a, b, eps) \ - do { \ - const double _a = (a), _b = (b); \ - if (!(std::fabs(_a - _b) <= (eps))) \ - throw ::plottest::Failure{std::string(__FILE__) + ":" + std::to_string(__LINE__) + \ - " CHECK_NEAR(" #a ", " #b ") got " + std::to_string(_a) + \ - " vs " + std::to_string(_b)}; \ - } while (0) +#define PLOT_CHECK_NEAR(a, b, ...) \ + ::plottest::check_near((a), (b) __VA_OPT__(,) __VA_ARGS__, __FILE__, __LINE__, #a, #b) diff --git a/tests/test_minmax_pyramid.cpp b/tests/test_minmax_pyramid.cpp new file mode 100644 index 0000000..432a5d5 --- /dev/null +++ b/tests/test_minmax_pyramid.cpp @@ -0,0 +1,106 @@ +#include "harness.hpp" +#include +#include + +using plot::data::kBucket; +using plot::data::MinMaxPyramid; +using plot::Range; + +PLOT_TEST(pyramid_bucket_extremes_preserved) { + MinMaxPyramid p; + p.reserve(1000); + for (int i = 0; i < 1000; ++i) p.append(i == 617 ? -50.0f : (i == 42 ? 900.0f : 1.0f)); + PLOT_CHECK(p.size() == 1000); + const Range r = p.reduce_buckets(0, 960); // exactly 15 whole buckets + PLOT_CHECK_NEAR(r.min, -50.0, 1e-6); + PLOT_CHECK_NEAR(r.max, 900.0, 1e-6); +} + +PLOT_TEST(pyramid_reduce_buckets_matches_bruteforce_on_aligned_ranges) { + std::mt19937 rng(1234); + std::uniform_real_distribution d(-100.0f, 100.0f); + std::vector v(5000); + MinMaxPyramid p; + for (auto& x : v) { x = d(rng); p.append(x); } + for (u64 a = 0; a < 5000; a += kBucket * 3) { + for (u64 b = a + kBucket; b <= 5000; b += kBucket * 7) { + const u64 ahi = (b / kBucket) * kBucket; // whole buckets only + if (ahi <= a) continue; + float lo = v[a], hi = v[a]; + for (u64 i = a; i < ahi; ++i) { lo = std::min(lo, v[i]); hi = std::max(hi, v[i]); } + const Range r = p.reduce_buckets(a, ahi); + PLOT_CHECK_NEAR(r.min, lo, 1e-4); + PLOT_CHECK_NEAR(r.max, hi, 1e-4); + } + } +} + +PLOT_TEST(pyramid_reduce_buckets_excludes_partial_edges) { + // Bucket 1 covers samples [64,128). reduce_buckets answers in whole buckets + // only, so the spike must show up in its own bucket and nowhere else. + MinMaxPyramid p; + for (int i = 0; i < 128; ++i) p.append(i == 70 ? 999.0f : 1.0f); + PLOT_CHECK_NEAR(p.reduce_buckets(64, 128).max, 999.0, 1e-6); + PLOT_CHECK_NEAR(p.reduce_buckets(0, 64).max, 1.0, 1e-6); + PLOT_CHECK(p.reduce_buckets(64, 64).span() == 0.0); // empty -> degenerate +} + +PLOT_TEST(pyramid_aligned_queries_never_include_a_partial_bucket) { + // Corpus of 100 samples: bucket 0 is whole, bucket 1 (samples 64..127) is + // partial. A whole-bucket query must not reach into the partial one, in + // either direction. + MinMaxPyramid p; + for (int i = 0; i < 100; ++i) p.append(i == 99 ? 777.0f : 1.0f); + PLOT_CHECK(p.reduce_buckets(64, 100).span() == 0.0); // no whole bucket inside + PLOT_CHECK_NEAR(p.reduce_buckets(0, 100).max, 1.0, 1e-6); // bucket 0 only + PLOT_CHECK(p.sum_buckets(64, 100) == 0.0); + PLOT_CHECK(p.sumsq_buckets(64, 100) == 0.0); +} + +PLOT_TEST(pyramid_level_selection) { + MinMaxPyramid p; + p.reserve(1u << 20); + for (u64 i = 0; i < (1u << 20); ++i) p.append(0.0f); + // <= kBucket samples per pixel must use level 0 (raw path handles it anyway). + PLOT_CHECK(p.level_for(1u << 20, 1.0) == 0); + // 512 samples per pixel with kBucket == 64 → bucket size 512 → level 3. + PLOT_CHECK(p.level_for(1u << 20, 512.0) == 3); + // Level selection must never exceed the deepest level that has data. + PLOT_CHECK(p.level_for(1000, 1e9) <= p.level_count() - 1); +} + +PLOT_TEST(pyramid_snap_returns_nearest_sample) { + MinMaxPyramid p; + for (int i = 0; i < 100; ++i) p.append(float(i)); + PLOT_CHECK(p.snap(10.4, 0.0, 1.0) == 10); + PLOT_CHECK(p.snap(10.6, 0.0, 1.0) == 11); + PLOT_CHECK(p.snap(-5.0, 0.0, 1.0) == 0); // clamped + PLOT_CHECK(p.snap(1e9, 0.0, 1.0) == 99); // clamped +} + +PLOT_TEST(pyramid_sums_match_bruteforce_on_aligned_ranges) { + std::mt19937 rng(99); + std::uniform_real_distribution d(-3.0f, 7.0f); + std::vector v(3000); + MinMaxPyramid p; + for (auto& x : v) { x = d(rng); p.append(x); } + for (u64 a = 0; a < 3000; a += kBucket * 5) { + for (u64 b = a + kBucket; b <= 3000; b += kBucket * 9) { + double s = 0.0, sq = 0.0; + for (u64 i = a; i < b; ++i) { s += v[i]; sq += double(v[i]) * double(v[i]); } + PLOT_CHECK_NEAR(p.sum_buckets(a, b), s, 1e-3); + PLOT_CHECK_NEAR(p.sumsq_buckets(a, b), sq, 1e-2); + } + } +} + +PLOT_TEST(pyramid_clear_and_rebuild) { + MinMaxPyramid p; + for (int i = 0; i < 500; ++i) p.append(1.0f); + p.clear(); + PLOT_CHECK(p.size() == 0); + PLOT_CHECK(p.level_count() <= 1); + p.append(7.0f); + PLOT_CHECK(p.size() == 1); + PLOT_CHECK_NEAR(p.reduce(0, 1).max, 7.0, 1e-6); +} diff --git a/tests/test_ring_buffer.cpp b/tests/test_ring_buffer.cpp new file mode 100644 index 0000000..4635992 --- /dev/null +++ b/tests/test_ring_buffer.cpp @@ -0,0 +1,19 @@ +#include "harness.hpp" +#include + +PLOT_TEST(ring_wraps_and_keeps_latest) { + plot::data::RingBuffer r(8); // capacity rounded up to a power of two + for (int i = 0; i < 20; ++i) r.push(f32(i)); + PLOT_CHECK(r.size() == 8); + PLOT_CHECK_NEAR(r[0], 12.0); // oldest retained + PLOT_CHECK_NEAR(r[7], 19.0); // newest +} +PLOT_TEST(ring_clear_resets) { + plot::data::RingBuffer r(8); + for (int i = 0; i < 5; ++i) r.push(f32(i)); + r.clear(); + PLOT_CHECK(r.size() == 0); + r.push(42.0f); + PLOT_CHECK(r.size() == 1); + PLOT_CHECK_NEAR(r[0], 42.0); +} diff --git a/tests/test_series.cpp b/tests/test_series.cpp new file mode 100644 index 0000000..93a5b0d --- /dev/null +++ b/tests/test_series.cpp @@ -0,0 +1,78 @@ +#include "harness.hpp" +#include +#include +#include +#include +#include +#include + +PLOT_TEST(array_series_bounds_and_gather) { + std::vector y{1, 5, -2, 3}; + plot::data::ArraySeries s(std::move(y), /*x0=*/10.0, /*dx=*/0.5); + PLOT_CHECK(s.size() == 4); + PLOT_CHECK_NEAR(s.x_bounds().min, 10.0, 1e-9); + PLOT_CHECK_NEAR(s.x_bounds().max, 11.5, 1e-9); + PLOT_CHECK_NEAR(s.y_bounds().min, -2.0, 1e-6); + std::vector out; + s.gather(0, 0, 4, out); + PLOT_CHECK(out.size() == 8); // interleaved x,y + PLOT_CHECK_NEAR(out[0], 10.0, 1e-6); + PLOT_CHECK_NEAR(out[1], 1.0, 1e-6); + PLOT_CHECK_NEAR(out[6], 11.5, 1e-6); +} + +PLOT_TEST(stream_series_acquire_publishes_snapshot) { + plot::data::StreamSeries s(1 << 12, /*x0=*/0.0, /*dx=*/1e-6); + for (int i = 0; i < 100; ++i) s.push(f32(i)); + PLOT_CHECK(s.produced() == 100); + s.acquire(); + PLOT_CHECK(s.size() == 100); + // Pyramid must track the published samples. + PLOT_CHECK_NEAR(s.pyramid().reduce(0, 100).max, 99.0, 1e-4); + for (int i = 100; i < 250; ++i) s.push(f32(i)); + s.acquire(); + PLOT_CHECK(s.size() == 250); +} + +PLOT_TEST(stream_series_drops_oldest_when_full) { + plot::data::StreamSeries s(256, 0.0, 1.0); + for (int i = 0; i < 1000; ++i) s.push(f32(i)); + s.acquire(); + PLOT_CHECK(s.size() <= 256); + PLOT_CHECK_NEAR(s.y_bounds().max, 999.0, 1e-4); +} + +PLOT_TEST(array_series_exact_range_queries_match_bruteforce) { + // The strong test for the exactness split: arbitrary, deliberately + // unaligned ranges must match a brute-force scan of the raw samples. + std::mt19937 rng(7); + std::uniform_real_distribution d(-5.0f, 5.0f); + std::vector y(4000); + for (auto& v : y) v = d(rng); + plot::data::ArraySeries s(y, 0.0, 1.0); + for (u64 a = 0; a < 4000; a += 173) { + for (u64 b = a + 1; b <= 4000; b += 331) { + float lo = y[a], hi = y[a]; + double sum = 0.0, sq = 0.0; + for (u64 i = a; i < b; ++i) { + lo = std::min(lo, y[i]); hi = std::max(hi, y[i]); + sum += y[i]; sq += double(y[i]) * double(y[i]); + } + const plot::Range r = s.reduce(a, b); + PLOT_CHECK_NEAR(r.min, lo, 1e-4); + PLOT_CHECK_NEAR(r.max, hi, 1e-4); + PLOT_CHECK_NEAR(s.sum(a, b), sum, 1e-2); + PLOT_CHECK_NEAR(s.sumsq(a, b), sq, 1e-1); + } + } +} + +PLOT_TEST(stream_series_exact_queries_survive_ring_wrap) { + plot::data::StreamSeries s(1u << 10, 0.0, 1.0); + for (int i = 0; i < 3000; ++i) s.push(f32(i % 97)); // wraps the ring twice + s.acquire(); + const u64 n = s.size(); + const plot::Range r = s.reduce(0, n); // retained window, wrapped + PLOT_CHECK_NEAR(r.max, 96.0, 1e-4); + PLOT_CHECK_NEAR(r.min, 0.0, 1e-4); +}