plot.py

#!/usr/bin/env python3
"""Plot mwbench output - samples.csv, amplification.csv, summary.csv,
delete_experiment.csv >> PNG + PDF figures.

Multi-CF aware: db-wide stats are plotted once (read from cf_index == 0 rows
since they are duplicated). Per-CF stats are overlaid across CFs or
faceted into one subplot per CF.

usage: plot.py [out_dir]   (defaults to newest ./out/run_* next to this script)
"""

import csv
import os
import sys

import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
from matplotlib.ticker import FuncFormatter

# ============================================================ palette
#
# pastel scheme derived from the seaborn "pastel" palette
# (which itself extends ColorBrewer Pastel1) -- chosen because it (a) reads
# softly enough to compose multiple overlaid series without visual fatigue,
# (b) photocopies / black-and-white reduces to distinguishable greys, and
# (c) is colorblind-friendly for the most common deficiencies (deuteran/
# protanopia preserve hue spread; tritanopia loses some but values remain
# distinguishable by lightness).
#
# the palette is two-tier:
#   PASTEL_*       -- soft fills + line strokes for the main signal
#   DEEP_*         -- saturated counterparts for emphasis (titles, refs,
#                     annotation borders). darker so they read against
#                     pastel-shaded backgrounds.

PASTEL_BLUE     = "#A1C9F4"
PASTEL_PEACH    = "#FFB482"
PASTEL_MINT     = "#8DE5A1"
PASTEL_SALMON   = "#FF9F9B"
PASTEL_LAVENDER = "#D0BBFF"
PASTEL_TAN      = "#DEBB9B"
PASTEL_PINK     = "#FAB0E4"
PASTEL_GREY     = "#CFCFCF"
PASTEL_CREAM    = "#FFFEA3"
PASTEL_CYAN     = "#B9F2F0"

DEEP_BLUE   = "#4C72B0"
DEEP_PEACH  = "#DD8452"
DEEP_MINT   = "#55A868"
DEEP_RED    = "#C44E52"
DEEP_PURPLE = "#8172B2"
DEEP_BROWN  = "#937860"
DEEP_PINK   = "#DA8BC3"
DEEP_GREY   = "#797979"
DEEP_CYAN   = "#64B5CD"  # rounding out the deep set

INK         = "#2C2C2C"  # primary axis / title color
PAPER       = "#FCFCFC"  # very light off-white figure background
RULE        = "#888888"  # mid-grey for reference lines / SLO dashes

# ============================================================ rcParams

plt.rcParams.update({
    # typography 
    "font.family":     "serif",
    "font.size":       10.5,
    "axes.titlesize":  12,
    "axes.labelsize":  10.5,
    "axes.titleweight": "normal",
    "axes.titlecolor":  INK,
    "axes.labelcolor":  INK,
    "xtick.color":      INK,
    "ytick.color":      INK,
    "xtick.labelsize":  9,
    "ytick.labelsize":  9,
    "legend.fontsize":  9,
    # in-axes legends get a translucent paper background so that when a
    # legend unavoidably sits over plotted lines it stays readable instead
    # of the lines bleeding through the text. margin legends
    # (add_phase_legend, the per-CF figure legends) opt back out with an
    # explicit frameon=False since they sit in empty figure space.
    "legend.frameon":     True,
    "legend.framealpha":  0.85,
    "legend.facecolor":   PAPER,
    "legend.edgecolor":   PASTEL_GREY,
    # headroom above the data so corner legends / explainer boxes have empty
    # space to sit in rather than landing on top of the topmost samples
    "axes.ymargin":       0.14,
    # canvas
    "figure.facecolor":  PAPER,
    "axes.facecolor":    PAPER,
    "savefig.facecolor": PAPER,
    "axes.edgecolor":    INK,
    "axes.linewidth":    0.7,
    "axes.spines.top":   False,   # clean modern frame: left+bottom only
    "axes.spines.right": False,
    "axes.grid":         True,
    "grid.color":        DEEP_GREY,
    "grid.alpha":        0.18,
    "grid.linewidth":    0.4,
    "grid.linestyle":    "-",
    # lines + markers 
    "lines.linewidth":   1.6,
    "lines.solid_capstyle": "round",
    "lines.dash_capstyle":  "round",
    # default cycler: pastel set used when nothing explicit is provided
    "axes.prop_cycle":   plt.cycler(color=[
        PASTEL_BLUE, PASTEL_PEACH, PASTEL_MINT, PASTEL_SALMON,
        PASTEL_LAVENDER, PASTEL_TAN, PASTEL_PINK, PASTEL_CYAN,
    ]),
    # output 
    "savefig.dpi":       200,    # bumped for print
    "savefig.bbox":      "tight",
    "figure.dpi":        140,
    "pdf.fonttype":      42,     # embedded TrueType (editable in vector tools)
    "ps.fonttype":       42,
})

# ============================================================ constants

# phase shade colors -- pastel set tuned so the 6 phases sit harmoniously
# even when several appear in the same figure. alpha at draw time is 0.40
# (set in shade_phases) so lines on top remain legible.
PHASES = {
    0: ("ingest",       PASTEL_BLUE),     # writing data
    1: ("cooldown",     PASTEL_GREY),     # quiescent
    2: ("delete",       PASTEL_SALMON),   # tombstoning -- warm warning
    3: ("post-delete",  PASTEL_PEACH),    # tombstones present, no compact
    4: ("compaction",   PASTEL_LAVENDER), # reclamation in flight
    5: ("post-compact", PASTEL_MINT),     # settled, clean
}

OP_TYPES   = ("point", "seek", "range")
# per-op stable colors so the same op reads the same hue everywhere
OP_COLORS  = {"point": DEEP_BLUE, "seek": DEEP_PEACH, "range": DEEP_MINT}
# per-percentile colors. p50 = calm green, p95 = warning orange,
# p99 = alarm red -- mirrors the visual "severity" intuition
PCT_COLORS = {"p50": DEEP_MINT,  "p95": DEEP_PEACH, "p99": DEEP_RED}

GIB = 1024 ** 3
MIB = 1024 ** 2

# common latency SLO references (µs) drawn as horizontal dashes
SLO_REFS_US = [(100, "100 µs"), (1000, "1 ms"), (10000, "10 ms")]

# ============================================================ io helpers

def load_csv(path):
    rows = []
    if not os.path.exists(path):
        return rows
    with open(path) as f:
        rdr = csv.DictReader(f)
        for r in rdr:
            for k, v in list(r.items()):
                if v is None or v == "":
                    r[k] = 0.0
                else:
                    try:
                        r[k] = float(v)
                    except ValueError:
                        pass
            rows.append(r)
    return rows


def col(rows, name):
    return [r.get(name, 0.0) for r in rows]


def partition_by_cf(rows):
    out = {}
    for r in rows:
        idx = int(r.get("cf_index", 0))
        out.setdefault(idx, []).append(r)
    return out


def load_config(run_dir):
    rows = []
    path = os.path.join(run_dir, "config.csv")
    if not os.path.exists(path):
        return rows
    with open(path) as f:
        rdr = csv.DictReader(f)
        for r in rdr:
            rows.append(r)
    return rows


def load_run_meta(run_dir):
    """run_meta.csv is a flat key,value table"""
    path = os.path.join(run_dir, "run_meta.csv")
    meta = {}
    if not os.path.exists(path):
        return meta
    with open(path) as f:
        rdr = csv.reader(f)
        for i, row in enumerate(rdr):
            if i == 0 or len(row) < 2: continue
            meta[row[0]] = row[1]
    return meta


def load_phase_bounds(run_dir):
    """build absolute (start_s, end_s, phase_id) bounds for the run from
    summary.csv. this is more reliable than inferring phase changes from
    per-row samples.csv inspection -- a phase that finished faster than
    one sample interval (e.g. a sub-second ingest on a small target) does
    not land any sample rows, so per-row inference would shade it out of
    existence. summary's per-phase duration_s captures those phases
    truthfully because the bench's per-phase aggregator measures wall
    time at the actual transition, not at the next sampler tick."""
    path = os.path.join(run_dir, "summary.csv")
    if not os.path.exists(path): return None
    bounds = []
    t = 0.0
    with open(path) as f:
        rdr = csv.DictReader(f)
        for r in rdr:
            try:
                ph  = int(float(r["phase"]))
                dur = float(r["duration_s"])
            except (KeyError, ValueError):
                continue
            if dur <= 0: continue
            bounds.append((t, t + dur, ph))
            t += dur
    return bounds or None


# module-global, set by main() once summary.csv is known. shade_phases
# prefers this when available; falls back to per-row inference otherwise.
PHASE_BOUNDS = None


def is_unified_mode(meta):
    """unified_memtable=1 in run_meta.csv (or implicit-on via objstore_mode).
    decides whether per-CF memtable stats are meaningful (they are not in
    unified mode -- the memtable lives in the db-wide aggregate only)."""
    if not meta: return False
    if meta.get("unified_memtable", "0") == "1": return True
    if meta.get("objstore_mode", "none") not in ("", "none"): return True
    return False


