Files
CoopAllTheThings/common/include/coop/tone_analysis.hpp
BlackMark 05039ab104 Clean up common/ comments and the log-ring strncpy warning
Comments must not document the past or reference plan circumstances:
drop the stale pointer to a never-created audio_correlate_layout.hpp,
the "step b" plan labels, the "carved out of reserved space" history
note, and a README pointer; reword a past-tense seqlock comment to
describe the failure mode in the present.

Replace the strncpy in log_ring_push with a bounded memcpy: same
semantics (truncate + NUL), but without the C4996 deprecation warning
on every host build.
2026-07-12 08:31:48 +02:00

367 lines
12 KiB
C++

// Quantitative fidelity analysis of captured audio, built to diagnose the exact
// failure modes the audio-mirror path can introduce: a wrong/guessed sample rate
// (pitch shift), periodic under-run/re-prime gaps (a "metallic" / choppy timbre),
// and dropped packets (clicks). It turns "the audio sounds slightly off" into
// numbers: pitch error in cents, SNR/THD, a click count, and a dropout count.
//
// The pitch/SNR/THD metrics assume the input is a single sine tone of a known
// frequency (use coop_tone as the source) -- that's the controlled signal that
// makes pitch error measurable. The discontinuity (click) and dropout metrics are
// content-agnostic, so they still mean something on a recording of real game audio.
//
// Header-only, no Windows / no audio-device dependency, so it is unit-tested with
// synthesized adversarial signals (tests/tone_analysis_test.cpp) and reused by the
// coop_audio_validate tool.
#pragma once
#include <algorithm>
#include <cmath>
#include <complex>
#include <cstddef>
#include <cstdint>
#include <vector>
namespace coop
{
// WAVE_FORMAT_* tags we decode (kept local to avoid an mmreg.h dependency, matching audio_mix.hpp).
inline constexpr std::uint32_t kToneFormatPcm = 1;
inline constexpr std::uint32_t kToneFormatFloat = 3;
// One channel's worth of measured fidelity. Fields are NaN/0 when not applicable
// (e.g. pitch metrics need a known expected_hz > 0).
struct ToneReport
{
bool valid = false; // enough samples to analyze
unsigned sample_rate = 0; // the rate the samples are interpreted at (the *declared* rate)
std::size_t frames = 0; // mono frames analyzed
double duration_sec = 0.0;
// --- Level ---
double rms = 0.0; // 0..1
double peak = 0.0; // 0..1
double clipped_fraction = 0.0; // fraction of samples at >= 0.999 full-scale
// --- Pitch (needs a known input tone frequency) ---
double expected_hz = 0.0; // the tone frequency that was played
double dominant_hz = 0.0; // the fundamental we recovered
double pitch_error_ratio = 0.0; // dominant / expected (1.0 = perfect)
double pitch_error_cents = 0.0; // 1200*log2(ratio); +/- ~10 cents starts to be audible
// --- Spectral purity (tone mode) ---
double snr_db = 0.0; // fundamental power vs everything else (DC + harmonics excluded from "signal")
double thd_percent = 0.0; // harmonics 2..6 vs fundamental
// --- Time-domain defects (content-agnostic) ---
unsigned glitch_count = 0; // discontinuity events (clicks): big isolated sample jumps
double glitch_rate_per_sec = 0.0;
unsigned dropout_count = 0; // gaps: stretches that fall near-silent mid-signal
double dropout_ms = 0.0; // total duration of those gaps
};
namespace detail
{
inline constexpr double kPi = 3.14159265358979323846;
// In-place iterative radix-2 Cooley-Tukey FFT. `a.size()` must be a power of two.
inline void fft(std::vector<std::complex<double>>& a)
{
const std::size_t n = a.size();
for (std::size_t i = 1, j = 0; i < n; ++i)
{
std::size_t bit = n >> 1;
for (; (j & bit) != 0; bit >>= 1)
{
j ^= bit;
}
j ^= bit;
if (i < j)
{
std::swap(a[i], a[j]);
}
}
for (std::size_t len = 2; len <= n; len <<= 1)
{
const double ang = -2.0 * kPi / static_cast<double>(len);
const std::complex<double> wlen(std::cos(ang), std::sin(ang));
for (std::size_t i = 0; i < n; i += len)
{
std::complex<double> w(1.0, 0.0);
for (std::size_t k = 0; k < len / 2; ++k)
{
const std::complex<double> u = a[i + k];
const std::complex<double> v = a[i + k + len / 2] * w;
a[i + k] = u + v;
a[i + k + len / 2] = u - v;
w *= wlen;
}
}
}
}
// Largest power of two <= n (0 for n==0).
inline std::size_t floor_pow2(std::size_t n)
{
std::size_t p = 1;
while ((p << 1) != 0 && (p << 1) <= n)
{
p <<= 1;
}
return n == 0 ? 0 : p;
}
} // namespace detail
// Decode one channel (default: channel 0) of interleaved PCM into normalized [-1,1]
// floats. Supports float32 and int16 (the formats the mirror's mixer handles); returns
// empty for anything else. `bytes` is the byte length of `pcm`.
inline std::vector<float> decode_channel(const std::uint8_t* pcm, std::size_t bytes, std::uint32_t format_tag,
std::uint32_t bits, std::uint32_t channels, std::uint32_t channel = 0)
{
std::vector<float> out;
if (pcm == nullptr || channels == 0 || channel >= channels)
{
return out;
}
if (format_tag == kToneFormatFloat && bits == 32)
{
const std::size_t frames = bytes / (channels * 4);
out.reserve(frames);
const auto* f = reinterpret_cast<const float*>(pcm);
for (std::size_t i = 0; i < frames; ++i)
{
out.push_back(f[i * channels + channel]);
}
}
else if (format_tag == kToneFormatPcm && bits == 16)
{
const std::size_t frames = bytes / (channels * 2);
out.reserve(frames);
const auto* s = reinterpret_cast<const std::int16_t*>(pcm);
for (std::size_t i = 0; i < frames; ++i)
{
out.push_back(static_cast<float>(s[i * channels + channel]) / 32768.0f);
}
}
return out;
}
// Analyze a single channel of normalized float samples. `expected_hz` is the known
// input tone frequency (pass 0 to skip the pitch/SNR/THD metrics for non-tone audio;
// the click/dropout/level metrics still apply).
inline ToneReport analyze_tone(const float* samples, std::size_t frames, unsigned sample_rate,
double expected_hz)
{
ToneReport r;
r.sample_rate = sample_rate;
r.frames = frames;
r.expected_hz = expected_hz;
if (samples == nullptr || frames < 64 || sample_rate == 0)
{
return r;
}
r.valid = true;
r.duration_sec = static_cast<double>(frames) / sample_rate;
// --- Level: RMS, peak, clipping ---
double sumsq = 0.0;
double peak = 0.0;
std::size_t clipped = 0;
for (std::size_t i = 0; i < frames; ++i)
{
const double x = samples[i];
sumsq += x * x;
const double a = std::fabs(x);
peak = std::max(peak, a);
if (a >= 0.999)
{
++clipped;
}
}
r.rms = std::sqrt(sumsq / static_cast<double>(frames));
r.peak = peak;
r.clipped_fraction = static_cast<double>(clipped) / static_cast<double>(frames);
// --- Discontinuity (click) detection: content-agnostic ---
// A click is an isolated sample-to-sample jump far larger than the signal's typical
// step. Use the median |first difference| as a robust scale (immune to the tone's own
// slope and to a few outliers), and flag steps beyond 8x it. Group samples within a
// short refractory window into one event so a single click isn't counted many times.
if (frames >= 3)
{
std::vector<float> diff(frames - 1);
for (std::size_t i = 1; i < frames; ++i)
{
diff[i - 1] = std::fabs(samples[i] - samples[i - 1]);
}
std::vector<float> sorted(diff);
std::nth_element(sorted.begin(), sorted.begin() + sorted.size() / 2, sorted.end());
const double median = sorted[sorted.size() / 2];
const double thresh = std::max(8.0 * median, 0.02 * std::max(peak, 1e-6));
const std::size_t refractory = std::max<std::size_t>(sample_rate / 1000, 8); // ~1 ms
std::size_t last_event_end = 0;
bool have_event = false;
for (std::size_t i = 0; i < diff.size(); ++i)
{
if (diff[i] > thresh)
{
if (!have_event || i > last_event_end)
{
++r.glitch_count;
}
have_event = true;
last_event_end = i + refractory;
}
}
r.glitch_rate_per_sec = r.glitch_count / std::max(r.duration_sec, 1e-9);
}
// --- Dropout detection: stretches that fall near-silent in an otherwise active signal ---
// Slide a ~5 ms window; flag windows whose RMS drops below 8% of the global RMS. Only
// meaningful when the signal is actually present (global RMS above a small floor).
if (r.rms > 1e-4)
{
const std::size_t win = std::max<std::size_t>(sample_rate * 5 / 1000, 16); // ~5 ms
const std::size_t hop = std::max<std::size_t>(win / 2, 1);
const double silence_thresh = 0.08 * r.rms;
bool in_gap = false;
std::size_t gap_first = 0; // first silent window's start sample
std::size_t gap_last = 0; // last silent window's end sample
std::size_t total_silent_samples = 0;
auto close_gap = [&]() {
if (in_gap)
{
total_silent_samples += (gap_last - gap_first);
in_gap = false;
}
};
for (std::size_t start = 0; start + win <= frames; start += hop)
{
double ws = 0.0;
for (std::size_t i = 0; i < win; ++i)
{
const double x = samples[start + i];
ws += x * x;
}
const double wr = std::sqrt(ws / static_cast<double>(win));
if (wr < silence_thresh)
{
if (!in_gap)
{
++r.dropout_count;
gap_first = start;
in_gap = true;
}
gap_last = start + win;
}
else
{
close_gap();
}
}
close_gap();
// Gap extent (first silent window start .. last silent window end) -- more faithful
// than counting interior windows, which drops the partially-silent edge windows.
r.dropout_ms = static_cast<double>(total_silent_samples) * 1000.0 / sample_rate;
}
// --- Spectral analysis (pitch / SNR / THD), Hann-windowed FFT ---
if (expected_hz > 0.0)
{
std::size_t n = detail::floor_pow2(frames);
n = std::min<std::size_t>(n, std::size_t(1) << 18); // cap cost (~5 s @ 48k)
if (n >= 1024)
{
std::vector<std::complex<double>> buf(n);
for (std::size_t i = 0; i < n; ++i)
{
const double w = 0.5 - 0.5 * std::cos(2.0 * detail::kPi * i / (n - 1)); // Hann
buf[i] = std::complex<double>(samples[i] * w, 0.0);
}
detail::fft(buf);
const std::size_t half = n / 2;
std::vector<double> mag(half);
for (std::size_t i = 0; i < half; ++i)
{
mag[i] = std::abs(buf[i]);
}
const double bin_hz = static_cast<double>(sample_rate) / static_cast<double>(n);
// Peak bin, ignoring DC/very low bins (skip < 20 Hz).
std::size_t lo = std::max<std::size_t>(static_cast<std::size_t>(20.0 / bin_hz), 1);
std::size_t peak_bin = lo;
for (std::size_t i = lo; i < half; ++i)
{
if (mag[i] > mag[peak_bin])
{
peak_bin = i;
}
}
// Quadratic (parabolic) interpolation on log-magnitude for a sub-bin estimate
// (accurate for a Hann-windowed peak).
double delta = 0.0;
if (peak_bin > 0 && peak_bin + 1 < half)
{
const double a = std::log(mag[peak_bin - 1] + 1e-30);
const double b = std::log(mag[peak_bin] + 1e-30);
const double c = std::log(mag[peak_bin + 1] + 1e-30);
const double denom = (a - 2.0 * b + c);
if (std::fabs(denom) > 1e-30)
{
delta = 0.5 * (a - c) / denom;
delta = std::max(-0.5, std::min(0.5, delta));
}
}
r.dominant_hz = (static_cast<double>(peak_bin) + delta) * bin_hz;
if (r.dominant_hz > 0.0)
{
r.pitch_error_ratio = r.dominant_hz / expected_hz;
r.pitch_error_cents = 1200.0 * std::log2(r.pitch_error_ratio);
}
// Power in a +/-half_w bin lobe around a target frequency. A bin-centered Hann
// peak is 3 bins wide, but an off-bin tone leaks into a wider skirt, so the
// fundamental lobe is generous (+/-8) to contain that skirt -- otherwise the
// leakage masquerades as noise and floors a clean tone's SNR (~35 dB instead of
// >60). Harmonic lobes stay tight (+/-3).
auto lobe_power = [&](double hz, long half_w) {
const long center = static_cast<long>(std::lround(hz / bin_hz));
double p = 0.0;
for (long k = center - half_w; k <= center + half_w; ++k)
{
if (k >= 0 && static_cast<std::size_t>(k) < half)
{
p += mag[k] * mag[k];
}
}
return p;
};
double total_power = 0.0;
for (std::size_t i = lo; i < half; ++i)
{
total_power += mag[i] * mag[i];
}
const double fund_power = lobe_power(r.dominant_hz, 8);
const double residual = std::max(total_power - fund_power, 1e-30);
r.snr_db = 10.0 * std::log10(std::max(fund_power, 1e-30) / residual);
double harm_power = 0.0;
for (int h = 2; h <= 6; ++h)
{
const double hz = r.dominant_hz * h;
if (hz < (sample_rate / 2.0))
{
harm_power += lobe_power(hz, 3);
}
}
r.thd_percent = 100.0 * std::sqrt(harm_power / std::max(fund_power, 1e-30));
}
}
return r;
}
} // namespace coop