def caption_from_config(cfg_rows, meta=None):
    """two-line caption: line 1 = workload (patterns + cf count + memtable
    mode + compression), line 2 = environment (tidesdb version, OS, disk).
    uses run_meta.csv when available, falls back to config.csv only."""
    if not cfg_rows and not meta:
        return ""
    line1_bits = []
    have_num_cfs = False
    if meta:
        wp = meta.get("write_pattern", "")
        rp = meta.get("read_pattern", "")
        if wp: line1_bits.append(f"write={wp}")
        if rp: line1_bits.append(f"read={rp}")
        if wp == "zipfian" or rp == "zipfian":
            line1_bits.append(f"theta={meta.get('zipf_skew', '')}")
        if "num_cfs" in meta:
            line1_bits.append(f"num_cfs={meta['num_cfs']}")
            have_num_cfs = True
        if "value_size" in meta:
            line1_bits.append(f"value={meta['value_size']}B")
        # surface memtable mode -- the single most load-bearing piece of
        # context for interpreting per-CF memtable / immutable panels
        line1_bits.append(
            "memtable=unified" if is_unified_mode(meta) else "memtable=per-cf")
    if cfg_rows:
        n = len(cfg_rows)
        def uniq(field):
            return sorted({r.get(field, "") for r in cfg_rows})
        comps = uniq("compression")
        btree = uniq("use_btree")
        if not have_num_cfs:
            line1_bits.append(f"num_cfs={n}")
        line1_bits.append("compression=" + ("/".join(comps) if len(comps) > 1 else comps[0]))
        line1_bits.append("btree=" + ("mixed" if len(btree) > 1 else btree[0]))
    line1 = "  ·  ".join(line1_bits)

    if not meta:
        return line1
    tdb = meta.get("tidesdb_version", "?")
    if meta.get("tidesdb_has_s3", "0") == "1": tdb += "+s3"
    os_str = f"{meta.get('os_name','?')} {meta.get('os_release','?')}"
    cpu = meta.get("cpu_model", "")
    if len(cpu) > 35: cpu = cpu[:32] + "..."
    disk = f"{meta.get('data_dir_device','?')} ({meta.get('data_dir_filesystem','?')})"
    line2_bits = [
        f"tidesdb {tdb}",
        os_str,
        cpu,
        f"disk={disk}",
        f"compiler={meta.get('compiler', '?').split()[0]}",
    ]
    line2 = "  ·  ".join(b for b in line2_bits if b and not b.endswith("?"))
    return line1 + "\n" + line2


def save(fig, out_dir, name):
    fig.savefig(os.path.join(out_dir, f"{name}.png"))
    fig.savefig(os.path.join(out_dir, f"{name}.pdf"))
    plt.close(fig)


# ============================================================ shared decorators

def shade_phases(ax, t, phase):
    """draw phase bands behind the lines.

    if module-level PHASE_BOUNDS is set (loaded once from summary.csv),
    shade by absolute time bands so phases that finished faster than one
    sample interval (a sub-second ingest, for example) still appear.

    fallback: walk the (t, phase) column pair and shade by per-row
    inference. this only sees phases that actually landed a sample row."""
    # the bands only span the data range; matplotlib's default 5% x-padding
    # would otherwise show as an unshaded white strip at the frame edges that
    # reads like a phantom phase before the run started. drop the x-margin so
    # the shading reaches the left and right spines.
    ax.set_xmargin(0)
    if PHASE_BOUNDS:
        for s, e, ph in PHASE_BOUNDS:
            _, color = PHASES.get(ph, ("?", "#eeeeee"))
            ax.axvspan(s, e, color=color, alpha=0.38, lw=0, zorder=0)
        return
    if not t: return
    start = 0
    cur = int(phase[0])
    for i in range(1, len(t)):
        if int(phase[i]) != cur:
            _, color = PHASES.get(cur, ("?", "#eeeeee"))
            ax.axvspan(t[start], t[i], color=color, alpha=0.38, lw=0, zorder=0)
            start = i
            cur = int(phase[i])
    _, color = PHASES.get(cur, ("?", "#eeeeee"))
    ax.axvspan(t[start], t[-1], color=color, alpha=0.38, lw=0, zorder=0)


def shade_phases_x(ax, x_vals, phase):
    """phase shading using an arbitrary x-axis (e.g., bytes_written)"""
    if not x_vals: return
    # no x-padding so the bands reach the frame edges -- see shade_phases. this
    # matters most here: the x-axis is bytes_written and the first sample is
    # already several GiB in, so the default margin left a wide white gap.
    ax.set_xmargin(0)
    start = 0
    cur = int(phase[0])
    for i in range(1, len(x_vals)):
        if int(phase[i]) != cur:
            _, color = PHASES.get(cur, ("?", "#eeeeee"))
            ax.axvspan(x_vals[start], x_vals[i], color=color, alpha=0.38,
                       lw=0, zorder=0)
            start = i
            cur = int(phase[i])
    _, color = PHASES.get(cur, ("?", "#eeeeee"))
    ax.axvspan(x_vals[start], x_vals[-1], color=color, alpha=0.38, lw=0,
               zorder=0)


def add_phase_legend(fig, y=0.07):
    """phase legend at the bottom; the run-meta caption sits below it"""
    handles = [plt.Rectangle((0, 0), 1, 1, color=c, alpha=0.5)
               for _, (_, c) in sorted(PHASES.items())]
    labels = [n for _, (n, _) in sorted(PHASES.items())]
    fig.legend(handles, labels, loc="lower center", ncol=len(PHASES),
               fontsize=8, frameon=False, bbox_to_anchor=(0.5, y))


def add_caption(fig, caption, y=0.005):
    """caption may be multi-line. tight_layout's rect bottom should leave
    enough room for both lines + the phase legend above it"""
    if not caption: return
    fig.text(0.5, y, caption, ha="center", va="bottom", fontsize=7.5,
             color=DEEP_GREY, family="monospace", linespacing=1.4)


def finalize(fig, out_dir, name, caption, top=0.95, phase=True):
    """lay out the bottom band -- phase legend stacked above the run caption
    -- and run tight_layout so the lowest subplot's x-axis labels never land
    on either of them, then save.

    the band budget is reckoned in INCHES, not figure fractions, so the same
    physical strip is reserved whether the figure is a 5-inch single panel or
    a 10-inch three-panel stack. fixed-fraction placement was the source of
    the x-axis-over-legend and caption-over-legend collisions: on a short
    figure a 0.11 fraction is too few inches to hold xlabels + legend +
    a two-line caption without overlap."""
    fig_h = fig.get_figheight()
    cap_lines = (caption.count("\n") + 1) if caption else 0
    cap_base_in = 0.10                 # caption baseline above the bottom edge
    cap_block_in = 0.15 * cap_lines    # caption text height
    if phase:
        # full band: caption, gap, phase legend, then x-label clearance above
        gap_in = 0.16
        leg_in = 0.20
        xlabel_gap_in = 0.42
    else:
        # caption only -- the axes' own x-decorations are kept above the rect
        # bottom by tight_layout, so we just need a small gap over the caption
        gap_in = 0.0
        leg_in = 0.0
        xlabel_gap_in = 0.10
    leg_base_in = cap_base_in + cap_block_in + gap_in
    bottom_in = leg_base_in + leg_in + xlabel_gap_in
    if phase:
        add_phase_legend(fig, y=leg_base_in / fig_h)
    if caption:
        add_caption(fig, caption, y=cap_base_in / fig_h)
    bottom = min(bottom_in / fig_h, top - 0.05)
    fig.tight_layout(rect=(0, bottom, 1, top))
    save(fig, out_dir, name)


def add_slo_lines(ax, label_only=True):
    """horizontal SLO dashes on log-y latency plots. label_only=True writes
    the µs/ms labels at the right edge so they don't stack at the leftmost x;
    label_only=False draws the lines without text (used on non-first panels
    sharing a y-axis with the labeled one)"""
    for v, lab in SLO_REFS_US:
        ax.axhline(v, color=RULE, linewidth=0.5, linestyle=":",
                   alpha=0.7, zorder=0)
        if label_only:
            ax.text(0.995, v, lab + " ", va="bottom", ha="right",
                    fontsize=7.5, color=RULE,
                    transform=ax.get_yaxis_transform())


def add_explainer(ax, text, xy=(0.99, 0.97), ha="right", va="top"):
    """small italic text box pinned to a corner of the axes; explains the
    panel's metric in one or two lines.

    the box must read on top of every plotted line. within one axes artist
    zorder handles that (text=3 > line=2), but matplotlib renders overlapping
    axes -- a base axes and its twinx -- in creation order regardless of
    zorder, so a box on the base axes would sit under any line on the later
    twin. anchor the text on the topmost axes sharing this position (the twin
    when there is one) so it always draws last. transAxes fractions are shared
    across the twin, so the placement is identical."""
    same_pos = [a for a in ax.figure.axes
                if a.get_position().bounds == ax.get_position().bounds]
    target = same_pos[-1] if same_pos else ax  # fig.axes is creation order
    target.text(xy[0], xy[1], text, transform=target.transAxes,
                ha=ha, va=va, fontsize=8, style="italic", color=INK, zorder=6,
                bbox=dict(boxstyle="round,pad=0.3", facecolor=PAPER,
                          edgecolor=PASTEL_GREY, linewidth=0.5, alpha=0.92))


def mask_zeros(xs):
    """replace exact 0.0 values with NaN so matplotlib draws gaps for
    'no measurement' rows instead of a misleading flat-zero baseline"""
    return [float("nan") if v == 0 else v for v in xs]


def fmt_bytes_gib(x, _pos=None):
    """matplotlib FuncFormatter signature is (x, pos). we don't use pos but
    the kw default keeps Pylance from flagging it as unused -- the
    formatter is still called with the positional pos by mpl."""
    del _pos
    if x >= 1024: return f"{x/1024:.1f} T"
    if x >= 1:    return f"{x:.1f} G"
    if x >= 1/1024: return f"{x*1024:.0f} M"
    return f"{x*1024*1024:.0f} K"


# ============================================================ figures

def plot_ingest(rows0, out_dir, caption):
    t = col(rows0, "elapsed_s")
    bytes_w = [r["bytes_written"] / GIB for r in rows0]
    mibs = col(rows0, "write_mibs")
    fig, ax1 = plt.subplots(figsize=(11, 5))
    shade_phases(ax1, t, col(rows0, "phase"))
    ax1.set_xlabel("elapsed (s)")
    ax1.set_ylabel("written (GiB)", color=DEEP_BLUE)
    ax1.plot(t, bytes_w, color=DEEP_BLUE, linewidth=1.6, label="written GiB")
    ax1.tick_params(axis="y", labelcolor=DEEP_BLUE)
    ax2 = ax1.twinx()
    ax2.set_ylabel("write throughput (MiB/s)", color=DEEP_RED)
    ax2.plot(t, mibs, color=DEEP_RED, alpha=0.85, label="MiB/s")
    ax2.tick_params(axis="y", labelcolor=DEEP_RED)
    ax1.set_title("ingest progress and write throughput")
    # cumulative bytes plateau pins blue to the top of the frame and the rate
    # collapses red to the bottom once ingest ends, leaving the mid-right band
    # empty -- park the explainer there so neither line runs through it
    add_explainer(ax1, "cumulative bytes (left, blue)\ninstantaneous rate (right, red)",
                  xy=(0.99, 0.5), va="center")
    finalize(fig, out_dir, "ingest", caption)


def plot_latency(rows0, out_dir, caption):
    """pooled across CFs. log y, SLO reference lines, p50/95/99 colors."""
    t = col(rows0, "elapsed_s")
    phase = col(rows0, "phase")
    fig, axes = plt.subplots(3, 1, figsize=(11, 9), sharex=True)
    for i, (ax, op) in enumerate(zip(axes, OP_TYPES)):
        shade_phases(ax, t, phase)
        ax.plot(t, col(rows0, f"{op}_p50_us"), label="p50",
                color=PCT_COLORS["p50"])
        ax.plot(t, col(rows0, f"{op}_p95_us"), label="p95",
                color=PCT_COLORS["p95"])
        ax.plot(t, col(rows0, f"{op}_p99_us"), label="p99",
                color=PCT_COLORS["p99"])
        ax.set_yscale("log")
        ax.set_ylabel(f"{op} (µs, log)")
        ax.legend(loc="upper right")
        add_slo_lines(ax, label_only=(i == 0))
    axes[-1].set_xlabel("elapsed (s)")
    axes[0].set_title("per-op read latency (pooled across all CFs)")
    add_explainer(axes[0], "p50/p95/p99 of probe latency.\n"
                           "dotted = common SLO refs.",
                  xy=(0.01, 0.97), ha="left")
    finalize(fig, out_dir, "latency", caption, top=0.97)


def plot_read_ops(rows0, out_dir, caption):
    """three per-op lines plus a `total` aggregate. since readers pick op
    uniformly the three per-op lines often coincide; the total trace makes
    the aggregate read rate visible at a glance, and per-op lines use
    different dash patterns so they remain distinguishable when overlapping"""
    t = col(rows0, "elapsed_s")
    pt = col(rows0, "point_ops_s")
    sk = col(rows0, "seek_ops_s")
    rg = col(rows0, "range_ops_s")
    tot = [a + b + c for a, b, c in zip(pt, sk, rg)]

    fig, ax = plt.subplots(figsize=(11, 5))
    shade_phases(ax, t, col(rows0, "phase"))
    # total first so the per-op lines draw on top
    ax.plot(t, tot, label="total", color=INK, linewidth=2.0,
            linestyle="-", alpha=0.9)
    op_styles = {"point": "-", "seek": (0, (4, 2)), "range": (0, (1, 2))}
    for op, ys in zip(OP_TYPES, [pt, sk, rg]):
        ax.plot(t, ys, label=op, color=OP_COLORS[op],
                linestyle=op_styles[op], linewidth=1.4, alpha=0.95)
    ax.set_yscale("log")
    ax.set_xlabel("elapsed (s)")
    ax.set_ylabel("ops / s (log)")
    ax.set_title("read probe throughput (pooled across CFs)")
    ax.legend(loc="lower right", ncol=4, fontsize=8)
    add_explainer(
        ax,
        "uniform op selection -> per-op lines roughly coincide.\n"
        "point = get  ·  seek = iter.seek  ·  range = scan N keys",
        xy=(0.01, 0.97), ha="left")
    finalize(fig, out_dir, "read_ops", caption)


def plot_db_stats(rows0, meta, out_dir, caption):
    """3 rows: storage state, memtable + worker queues, and engine health
    (memory pressure / pending flushes / unified WAL generation)"""
    unified = is_unified_mode(meta)
    mode_tag = "unified buffer" if unified else "Σ per-CF"
    t = col(rows0, "elapsed_s")
    phase = col(rows0, "phase")
    fig, (ax1, ax2, ax3) = plt.subplots(3, 1, figsize=(11, 10.5), sharex=True)
    for ax in (ax1, ax2, ax3):
        shade_phases(ax, t, phase)

    # row 1: sstable_count + immutable_count + open file handles
    ax1.plot(t, col(rows0, "sstable_count"), color=DEEP_BLUE,
             label="sstables (db-wide)")
    ax1.plot(t, col(rows0, "open_sstables"), color=DEEP_CYAN,
             linestyle="--", linewidth=1.0, alpha=0.85,
             label="open sstable handles")
    ax1.set_ylabel("sstables", color=DEEP_BLUE)
    ax1.tick_params(axis="y", labelcolor=DEEP_BLUE)
    ax1.legend(loc="lower right", fontsize=8)
    ax1r = ax1.twinx()
    ax1r.plot(t, col(rows0, "immutable_count"), color=DEEP_PEACH,
              alpha=0.85, label="immutable memtables")
    ax1r.set_ylabel(f"immutable memtables ({mode_tag})", color=DEEP_PEACH)
    ax1r.tick_params(axis="y", labelcolor=DEEP_PEACH)
    ax1.set_title("tidesdb storage state")
    add_explainer(ax1, "sstables = on-disk runs (blue solid)\n"
                       "open handles = file descriptors (cyan dashed)\n"
                       f"immutables = memtables being flushed (orange, {mode_tag})",
                  xy=(0.01, 0.97), ha="left")

    # row 2: memtable bytes + queue depths
    ax2.plot(t, [r["memtable_bytes"] / MIB for r in rows0],
             color=DEEP_MINT, label="memtable (MiB)")
    ax2.set_ylabel(f"memtable (MiB, {mode_tag})", color=DEEP_MINT)
    ax2.tick_params(axis="y", labelcolor=DEEP_MINT)
    ax2r = ax2.twinx()
    ax2r.plot(t, col(rows0, "flush_qsize"), color=DEEP_PURPLE,
              alpha=0.85, label="flush_q")
    ax2r.plot(t, col(rows0, "compact_qsize"), color=DEEP_BROWN,
              alpha=0.85, label="compact_q")
    ax2r.set_ylabel("queue depth")
    ax2r.legend(loc="upper right", fontsize=8)
    add_explainer(ax2,
                  f"memtable bytes (green, {mode_tag})\n"
                  "worker queues (right axis): rising = backlog",
                  xy=(0.01, 0.97), ha="left")

    # row 3: engine health -- memory pressure (0..3) + pending flushes,
    # plus the unified WAL generation when in unified mode
    ax3.plot(t, col(rows0, "flush_pending_count"), color=DEEP_PURPLE,
             label="flush_pending (queued+in-flight)")
    ax3.plot(t, col(rows0, "flush_qsize"), color=DEEP_BROWN, linestyle="--",
             linewidth=1.0, alpha=0.8, label="flush_q (queued only)")
    ax3.set_ylabel("pending flushes", color=DEEP_PURPLE)
    ax3.tick_params(axis="y", labelcolor=DEEP_PURPLE)
    # flush_pending spikes during ingest sit at the upper left, so park this
    # legend along the top-center where the post-ingest trace is flat near zero
    ax3.legend(loc="upper center", fontsize=8)
    ax3r = ax3.twinx()
    ax3r.step(t, col(rows0, "memory_pressure"), color=DEEP_RED, where="post",
              linewidth=1.4, label="memory_pressure")
    if unified:
        ax3r.plot(t, col(rows0, "unified_wal_generation"), color=DEEP_GREY,
                  linestyle=":", linewidth=1.0, alpha=0.8,
                  label="unified WAL gen")
    ax3r.set_ylim(-0.2, 3.4)
    ax3r.set_ylabel("memory_pressure (0-3)", color=DEEP_RED)
    ax3r.tick_params(axis="y", labelcolor=DEEP_RED)
    ax3r.legend(loc="upper right", fontsize=8)
    ax3.set_xlabel("elapsed (s)")
    # legends sit along the top, flush_pending spikes at the upper left, and
    # the memory_pressure step rides the bottom -- the empty band is mid-right,
    # so the explainer goes there to clear all three (add_explainer anchors it
    # on the topmost twin so the red ax3r line can't draw over it)
    add_explainer(ax3,
                  "flush_pending - flush_q = flushes actively writing.\n"
                  "memory_pressure (red step): 0 normal -> 3 critical.",
                  xy=(0.99, 0.5), va="center")

    finalize(fig, out_dir, "db_stats", caption, top=0.96)


LEVEL_CAP = 16  # match the C side's lvl0_ssts..lvl15_ssts schema

def detect_active_levels(by_cf, suffix="ssts"):
    """return the subset of [0..LEVEL_CAP) that has any non-zero value in any
    CF across the whole run for the given per-level field suffix (ssts /
    bytes / keys / tomb). levels that never filled are excluded so the
    stacked area + legend don't show empty L11..L15 padding"""
    active = []
    for j in range(LEVEL_CAP):
        for rs in by_cf.values():
            if any(r.get(f"lvl{j}_{suffix}", 0) > 0 for r in rs):
                active.append(j); break
    return active or [0]  # at least L0 even if no flushes happened


def plot_per_cf_levels(by_cf, meta, out_dir, caption):
    """one subplot per CF. left axis: stacked level counts (only levels that
    saw activity in this run -- tidesdb grows levels organically up to 32).
    right axis: cf_memtable_size so CFs that haven't flushed yet still show
    activity (full memtable, empty disk).

    in unified-memtable mode tidesdb collapses every CF's active memtable
    into one shared buffer -- cf_memtable_size reads 0 for every CF -- so
    we drop the per-CF overlay and skip the right axis entirely.

    a single shared legend lives in the figure margin so it doesn't eat
    subplot area"""
    unified = is_unified_mode(meta)
    cfs = sorted(by_cf.keys())
    n = len(cfs)
    cols_n = min(n, 2)
    rows_n = (n + cols_n - 1) // cols_n
    fig, axes = plt.subplots(rows_n, cols_n, figsize=(11, 3.0 * rows_n + 0.6),
                             sharex=True, squeeze=False)

    active_levels = detect_active_levels(by_cf)
    stack_handles = None
    for i, cf in enumerate(cfs):
        ax = axes[i // cols_n][i % cols_n]
        rs = by_cf[cf]
        t = col(rs, "elapsed_s")
        levels = [col(rs, f"lvl{j}_ssts") for j in active_levels]
        shade_phases(ax, t, col(rs, "phase"))
        polys = ax.stackplot(t, *levels,
                             labels=[f"L{j}" for j in active_levels],
                             alpha=0.85)
        if stack_handles is None: stack_handles = polys
        ax.set_title(f"cf_{cf}", fontsize=10)
        ax.set_ylabel("sstables")
        if not unified:
            axr = ax.twinx()
            mtab = [r["cf_memtable_size"] / MIB for r in rs]
            axr.plot(t, mtab, color=RULE, linewidth=1.0,
                     linestyle="--", alpha=0.85, label="memtable")
            axr.set_ylabel("memtable (MiB)", color=RULE, fontsize=9)
            axr.tick_params(axis="y", labelcolor=RULE, labelsize=8)
    for i in range(n, rows_n * cols_n):
        axes[i // cols_n][i % cols_n].set_visible(False)
    for ax in axes[-1]: ax.set_xlabel("elapsed (s)")
    title = "per-level sstable distribution, per CF"
    if unified:
        title += "  (unified memtable: per-CF memtable overlay omitted)"
    fig.suptitle(title, fontsize=12)
    if stack_handles is not None:
        labels = [f"L{j}" for j in active_levels]
        handles = list(stack_handles)
        if not unified:
            labels.append("memtable (dashed)")
            handles.append(plt.Line2D([0], [0], color=RULE, linestyle="--",
                                      linewidth=1.0))
        ncol = min(len(labels), 8)
        fig.legend(handles, labels, loc="upper right", ncol=ncol,
                   fontsize=8, frameon=False, bbox_to_anchor=(0.99, 0.97))
    finalize(fig, out_dir, "per_cf_levels", caption, top=0.93)


def plot_per_cf_level_bytes(by_cf, out_dir, caption):
    """companion to plot_per_cf_levels: stacked per-level BYTES (klog+vlog)
    per CF. the count view shows how many runs sit at each level; this shows
    where the data volume actually lives, which is the figure that reveals
    the LSM's size pyramid (each level ~level_size_ratio larger than the one
    above). y-axis auto-scales via the GiB formatter."""
    cfs = sorted(by_cf.keys())
    n = len(cfs)
    cols_n = min(n, 2)
    rows_n = (n + cols_n - 1) // cols_n
    fig, axes = plt.subplots(rows_n, cols_n, figsize=(11, 3.0 * rows_n + 0.6),
                             sharex=True, squeeze=False)

    active_levels = detect_active_levels(by_cf, "bytes")
    stack_handles = None
    for i, cf in enumerate(cfs):
        ax = axes[i // cols_n][i % cols_n]
        rs = by_cf[cf]
        t = col(rs, "elapsed_s")
        # stack in GiB so the shared formatter reads naturally
        levels = [[r.get(f"lvl{j}_bytes", 0) / GIB for r in rs]
                  for j in active_levels]
        shade_phases(ax, t, col(rs, "phase"))
        polys = ax.stackplot(t, *levels,
                             labels=[f"L{j}" for j in active_levels],
                             alpha=0.85)
        if stack_handles is None: stack_handles = polys
        ax.set_title(f"cf_{cf}", fontsize=10)
        ax.set_ylabel("level size")
        ax.yaxis.set_major_formatter(FuncFormatter(fmt_bytes_gib))
    for i in range(n, rows_n * cols_n):
        axes[i // cols_n][i % cols_n].set_visible(False)
    for ax in axes[-1]: ax.set_xlabel("elapsed (s)")
    fig.suptitle("per-level data size distribution, per CF", fontsize=12)
    if stack_handles is not None:
        labels = [f"L{j}" for j in active_levels]
        ncol = min(len(labels), 8)
        fig.legend(list(stack_handles), labels, loc="upper right", ncol=ncol,
                   fontsize=8, frameon=False, bbox_to_anchor=(0.99, 0.97))
    finalize(fig, out_dir, "per_cf_level_bytes", caption, top=0.93)


def plot_per_cf_keys(by_cf, out_dir, caption):
    cfs = sorted(by_cf.keys())
    fig, ax = plt.subplots(figsize=(11, 5))
    rs0 = by_cf[cfs[0]]
    shade_phases(ax, col(rs0, "elapsed_s"), col(rs0, "phase"))
    for cf in cfs:
        rs = by_cf[cf]
        ax.plot(col(rs, "elapsed_s"), col(rs, "cf_total_keys"),
                label=f"cf_{cf}", linewidth=1.5)
    ax.set_xlabel("elapsed (s)")
    ax.set_ylabel("live keys (cf_total_keys)")
    ax.set_title("per-CF live key count over time")
    ax.legend(loc="lower right")
    # legend sits lower-right; keep the explainer upper-left where the early
    # ramp leaves whitespace, so the two never overlap
    add_explainer(ax,
                  "writers fan out every commit to all CFs,\n"
                  "so curves overlap until the delete phase\n"
                  "deducts from the target CF(s).",
                  xy=(0.01, 0.97), ha="left")
    finalize(fig, out_dir, "per_cf_keys", caption)


def plot_per_cf_tombstones(by_cf, out_dir, caption):
    cfs = sorted(by_cf.keys())
    fig, axes = plt.subplots(2, 1, figsize=(11, 8), sharex=True)
    rs0 = by_cf[cfs[0]]
    for ax in axes: shade_phases(ax, col(rs0, "elapsed_s"), col(rs0, "phase"))
    for cf in cfs:
        rs = by_cf[cf]
        axes[0].plot(col(rs, "elapsed_s"),
                     [r["cf_tombstone_ratio"] * 100 for r in rs],
                     label=f"cf_{cf}", linewidth=1.4)
        axes[1].plot(col(rs, "elapsed_s"),
                     [r["cf_max_density"] * 100 for r in rs],
                     label=f"cf_{cf}", linewidth=1.4)
    axes[0].set_ylabel("tombstone ratio (%)")
    axes[1].set_ylabel("max per-sst tombstone density (%)")
    axes[1].set_xlabel("elapsed (s)")
    axes[0].legend(loc="upper right")
    axes[0].set_title("per-CF tombstone observability")
    add_explainer(axes[0], "tombstones / total keys.\n"
                            "rises on delete, drops to 0 after\n"
                            "compaction merges them out.",
                  xy=(0.01, 0.97), ha="left")
    add_explainer(axes[1], "worst single-sst tombstone density.\n"
                            "drives tombstone_density_trigger.",
                  xy=(0.01, 0.97), ha="left")
    finalize(fig, out_dir, "per_cf_tombstones", caption, top=0.96)


def plot_per_cf_latency(by_cf, out_dir, caption):
    cfs = sorted(by_cf.keys())
    fig, axes = plt.subplots(3, 1, figsize=(11, 9), sharex=True)
    rs0 = by_cf[cfs[0]]
    for ax in axes: shade_phases(ax, col(rs0, "elapsed_s"), col(rs0, "phase"))
    for i, (ax, op) in enumerate(zip(axes, OP_TYPES)):
        for cf in cfs:
            rs = by_cf[cf]
            ax.plot(col(rs, "elapsed_s"),
                    mask_zeros(col(rs, f"cf_{op}_p99_us")),
                    label=f"cf_{cf}", alpha=0.9, linewidth=1.3)
        ax.set_yscale("log")
        ax.set_ylabel(f"{op} p99 (µs, log)")
        add_slo_lines(ax, label_only=(i == 0))
        if i == 0: ax.legend(loc="upper right", ncol=min(len(cfs), 4))
    axes[-1].set_xlabel("elapsed (s)")
    axes[0].set_title("per-CF read p99 latency")
    add_explainer(axes[0], "compare CFs with different configs\n"
                           "(e.g. btree vs skip-list, compression).\n"
                           "lines drop where no probes hit that CF.",
                  xy=(0.99, 0.04), va="bottom")
    finalize(fig, out_dir, "per_cf_latency", caption, top=0.97)


def plot_per_cf_amp(by_cf, out_dir, caption):
    """per-CF read amplification + on-disk data size over time. read_amp is
    tidesdb's point-lookup cost multiplier (roughly the bloom-adjusted number
    of levels a get must consult); it should fall as compaction merges runs.
    cf_data_size is each CF's klog+vlog bytes -- in mirror mode the CFs track
    each other; configs that differ (compression, btree) diverge here."""
    cfs = sorted(by_cf.keys())
    fig, axes = plt.subplots(2, 1, figsize=(11, 8), sharex=True)
    rs0 = by_cf[cfs[0]]
    for ax in axes:
        shade_phases(ax, col(rs0, "elapsed_s"), col(rs0, "phase"))
    for cf in cfs:
        rs = by_cf[cf]
        axes[0].plot(col(rs, "elapsed_s"), col(rs, "cf_read_amp"),
                     label=f"cf_{cf}", linewidth=1.4)
        axes[1].plot(col(rs, "elapsed_s"),
                     [r.get("cf_data_size_bytes", 0) / GIB for r in rs],
                     label=f"cf_{cf}", linewidth=1.4)
    axes[0].set_ylabel("read amp (× levels probed)")
    axes[1].set_ylabel("data size")
    axes[1].yaxis.set_major_formatter(FuncFormatter(fmt_bytes_gib))
    axes[1].set_xlabel("elapsed (s)")
    axes[0].legend(loc="upper right", ncol=min(len(cfs), 4))
    axes[0].set_title("per-CF read amplification and data size")
    add_explainer(axes[0], "lower read_amp = fewer runs per lookup.\n"
                           "spikes during ingest, settles after compaction.",
                  xy=(0.01, 0.97), ha="left")
    add_explainer(axes[1], "per-CF klog+vlog bytes (db-reported).\n"
                           "diverges when CFs differ in compression / format.",
                  xy=(0.01, 0.97), ha="left")
    finalize(fig, out_dir, "per_cf_amp", caption, top=0.96)


def plot_per_cf_btree(by_cf, out_dir, caption):
    """b+tree klog structure, per CF -- only emitted when at least one CF was
    built with use_btree=1 (otherwise every value is 0 and the figure is
    noise). tree height bounds the per-sstable seek cost; node count tracks
    total index size."""
    cfs = sorted(by_cf.keys())
    has_btree = any(
        any(r.get("cf_btree_total_nodes", 0) > 0
            or r.get("cf_btree_max_height", 0) > 0 for r in by_cf[cf])
        for cf in cfs)
    if not has_btree:
        return
    fig, axes = plt.subplots(2, 1, figsize=(11, 8), sharex=True)
    rs0 = by_cf[cfs[0]]
    for ax in axes:
        shade_phases(ax, col(rs0, "elapsed_s"), col(rs0, "phase"))
    for cf in cfs:
        rs = by_cf[cf]
        # skip CFs that never built a btree (skip-list CFs in a mixed run)
        if not any(r.get("cf_btree_total_nodes", 0) > 0 for r in rs):
            continue
        axes[0].plot(col(rs, "elapsed_s"), col(rs, "cf_btree_avg_height"),
                     label=f"cf_{cf} avg", linewidth=1.4)
        axes[0].plot(col(rs, "elapsed_s"), col(rs, "cf_btree_max_height"),
                     label=f"cf_{cf} max", linewidth=1.0, linestyle="--",
                     alpha=0.8)
        axes[1].plot(col(rs, "elapsed_s"), col(rs, "cf_btree_total_nodes"),
                     label=f"cf_{cf}", linewidth=1.4)
    axes[0].set_ylabel("b+tree height (levels)")
    axes[1].set_ylabel("total b+tree nodes")
    axes[1].set_xlabel("elapsed (s)")
    axes[0].legend(loc="upper right", ncol=min(len(cfs), 4), fontsize=8)
    axes[0].set_title("per-CF b+tree klog structure")
    add_explainer(axes[0], "avg (solid) / max (dashed) tree height across\n"
                           "sstables. height bounds per-sst seek cost.",
                  xy=(0.01, 0.97), ha="left")
    add_explainer(axes[1], "total index nodes across all sstables.\n"
                           "grows with data, drops as compaction merges.",
                  xy=(0.01, 0.97), ha="left")
    finalize(fig, out_dir, "per_cf_btree", caption, top=0.96)


def plot_cache(rows0, out_dir, caption):
    """2 rows: hit rate (cumulative vs windowed) over cumulative hit/miss
    counts, and cache occupancy (bytes used + entry count)."""
    t = col(rows0, "elapsed_s")
    fig, (ax, axf) = plt.subplots(2, 1, figsize=(11, 8), sharex=True)
    shade_phases(ax, t, col(rows0, "phase"))
    shade_phases(axf, t, col(rows0, "phase"))

    # row 1: hit rate (cumulative solid, windowed dashed) + cumulative ops
    ln1, = ax.plot(t, [r["cache_hit_rate"] * 100 for r in rows0],
                   color=DEEP_BLUE, label="hit rate cum (%)")
    ln_w, = ax.plot(t, [r.get("cache_hit_rate_window", 0) * 100 for r in rows0],
                    color=DEEP_CYAN, linestyle="--", linewidth=1.1, alpha=0.9,
                    label="hit rate window (%)")
    ax.set_ylim(0, 100)
    ax.set_ylabel("cache hit rate (%)", color=DEEP_BLUE)
    ax2 = ax.twinx()
    ln2, = ax2.plot(t, col(rows0, "cache_hits"), color=DEEP_MINT,
                    label="hits (cum)", alpha=0.7)
    ln3, = ax2.plot(t, col(rows0, "cache_misses"), color=DEEP_RED,
                    label="misses (cum)", alpha=0.7)
    ax2.set_ylabel("cumulative ops")
    ax.set_title("block clock-cache hit rate and occupancy")
    ax.legend(handles=[ln1, ln_w, ln2, ln3], loc="center right", fontsize=8)
    add_explainer(ax, "cumulative (blue) vs this-window (cyan dashed) hit rate.\n"
                      "windowed reveals recent behaviour the cumulative masks.",
                  xy=(0.99, 0.04), va="bottom")

    # row 2: occupancy -- bytes used (left) + entry count (right)
    axf.plot(t, [r.get("cache_total_bytes", 0) / MIB for r in rows0],
             color=DEEP_BLUE, label="cache used (MiB)")
    axf.set_ylabel("cache used (MiB)", color=DEEP_BLUE)
    axf.tick_params(axis="y", labelcolor=DEEP_BLUE)
    axf.set_xlabel("elapsed (s)")
    axfr = axf.twinx()
    axfr.plot(t, col(rows0, "cache_total_entries"), color=DEEP_PEACH,
              alpha=0.85, label="entries")
    axfr.set_ylabel("cached entries", color=DEEP_PEACH)
    axfr.tick_params(axis="y", labelcolor=DEEP_PEACH)
    add_explainer(axf, "resident cache size (blue) and entry count (orange).\n"
                       "flat ceiling = cache saturated / evicting.",
                  xy=(0.01, 0.97), ha="left")

    finalize(fig, out_dir, "cache", caption, top=0.96)


def plot_disk(rows0, out_dir, caption):
    t = col(rows0, "elapsed_s")
    fig, ax = plt.subplots(figsize=(11, 5))
    shade_phases(ax, t, col(rows0, "phase"))
    ax.plot(t, [r["disk_bytes"] / GIB for r in rows0], label="disk (du)")
    ax.plot(t, [r["data_size_bytes"] / GIB for r in rows0],
            label="data (db-stats)")
    ax.plot(t, [r["memtable_bytes"] / GIB for r in rows0],
            label="memtable", alpha=0.7)
    ax.set_xlabel("elapsed (s)")
    ax.set_ylabel("size (GiB)")
    ax.yaxis.set_major_formatter(FuncFormatter(fmt_bytes_gib))
    ax.set_title("disk vs db-reported data size vs memtable")
    ax.legend(loc="upper left")
    add_explainer(ax, "disk = filesystem du walk\n"
                       "data = sum of klog+vlog from tidesdb\n"
                       "memtable = in-RAM active buffer",
                  xy=(0.99, 0.04), va="bottom")
    finalize(fig, out_dir, "disk", caption)


def plot_integrity(rows0, out_dir, caption):
    t = col(rows0, "elapsed_s")
    fig, ax = plt.subplots(figsize=(11, 5))
    shade_phases(ax, t, col(rows0, "phase"))
    for op in OP_TYPES:
        ax.plot(t, col(rows0, f"{op}_misses"),
                label=f"{op} unexpected misses (cum)", color=OP_COLORS[op])
    ax.plot(t, col(rows0, "mismatches"),
            label="value mismatches (cum)", linewidth=1.8, color=INK)
    ax.set_xlabel("elapsed (s)")
    ax.set_ylabel("cumulative count")
    ax.set_title("data integrity")
    ax.legend(loc="upper left")
    add_explainer(ax, "any non-zero = real data loss / corruption.\n"
                       "deliberate deletes are excluded from misses.",
                  xy=(0.99, 0.04), va="bottom")
    finalize(fig, out_dir, "integrity", caption)


def plot_p99_vs_size(rows0, out_dir, caption):
    """phase-shaded along the byte axis since x is not time"""
    bytes_g = [r["bytes_written"] / GIB for r in rows0]
    phase = col(rows0, "phase")
    fig, ax = plt.subplots(figsize=(11, 5))
    shade_phases_x(ax, bytes_g, phase)
    for op in OP_TYPES:
        ax.plot(bytes_g, col(rows0, f"{op}_p99_us"), label=f"{op} p99",
                color=OP_COLORS[op])
    ax.set_yscale("log")
    ax.set_xlabel("data written (GiB)")
    ax.set_ylabel("p99 latency (µs, log)")
    ax.set_title("read p99 latency vs ingested size")
    # SLO labels live at the right edge, explainer at bottom-right, so the
    # legend goes upper-left to stay clear of both
    ax.legend(loc="upper left")
    add_slo_lines(ax, label_only=True)
    add_explainer(ax, "shading = phase when each sample was taken.\n"
                       "left-to-right = early run >> late run.",
                  xy=(0.99, 0.04), va="bottom")
    finalize(fig, out_dir, "p99_vs_size", caption)


def plot_amplification(amp_rows, out_dir, caption):
    """three panels. each flow-amp panel draws the whole-run cumulative ratio
    bold -- the headline figure the run converges to -- with the per-sample-
    window ratios kept faint behind it for the transient detail (compaction
    bursts) the cumulative curve smooths away.

    write amp cumulative is recovered here by dividing the emitted *_total
    counters. read amp cumulative is read straight from ra_app_cum / ra_io_cum:
    the file can't reconstruct it because rd_logical (the denominator) is only
    ever a per-window quantity, so the bench accumulates it run-wide. space amp
    is already a whole-DB state ratio (on-disk size / live logical at the
    instant), so it has no separate cumulative form -- the single line is it."""
    if not amp_rows:
        return
    t = col(amp_rows, "elapsed_s")
    phase = col(amp_rows, "phase")
    fig, axes = plt.subplots(3, 1, figsize=(11, 10), sharex=True)
    for ax in axes: shade_phases(ax, t, phase)

    def ratio(num_col, den_col):
        """whole-run cumulative ratio of two cumulative counters; NaN until the
        denominator starts moving so the line begins at the first real datum"""
        return [n / d if d > 0 else float("nan")
                for n, d in zip(col(amp_rows, num_col), col(amp_rows, den_col))]

    def read_amp_cum(cum_col, num_window_col):
        """running cumulative read amp. newer runs carry the bench's exact
        ra_*_cum columns -- use them directly. runs predating those columns
        don't, but the cumulative is still recoverable from the per-window
        data: rd_logical (the denominator) isn't emitted per row, yet in any
        window where the windowed ratio is defined it equals
        rchar_window / ra_app, so the bytes-weighted running ratio is
        sum(num_window) / sum(rd_logical) over exactly the windows that passed
        the bench's write-dominated gate (ra_app > 0)."""
        cum = col(amp_rows, cum_col)
        if any(v > 0 for v in cum):
            return mask_zeros(cum)
        num_w  = col(amp_rows, num_window_col)
        rchar  = col(amp_rows, "rchar_window")
        ra_app = col(amp_rows, "ra_app")
        out, snum, sden = [], 0.0, 0.0
        for nw, rw, raa in zip(num_w, rchar, ra_app):
            if raa > 0 and rw > 0:
                snum += nw
                sden += rw / raa
            out.append(snum / sden if sden > 0 else float("nan"))
        return out

    FAINT = dict(linewidth=0.8, alpha=0.30)  # per-window context behind the cum
    BOLD = 2.0

    #  panel 1: write amp
    ax = axes[0]
    user_bw = col(amp_rows, "user_bytes_window")
    def mask_when_no_writes(xs):
        return [float("nan") if b == 0 else v for v, b in zip(xs, user_bw)]
    # faint per-window first so the cumulative lines draw on top
    ax.plot(t, mask_when_no_writes(col(amp_rows, "wa_app")),  color=DEEP_BLUE,  **FAINT)
    ax.plot(t, mask_when_no_writes(col(amp_rows, "wa_io")),   color=DEEP_PEACH, **FAINT)
    ax.plot(t, mask_when_no_writes(col(amp_rows, "wa_disk")), color=DEEP_MINT,  **FAINT)
    # bold whole-run cumulative (the headline)
    ax.plot(t, ratio("wchar_total", "user_bytes_total"),
            label="wa_app  (syscalls)", color=DEEP_BLUE, linewidth=BOLD)
    ax.plot(t, ratio("write_bytes_total", "user_bytes_total"),
            label="wa_io  (block layer)", color=DEEP_PEACH, linewidth=BOLD)
    ax.plot(t, ratio("disk_total_bytes", "user_bytes_total"),
            label="wa_disk  (du / ingested)", color=DEEP_MINT, linewidth=BOLD)
    ax.axhline(1.0, color=RULE, linewidth=0.6, linestyle=":", zorder=0)
    ax.set_ylabel("write amp (× bytes / ingested)")
    ax.set_title("amplification metrics  (bold = whole-run cumulative)")
    ax.legend(loc="upper right", fontsize=8, ncol=3)
    add_explainer(ax,
                  "write amp = bytes written / bytes ingested.\n"
                  "bold = cumulative over the run; faint = per-window.\n"
                  "wa > 1 >> LSM rewrote data (compaction).",
                  xy=(0.01, 0.97), ha="left")

    #  panel 2: read amp
    ax = axes[1]
    ax.plot(t, mask_zeros(col(amp_rows, "ra_app")), color=DEEP_BLUE,  **FAINT)
    ax.plot(t, mask_zeros(col(amp_rows, "ra_io")),  color=DEEP_PEACH, **FAINT)
    ax.plot(t, read_amp_cum("ra_app_cum", "rchar_window"),
            label="ra_app  (syscalls)", color=DEEP_BLUE, linewidth=BOLD)
    ax.plot(t, read_amp_cum("ra_io_cum", "read_bytes_window"),
            label="ra_io  (block layer)", color=DEEP_PEACH, linewidth=BOLD)
    ax.axhline(1.0, color=RULE, linewidth=0.6, linestyle=":", zorder=0)
    ax.set_ylabel("read amp (× bytes / returned)")
    ax.legend(loc="upper right", fontsize=8)
    add_explainer(ax,
                  "read amp = kernel-read bytes / bytes the reader got back.\n"
                  "bold = cumulative; faint = per-window (spikes in compaction).",
                  xy=(0.01, 0.97), ha="left")

    #  panel 3: space amp
    ax = axes[2]
    ax.plot(t, col(amp_rows, "sa_disk"),
            label="sa_disk  (du / live logical)", color=DEEP_BLUE)
    ax.plot(t, col(amp_rows, "sa_data"),
            label="sa_data  (db-stats / live)", color=DEEP_MINT)
    ax.axhline(1.0, color=RULE, linewidth=0.6, linestyle=":", zorder=0)
    ax.set_ylabel("space amp (× on-disk / live)")
    ax.set_xlabel("elapsed (s)")
    ax.legend(loc="upper right", fontsize=8)
    add_explainer(ax,
                  "space amp = on-disk size / live logical bytes.\n"
                  "in mirror mode > 1 reflects the num_cfs multiplier.",
                  xy=(0.01, 0.97), ha="left")

    finalize(fig, out_dir, "amplification", caption, top=0.97)


def plot_power(rows0, out_dir, caption):
    """three panels:
       1) watts per RAPL domain with cpu_share filled in the background
       2) cumulative attributed energy
       3) energy efficiency -- J/GiB written (left) and J/M reads (right)

    skipped entirely when RAPL is unavailable (every power sample is 0), so
    machines without powercap don't get an empty all-zero figure.
    """
    if not rows0:
        return
    if not any(r.get("power_pkg_w", 0) > 0 for r in rows0):
        return
    t = col(rows0, "elapsed_s")
    fig, axes = plt.subplots(3, 1, figsize=(11, 10), sharex=True)
    for ax in axes: shade_phases(ax, t, col(rows0, "phase"))

    # panel 1: power + cpu_share background
    ax = axes[0]
    cpu = col(rows0, "cpu_share")
    ax_bg = ax.twinx()
    ax_bg.fill_between(t, 0, cpu, color=RULE, alpha=0.18, zorder=0,
                       label="cpu_share")
    ax_bg.set_ylim(0, 1)
    ax_bg.set_ylabel("cpu_share (fill)", color=RULE)
    ax_bg.tick_params(axis="y", labelcolor=RULE)
    ax.plot(t, col(rows0, "power_pkg_w"),    label="pkg",    color=DEEP_RED,
            linewidth=1.6)
    ax.plot(t, col(rows0, "power_core_w"),   label="core",   color=DEEP_PEACH,
            linewidth=1.2)
    ax.plot(t, col(rows0, "power_uncore_w"), label="uncore", color=DEEP_PURPLE,
            linewidth=1.0, alpha=0.85)
    if any(col(rows0, "power_dram_w")):
        ax.plot(t, col(rows0, "power_dram_w"), label="dram",
                color=DEEP_CYAN, linewidth=1.0, alpha=0.85)
    ax.set_ylabel("power (W)")
    ax.set_title("RAPL package power and attributed energy")
    ax.legend(loc="upper right", fontsize=8, ncol=4)
    add_explainer(ax,
                  "lines = whole-package watts (RAPL).\n"
                  "grey fill = cpu_share = bench's slice of all-core time.",
                  xy=(0.01, 0.97), ha="left")

    # panel 2: cumulative attributed energy
    ax = axes[1]
    ax.plot(t, col(rows0, "proc_pkg_energy_j_cum"), color=DEEP_BLUE,
            label="process pkg energy (J, cum)")
    ax.set_ylabel("attributed energy (J)")
    ax.legend(loc="upper left", fontsize=8)
    add_explainer(ax,
                  "= Σ (pkg_W · Δt · cpu_share).\n"
                  "approximates the energy this bench burned on its slice "
                  "of the CPU.",
                  xy=(0.99, 0.04), va="bottom")

    # panel 3: efficiency, separate y-axes for distinct units
    ax = axes[2]
    eff_w = mask_zeros(col(rows0, "energy_per_gib_written"))
    eff_r = mask_zeros(col(rows0, "energy_per_mread"))
    ln1, = ax.plot(t, eff_w, label="J / GiB written", color=DEEP_RED)
    ax.set_ylabel("J / GiB written", color=DEEP_RED)
    ax.tick_params(axis="y", labelcolor=DEEP_RED)
    axr = ax.twinx()
    ln2, = axr.plot(t, eff_r, label="J / M read ops", color=DEEP_PURPLE,
                    alpha=0.85)
    axr.set_ylabel("J / M read ops", color=DEEP_PURPLE)
    axr.tick_params(axis="y", labelcolor=DEEP_PURPLE)
    ax.set_xlabel("elapsed (s)")
    ax.legend(handles=[ln1, ln2], loc="upper right", fontsize=8)
    add_explainer(ax,
                  "ingest efficiency on the left, read efficiency on the right.\n"
                  "lower = more work per joule.",
                  xy=(0.01, 0.97), ha="left")

    finalize(fig, out_dir, "power", caption, top=0.97)


def plot_commit_latency(rows0, num_cfs, out_dir, caption):
    t = col(rows0, "elapsed_s")
    fig, axes = plt.subplots(2, 1, figsize=(11, 7), sharex=True)
    for ax in axes: shade_phases(ax, t, col(rows0, "phase"))

    ax = axes[0]
    ax.plot(t, mask_zeros(col(rows0, "commit_p50_us")), label="p50",
            color=PCT_COLORS["p50"])
    ax.plot(t, mask_zeros(col(rows0, "commit_p95_us")), label="p95",
            color=PCT_COLORS["p95"])
    ax.plot(t, mask_zeros(col(rows0, "commit_p99_us")), label="p99",
            color=PCT_COLORS["p99"])
    ax.set_yscale("log")
    ax.set_ylabel("commit latency (µs, log)")
    ax.set_title("writer commit latency")
    ax.legend(loc="upper right")
    add_slo_lines(ax, label_only=True)
    explain = "p50/p95/p99 of every writer's txn_commit nanos."
    if num_cfs > 1:
        explain += f"\nmirror mode >> each commit fans out to {num_cfs} CFs."
    add_explainer(ax, explain, xy=(0.01, 0.97), ha="left")

    ax = axes[1]
    ax.plot(t, col(rows0, "commits_in_window"), color=DEEP_BLUE)
    ax.set_ylabel("commits per sample window")
    ax.set_xlabel("elapsed (s)")
    add_explainer(ax,
                  "throughput proxy: commits this window.\n"
                  "drops to 0 outside ingest.",
                  xy=(0.99, 0.97))

    finalize(fig, out_dir, "commit_latency", caption, top=0.96)


def plot_process_resources(rows0, ncpu, out_dir, caption):
    t = col(rows0, "elapsed_s")
    fig, axes = plt.subplots(2, 1, figsize=(11, 7), sharex=True)
    for ax in axes: shade_phases(ax, t, col(rows0, "phase"))

    ax = axes[0]
    rss = [r["ru_maxrss_kb"] / 1024.0 for r in rows0]
    ax.plot(t, rss, color=DEEP_BLUE, label="RSS (MiB, peak)")
    ax.set_ylabel("RSS (MiB)", color=DEEP_BLUE)
    ax.tick_params(axis="y", labelcolor=DEEP_BLUE)
    ax2 = ax.twinx()
    ax2.plot(t, col(rows0, "ru_majflt"), color=DEEP_RED, linewidth=0.9,
             alpha=0.7, label="major faults")
    ax2.set_ylabel("major page faults (cum)", color=DEEP_RED)
    ax2.tick_params(axis="y", labelcolor=DEEP_RED)
    ax.set_title("process resources")
    add_explainer(ax,
                  "peak RSS (blue) - never decreases (ru_maxrss).\n"
                  "major faults (red) - slow path: pages from disk.",
                  xy=(0.01, 0.97), ha="left")

    ax = axes[1]
    pct = []
    for i, r in enumerate(rows0):
        if i == 0:
            pct.append(0.0); continue
        dt = (t[i] - t[i - 1]) * 1e6  # us
        d_cpu = r["ru_utime_us_d"] + r["ru_stime_us_d"]
        pct.append((d_cpu / dt) * 100.0 if dt > 0 else 0.0)
    ax.plot(t, pct, color=DEEP_PURPLE)
    ax.axhline(ncpu * 100, color=RULE, linewidth=0.6, linestyle=":",
               zorder=0)
    ax.text(0.995, ncpu * 100, f"max ({ncpu} cores) ",
            transform=ax.get_yaxis_transform(),
            va="bottom", ha="right", fontsize=8, color=RULE)
    ax.set_ylabel(f"CPU % (max = {ncpu * 100})")
    ax.set_xlabel("elapsed (s)")
    add_explainer(ax,
                  f"100% = one core saturated;\n"
                  f"{ncpu * 100}% = every core saturated.",
                  xy=(0.01, 0.97), ha="left")

    finalize(fig, out_dir, "process_resources", caption, top=0.96)


def plot_objstore(rows0, meta, out_dir, caption):
    """object-store mode telemetry: cache fill + upload queue + cumulative
    uploads. only emits when objstore is actually enabled in this run.
    if objstore_replica_mode==1 in any row we tag the title as 'replica'."""
    if not rows0: return
    if not any(r.get("objstore_enabled", 0) > 0 for r in rows0):
        return
    t = col(rows0, "elapsed_s")
    enabled_any_replica = any(r.get("objstore_replica_mode", 0) > 0 for r in rows0)
    is_replica = enabled_any_replica or (meta and meta.get("role") == "replica")

    fig, axes = plt.subplots(3, 1, figsize=(11, 9), sharex=True)
    for ax in axes: shade_phases(ax, t, col(rows0, "phase"))

    # panel 1: local cache size (used vs max) + cache file count on right axis
    ax = axes[0]
    cache_used = [r["objstore_local_cache_bytes"] / MIB for r in rows0]
    cache_max  = [r["objstore_local_cache_max"]   / MIB for r in rows0]
    ax.plot(t, cache_used, color=DEEP_BLUE, label="cache used (MiB)")
    if any(v > 0 for v in cache_max):
        ax.plot(t, cache_max, color=DEEP_RED, linestyle="--", alpha=0.7,
                label="cache cap (MiB)")
    ax.set_ylabel("local cache (MiB)", color=DEEP_BLUE)
    ax.tick_params(axis="y", labelcolor=DEEP_BLUE)
    axr = ax.twinx()
    axr.plot(t, col(rows0, "objstore_local_cache_files"), color=RULE,
             linewidth=0.9, alpha=0.7, label="cache files")
    axr.set_ylabel("cache files", color=RULE)
    axr.tick_params(axis="y", labelcolor=RULE, labelsize=8)
    role_tag = "(replica)" if is_replica else "(primary)"
    backend = meta.get("objstore_mode", "?") if meta else "?"
    ax.set_title(f"object-store mode telemetry  {role_tag}  backend={backend}")
    ax.legend(loc="upper right", fontsize=8)
    add_explainer(ax,
                  "local cache holds recently-accessed sstable / vlog files.\n"
                  "cap=0 means unlimited; eviction kicks in only when bounded.",
                  xy=(0.01, 0.97), ha="left")

    # panel 2: upload pipeline -- queue depth + cumulative uploads
    ax = axes[1]
    ax.plot(t, col(rows0, "objstore_upload_queue_depth"), color=DEEP_PEACH,
            label="upload queue depth")
    ax.set_ylabel("upload queue depth", color=DEEP_PEACH)
    ax.tick_params(axis="y", labelcolor=DEEP_PEACH)
    axr = ax.twinx()
    axr.plot(t, col(rows0, "objstore_total_uploads"), color=DEEP_MINT,
             label="cumulative uploads")
    axr.plot(t, col(rows0, "objstore_total_upload_failures"), color=DEEP_RED,
             linewidth=1.4, label="upload failures")
    axr.set_ylabel("cumulative count")
    axr.legend(loc="upper right", fontsize=8)
    add_explainer(ax,
                  "queue depth = backlog of pending uploads.\n"
                  "cumulative uploads should climb steadily during writes;\n"
                  "any non-zero failure curve is a red flag.",
                  xy=(0.01, 0.97), ha="left")

    # panel 3: WAL upload generation -- how far behind primary the replica is
    ax = axes[2]
    ax.plot(t, col(rows0, "objstore_last_uploaded_gen"), color=DEEP_PURPLE)
    ax.set_xlabel("elapsed (s)")
    ax.set_ylabel("last uploaded WAL gen")
    add_explainer(ax,
                  "monotonic counter -- a primary advances it as it ships\n"
                  "WAL segments; a replica observes it via MANIFEST poll.",
                  xy=(0.01, 0.97), ha="left")

    finalize(fig, out_dir, "objstore", caption, top=0.96)


def plot_delete_experiment(out_dir, delete_path, caption):
    if not os.path.exists(delete_path):
        return
    drows = load_csv(delete_path)
    if not drows:
        return
    labels = [r["phase"] for r in drows]
    disk = [r["disk_bytes"] / GIB for r in drows]
    data = [r["data_bytes"] / GIB for r in drows]
    ssts = [r["sstable_count"] for r in drows]
    tomb = [r.get("tombstones", 0) for r in drows]
    total_del = drows[-1].get("total_deleted", 0)

    fig, axes = plt.subplots(1, 3, figsize=(15, 5))
    x = list(range(len(labels)))
    w = 0.35

    axes[0].bar([i - w / 2 for i in x], disk, width=w, label="disk (du)",
                color=DEEP_BLUE)
    axes[0].bar([i + w / 2 for i in x], data, width=w,
                label="data (db-stats)", color=DEEP_PEACH)
    axes[0].set_xticks(x); axes[0].set_xticklabels(labels, rotation=15)
    axes[0].set_ylabel("GiB")
    axes[0].set_title("disk + data size")
    axes[0].legend()

    axes[1].bar(x, ssts, color=DEEP_MINT)
    axes[1].set_xticks(x); axes[1].set_xticklabels(labels, rotation=15)
    axes[1].set_ylabel("sstables")
    axes[1].set_title("sstable count")

    axes[2].bar(x, tomb, color=DEEP_RED)
    axes[2].set_xticks(x); axes[2].set_xticklabels(labels, rotation=15)
    axes[2].set_ylabel("tombstones")
    axes[2].set_title(f"tombstones (target deleted: {int(total_del):,})")

    fig.suptitle("delete-reclaim experiment snapshots", y=0.99, fontsize=12)
    finalize(fig, out_dir, "delete_experiment", caption, top=0.94, phase=False)


def plot_summary_table(summary_path, out_dir, caption):
    if not os.path.exists(summary_path):
        return
    rows = load_csv(summary_path)
    if not rows:
        return
    headers = ["phase", "dur_s", "GiB", "MiB/s", "wa_io", "wa_disk",
               "sa_disk", "pkg_W", "core_W", "cpu%",
               "p99_pt_µs", "p99_cm_µs", "RSS_MiB", "tomb_dns%"]
    data = []
    for r in rows:
        data.append([
            r.get("phase_name", "?"),
            f"{r['duration_s']:.1f}",
            f"{r['ingest_gib']:.2f}",
            f"{r['mean_write_mibs']:.1f}",
            f"{r['mean_wa_io']:.2f}",
            f"{r['mean_wa_disk']:.2f}",
            f"{r['mean_sa_disk']:.2f}",
            f"{r['mean_pkg_w']:.1f}",
            f"{r['mean_core_w']:.1f}",
            f"{r['mean_cpu_share']*100:.1f}",
            f"{r['max_point_p99_us']:.1f}",
            f"{r['max_commit_p99_us']:.1f}",
            f"{r['max_rss_mib']:.0f}",
            f"{r.get('max_tombstone_density', 0)*100:.1f}",
        ])
    fig, ax = plt.subplots(figsize=(13.5, 0.5 + 0.35 * (len(data) + 1)))
    ax.axis("off")
    tbl = ax.table(cellText=data, colLabels=headers, loc="center",
                   cellLoc="right", colLoc="center")
    tbl.auto_set_font_size(False)
    tbl.set_fontsize(9)
    tbl.scale(1.0, 1.4)
    for r in range(1, len(data) + 1):
        tbl[(r, 0)].set_text_props(ha="left")
    ax.set_title("per-phase summary  (means over per-tick samples)", pad=10)
    finalize(fig, out_dir, "summary_table", caption, top=1.0, phase=False)


# ============================================================ driver

def latest_run(base):
    if not os.path.isdir(base): return None
    runs = [d for d in os.listdir(base)
            if d.startswith("run_") and os.path.isdir(os.path.join(base, d))]
    if not runs: return None
    runs.sort()
    return os.path.join(base, runs[-1])


def detect_ncpu(rows0):
    """back out ncpu from the relationship cpu_share = cpu_time / (wall * ncpu).
    falls back to 8 if there isn't enough signal in the data."""
    best = 1
    for r in rows0:
        cpu_share = r.get("cpu_share", 0)
        utime = r.get("ru_utime_us_d", 0)
        stime = r.get("ru_stime_us_d", 0)
        if cpu_share > 0.05 and (utime + stime) > 1000:
            # cpu_share = (utime+stime) / (dt_us * ncpu)
            # dt_us is unknown but constant per tick; we can't recover ncpu cleanly.
            break
    try:
        return os.cpu_count() or best
    except Exception:
        return best


def main():
    script_dir = os.path.dirname(os.path.abspath(__file__))
    if len(sys.argv) > 1:
        out_dir = sys.argv[1]
    else:
        base = os.path.join(script_dir, "out")
        out_dir = latest_run(base) or base
        print(f"using: {out_dir}")
    samples_path = os.path.join(out_dir, "samples.csv")
    amp_path     = os.path.join(out_dir, "amplification.csv")
    summary_path = os.path.join(out_dir, "summary.csv")
    delete_path  = os.path.join(out_dir, "delete_experiment.csv")
    if not os.path.exists(samples_path):
        print(f"missing: {samples_path}", file=sys.stderr)
        sys.exit(1)

    rows = load_csv(samples_path)
    if not rows:
        print("empty samples.csv", file=sys.stderr); sys.exit(1)

    by_cf = partition_by_cf(rows)
    cfs = sorted(by_cf.keys())
    num_cfs = len(cfs)
    rows0 = by_cf[cfs[0]]
    cfg_rows = load_config(out_dir)
    meta = load_run_meta(out_dir)
    caption = caption_from_config(cfg_rows, meta)
    amp_rows = load_csv(amp_path)
    ncpu = int(meta.get("ncpu", detect_ncpu(rows0))) if meta else detect_ncpu(rows0)

    # use absolute phase bounds from summary.csv so sub-sample-interval
    # phases (e.g. a 0.1 s ingest) still show up in shading
    global PHASE_BOUNDS
    PHASE_BOUNDS = load_phase_bounds(out_dir)

    plot_ingest(rows0, out_dir, caption)
    plot_latency(rows0, out_dir, caption)
    plot_read_ops(rows0, out_dir, caption)
    plot_db_stats(rows0, meta, out_dir, caption)
    plot_cache(rows0, out_dir, caption)
    plot_disk(rows0, out_dir, caption)
    plot_integrity(rows0, out_dir, caption)
    plot_p99_vs_size(rows0, out_dir, caption)

    plot_per_cf_levels(by_cf, meta, out_dir, caption)
    plot_per_cf_level_bytes(by_cf, out_dir, caption)
    plot_per_cf_keys(by_cf, out_dir, caption)
    plot_per_cf_tombstones(by_cf, out_dir, caption)
    plot_per_cf_latency(by_cf, out_dir, caption)
    plot_per_cf_amp(by_cf, out_dir, caption)
    plot_per_cf_btree(by_cf, out_dir, caption)  # no-op unless a CF uses btree

    plot_amplification(amp_rows, out_dir, caption)
    plot_power(rows0, out_dir, caption)
    plot_commit_latency(rows0, num_cfs, out_dir, caption)
    plot_process_resources(rows0, ncpu, out_dir, caption)
    plot_objstore(rows0, meta, out_dir, caption)

    plot_delete_experiment(out_dir, delete_path, caption)
    plot_summary_table(summary_path, out_dir, caption)

    print(f"wrote {num_cfs}-CF run figures to {out_dir}")


if __name__ == "__main__":
    main()