From 9223d124a276e8b1d36e3d2a1f8104efe92c0c4e Mon Sep 17 00:00:00 2001
From: micah-olivas <micah5051olivas@gmail.com>
Date: Fri, 24 Jul 2026 23:31:07 -0700
Subject: [PATCH] Fix correctness bugs in pileup alignment metrics

The pileup viewer's rendering was sound, but the metrics layered on top of the
alignment could misreport what the reads actually showed. Since this view is
used to judge whether a well is recoverable, each of these could push the call
the wrong way.

Deletions were encoded identically to uncovered positions, as ("-", True), so
they were excluded from the identity denominator: a read with a 30 bp deletion
and no substitutions reported 100% identity, scoring better than a read with
five substitutions. Deletions and uncovered positions are now distinct symbols
(PILEUP_DEL vs PILEUP_NOCOV) and only uncovered positions are excluded.

The consensus padded deleted positions with N, keeping the sequence at
reference length, so a frameshift translated in frame and looked like a single
point substitution. The insert is now sliced in reference coordinates and
deleted positions are then dropped, so the frame actually shifts, and the net
indel is labelled as a frameshift or an in-frame indel.

n_reads held the post-filter aligned-row count while frac came from the
pre-filter group, so a header could read "0 reads (75%)". The empty-state
message interpolated the same zeroed count and so always said "(0 reads
unaligned)" - the one number it existed to report. Sections now carry n_reads
(group size) and n_aligned separately.

Problem-position flagging applied a >10% threshold with no depth floor, so a
single read covering a single position flagged it at 100%. Flagging now
requires PILEUP_MIN_FLAG_DEPTH (3) covering reads.

Also:
- All three row builders now route through one _rows_from_aligned_bam, so the
  symbol encoding cannot drift between them. _build_pileup_from_bam was dead
  code still carrying the old encoding.
- _run_pileup_alignment closes the parent's copy of minimap2's stdout and
  checks both exit codes; previously, if samtools died early minimap2 never
  got EPIPE and wait() could block forever, and a failed pipeline silently
  looked like an empty pileup.
- Deletions render in their own colour, with legend and tooltip entries, and
  the legend states that insertions are not shown.
- Canvas backing store is clamped to a maximum area; large pileups exceeded
  iOS Safari's cap and rendered blank.
- Inline JSON escapes "</" so a payload cannot close the script block.
- Removed an unreachable JS fallback branch ([] is truthy).

Adds tests/test_pileup_html.py (30 tests) covering all of the above; the
render layer needs no aligner, and the alignment-failure tests skip when
minimap2/samtools are absent. Full suite: 202 passed.
---
 src/usortm/demux/streakout.py | 655 ++++++++++++++++++++--------------
 tests/test_pileup_html.py     | 312 ++++++++++++++++
 2 files changed, 703 insertions(+), 264 deletions(-)
 create mode 100644 tests/test_pileup_html.py

diff --git a/src/usortm/demux/streakout.py b/src/usortm/demux/streakout.py
index 563f405..c3d1304 100644
--- a/src/usortm/demux/streakout.py
+++ b/src/usortm/demux/streakout.py
@@ -46,6 +46,25 @@
 logger = logging.getLogger(__name__)
 
 
+# ---------------------------------------------------------------------------
+# Pileup row encoding
+# ---------------------------------------------------------------------------
+# A pileup row holds one (symbol, is_match) tuple per reference position.
+# Deletions and uncovered positions must stay distinct: a deletion is positive
+# evidence that the read disagrees with the reference, while an uncovered
+# position is merely absent data.  Collapsing the two makes a read with a large
+# deletion score 100% identical to the reference.
+PILEUP_NOCOV = "-"   # read does not span this reference position
+PILEUP_DEL = "*"     # position deleted within the read's aligned span
+
+# Minimum read depth before a column may be flagged as a problem position.
+# Without this, a single read covering a single position flags it at 100%.
+PILEUP_MIN_FLAG_DEPTH = 3
+
+# Fraction of covering reads that must disagree with the reference to flag.
+PILEUP_MISMATCH_THRESHOLD = 0.10
+
+
 # ---------------------------------------------------------------------------
 # CIGAR classification (mirrors utils.extract_matches logic)
 # ---------------------------------------------------------------------------
@@ -818,45 +837,13 @@ def _build_pileup_from_bam(
     """
     if min_overlap_pos < 0:
         min_overlap_pos = ref_len // 2
-    rows = []
-    try:
-        with pysam.AlignmentFile(bam_path, "rb", check_sq=False) as bf:
-            for read in bf.fetch(until_eof=True):
-                if read.query_name not in read_names:
-                    continue
-                if read.is_unmapped or read.is_secondary or read.is_supplementary:
-                    continue
-                if min_overlap_pos and (
-                    read.reference_end is None
-                    or read.reference_start is None
-                    or read.reference_end <= min_overlap_pos
-                    or read.reference_start >= min_overlap_pos
-                ):
-                    continue
-                row = [("-", True)] * ref_len
-                # Use get_aligned_pairs() without with_seq=True to avoid
-                # requiring the MD tag, then look up the reference base from
-                # ref_seq directly.
-                pairs = read.get_aligned_pairs()
-                for qpos, rpos in pairs:
-                    if rpos is None or rpos >= ref_len:
-                        continue
-                    if qpos is None:
-                        row[rpos] = ("-", True)
-                    else:
-                        qbase = read.query_sequence[qpos]
-                        is_match = qbase.upper() == ref_seq[rpos].upper()
-                        row[rpos] = (qbase, is_match)
-                rows.append(row)
-    except Exception as exc:
-        logger.warning("BAM pileup extraction failed: %s", exc)
-
-    # Cluster reads by mismatch pattern (see _build_pileup_grid).
-    if rows:
-        rows.sort(key=lambda row: "".join(
-            "." if m or b == "-" else b.upper() for b, m in row
-        ))
-
+    # with_seq=False avoids requiring an MD tag; the reference base is looked
+    # up from ref_seq instead.
+    rows = _rows_from_aligned_bam(
+        bam_path, ref_seq, ref_len, min_overlap_pos,
+        read_names=read_names, with_seq=False, check_sq=False,
+    )
+    _sort_rows_by_mismatch(rows)
     return rows
 
 
@@ -927,7 +914,9 @@ def _make_section(ginfo: dict, grp: pd.DataFrame) -> Optional[dict]:
         )
         return {
             "ref_id": ref_id,
-            "n_reads": len(pileup_rows),
+            # n_reads is the group size; n_aligned is how many produced a row.
+            "n_reads": len(grp),
+            "n_aligned": len(pileup_rows),
             "frac": frac,
             "status": status,
             "is_recoverable": is_recoverable,
@@ -951,19 +940,26 @@ def _make_section(ginfo: dict, grp: pd.DataFrame) -> Optional[dict]:
             status = "Clean" if _cigar_is_clean(ginfo.get("cigar")) else "Mutation"
 
             variant_fasta = os.path.join(single_ref_dir, f"{ginfo['variant']}.fasta")
-            if os.path.exists(variant_fasta) and minimap2_path and samtools_path:
-                pileup_rows = _build_pileup_from_bam_realign(
-                    bam_path, read_names, variant_fasta,
-                    minimap2_path, samtools_path,
+            if not (os.path.exists(variant_fasta) and minimap2_path and samtools_path):
+                # No reference to display against: skip the group, matching
+                # _make_section rather than emitting an empty section.
+                logger.warning(
+                    "Pileup: skipping group %s (missing FASTA or aligner)",
+                    ginfo.get("variant"),
                 )
-                ref_seq = str(next(SeqIO.parse(variant_fasta, "fasta")).seq)
-            else:
-                pileup_rows = []
-                ref_seq = ""
+                continue
+
+            pileup_rows = _build_pileup_from_bam_realign(
+                bam_path, read_names, variant_fasta,
+                minimap2_path, samtools_path,
+            )
+            ref_seq = str(next(SeqIO.parse(variant_fasta, "fasta")).seq)
 
             group_sections.append({
                 "ref_id": ginfo["variant"],
-                "n_reads": len(pileup_rows),
+                # n_reads is the group size; n_aligned is how many produced a row.
+                "n_reads": len(read_names),
+                "n_aligned": len(pileup_rows),
                 "frac": frac,
                 "status": status,
                 "is_recoverable": is_recoverable,
@@ -1002,6 +998,157 @@ def _make_section(ginfo: dict, grp: pd.DataFrame) -> Optional[dict]:
         f.write(html)
 
 
+def _run_pileup_alignment(
+    fq_path: str,
+    ref_target: str,
+    out_bam: str,
+    minimap2_path: str,
+    samtools_path: str,
+) -> bool:
+    """Align *fq_path* to *ref_target*, writing a sorted, indexed BAM.
+
+    Returns True on success.  The parent's copy of minimap2's stdout is closed
+    once samtools owns it, so minimap2 gets EPIPE if samtools exits early
+    instead of leaving ``wait()`` blocked forever.  Both exit codes are
+    checked: a silently failed pipeline used to look like an empty pileup.
+    """
+    try:
+        mm2 = subprocess.Popen(
+            [minimap2_path, "-a", "--MD", "--secondary=no", ref_target, fq_path],
+            stdout=subprocess.PIPE, stderr=subprocess.DEVNULL,
+        )
+    except Exception as exc:
+        logger.warning("Pileup alignment could not start minimap2: %s", exc)
+        return False
+
+    try:
+        sort_proc = subprocess.run(
+            [samtools_path, "sort", "-o", out_bam],
+            stdin=mm2.stdout, stderr=subprocess.DEVNULL, check=False,
+        )
+    except Exception as exc:
+        logger.warning("Pileup alignment failed during samtools sort: %s", exc)
+        mm2.kill()
+        mm2.wait()
+        return False
+    finally:
+        if mm2.stdout is not None:
+            mm2.stdout.close()
+
+    mm2_rc = mm2.wait()
+    if mm2_rc != 0:
+        logger.warning("minimap2 exited %s during pileup alignment", mm2_rc)
+        return False
+    if sort_proc.returncode != 0:
+        logger.warning(
+            "samtools sort exited %s during pileup alignment", sort_proc.returncode
+        )
+        return False
+
+    index_proc = subprocess.run(
+        [samtools_path, "index", out_bam],
+        stderr=subprocess.DEVNULL, check=False,
+    )
+    if index_proc.returncode != 0:
+        # Only sequential reads follow, which do not need the index.
+        logger.debug("samtools index exited %s; continuing", index_proc.returncode)
+    return True
+
+
+def _rows_from_aligned_bam(
+    bam_path: str,
+    ref_seq: str,
+    ref_len: int,
+    min_overlap_pos: int,
+    read_names: set = None,
+    with_seq: bool = True,
+    check_sq: bool = True,
+) -> list:
+    """Parse an aligned BAM into pileup rows, one row per aligned read.
+
+    Each row is a list of ``(symbol, is_match)`` tuples indexed by reference
+    position.  *symbol* is a base character, :data:`PILEUP_DEL` for a deletion
+    inside the read's aligned span, or :data:`PILEUP_NOCOV` where the read does
+    not reach.
+
+    Insertions relative to the reference cannot be represented in this
+    column-per-reference-position layout and are dropped.  The rendered page
+    says so, so that a clean-looking pileup is not misread as "no insertions".
+
+    This is the single place the symbol encoding is applied; every pileup
+    builder routes through it so the encoding cannot drift between them.
+
+    Args:
+        read_names: When given, only reads with these names are included, and
+            secondary/supplementary alignments are skipped.
+        with_seq: Pass ``with_seq=True`` to ``get_aligned_pairs``, which needs
+            an MD tag.  Set False to look the reference base up from *ref_seq*
+            instead, for BAMs written without MD.
+        check_sq: Passed to :class:`pysam.AlignmentFile`; False tolerates a
+            missing ``@SQ`` header.
+    """
+    rows = []
+    try:
+        with pysam.AlignmentFile(bam_path, "rb", check_sq=check_sq) as bf:
+            # until_eof=True does not require an index, so a failed
+            # `samtools index` still yields a readable pileup.
+            for read in bf.fetch(until_eof=True):
+                if read_names is not None:
+                    if read.query_name not in read_names:
+                        continue
+                    if read.is_secondary or read.is_supplementary:
+                        continue
+                if read.is_unmapped:
+                    continue
+                # Skip reads that don't span the midpoint of the reference.
+                # 5' concatemers end before the midpoint; 3' concatemers
+                # start after it.  Only full-length reads cross it from
+                # both sides and cover the variable region.
+                if min_overlap_pos and (
+                    read.reference_end is None
+                    or read.reference_start is None
+                    or read.reference_end <= min_overlap_pos
+                    or read.reference_start >= min_overlap_pos
+                ):
+                    continue
+                row = [(PILEUP_NOCOV, False)] * ref_len
+                pairs = (
+                    read.get_aligned_pairs(with_seq=True) if with_seq
+                    else read.get_aligned_pairs()
+                )
+                for pair in pairs:
+                    qpos, rpos = pair[0], pair[1]
+                    if rpos is None or rpos >= ref_len:
+                        continue  # insertion, or beyond the reference end
+                    if qpos is None:
+                        row[rpos] = (PILEUP_DEL, False)
+                    else:
+                        qbase = read.query_sequence[qpos]
+                        is_match = qbase.upper() == ref_seq[rpos].upper()
+                        row[rpos] = (qbase, is_match)
+                rows.append(row)
+    except Exception as exc:
+        logger.warning("BAM parsing failed: %s", exc)
+
+    return rows
+
+
+def _sort_rows_by_mismatch(rows: list) -> None:
+    """Sort *rows* in place so reads sharing a mismatch pattern are adjacent.
+
+    Key: "." for a match or for no coverage, the base for a mismatch, ``*`` for
+    a deletion.  Identical patterns sort together and dots sort before letters,
+    so the cleanest reads come first and subpopulations stack visibly.
+    """
+    def _mismatch_key(row):
+        return "".join(
+            "." if (is_match or sym == PILEUP_NOCOV) else sym.upper()
+            for sym, is_match in row
+        )
+
+    rows.sort(key=_mismatch_key)
+
+
 def _build_pileup_from_bam_realign(
     bam_path: str,
     read_names: set,
@@ -1055,59 +1202,15 @@ def _build_pileup_from_bam_realign(
             return []
 
         # Align to target variant reference
-        try:
-            mm2 = subprocess.Popen(
-                [minimap2_path, "-a", "--MD", "--secondary=no", target_fasta, fq_path],
-                stdout=subprocess.PIPE, stderr=subprocess.DEVNULL,
-            )
-            subprocess.run(
-                [samtools_path, "sort", "-o", out_bam],
-                stdin=mm2.stdout, stderr=subprocess.DEVNULL, check=False,
-            )
-            mm2.wait()
-            subprocess.run(
-                [samtools_path, "index", out_bam],
-                stderr=subprocess.DEVNULL, check=False,
-            )
-        except Exception as exc:
-            logger.warning("Re-alignment for pileup failed: %s", exc)
+        if not _run_pileup_alignment(
+            fq_path, target_fasta, out_bam, minimap2_path, samtools_path
+        ):
             return []
 
-        # Parse aligned BAM into pileup rows
-        rows = []
-        try:
-            with pysam.AlignmentFile(out_bam, "rb") as bf:
-                for read in bf:
-                    if read.is_unmapped:
-                        continue
-                    if min_overlap_pos and (
-                        read.reference_end is None
-                        or read.reference_start is None
-                        or read.reference_end <= min_overlap_pos
-                        or read.reference_start >= min_overlap_pos
-                    ):
-                        continue
-                    row = [("-", True)] * ref_len
-                    pairs = read.get_aligned_pairs(with_seq=True)
-                    for qpos, rpos, rbase in pairs:
-                        if rpos is None or rpos >= ref_len:
-                            continue
-                        if qpos is None:
-                            row[rpos] = ("-", True)
-                        else:
-                            qbase = read.query_sequence[qpos]
-                            is_match = qbase.upper() == ref_seq[rpos].upper()
-                            row[rpos] = (qbase, is_match)
-                    rows.append(row)
-        except Exception as exc:
-            logger.warning("Re-aligned BAM pileup parsing failed: %s", exc)
-
-    # Cluster reads by mismatch pattern (see _build_pileup_grid).
-    if rows:
-        rows.sort(key=lambda row: "".join(
-            "." if m or b == "-" else b.upper() for b, m in row
-        ))
+        rows = _rows_from_aligned_bam(out_bam, ref_seq, ref_len, min_overlap_pos)
 
+    # Cluster reads by mismatch pattern so subpopulations are adjacent.
+    _sort_rows_by_mismatch(rows)
     return rows
 
 
@@ -1129,7 +1232,8 @@ def _build_pileup_grid(
     (typically >1 kb).  Pass 0 to disable the filter entirely.
 
     Returns a list of rows, where each row is a list of
-    (base_char, is_match) tuples indexed by reference position.
+    (symbol, is_match) tuples indexed by reference position.  See
+    :func:`_rows_from_aligned_bam` for the symbol encoding.
     """
     if min_overlap_pos < 0:
         min_overlap_pos = ref_len // 2
@@ -1144,56 +1248,12 @@ def _build_pileup_grid(
 
         # Align (use pre-built .mmi index if available)
         mm2_ref = ref_index if ref_index else ref_fasta
-        try:
-            mm2 = subprocess.Popen(
-                [minimap2_path, "-a", "--MD", "--secondary=no", mm2_ref, fq_path],
-                stdout=subprocess.PIPE, stderr=subprocess.DEVNULL,
-            )
-            subprocess.run(
-                [samtools_path, "sort", "-o", bam_path],
-                stdin=mm2.stdout, stderr=subprocess.DEVNULL, check=False,
-            )
-            mm2.wait()
-            subprocess.run(
-                [samtools_path, "index", bam_path],
-                stderr=subprocess.DEVNULL, check=False,
-            )
-        except Exception as exc:
-            logger.warning("Pileup alignment failed: %s", exc)
+        if not _run_pileup_alignment(
+            fq_path, mm2_ref, bam_path, minimap2_path, samtools_path
+        ):
             return []
 
-        # Parse BAM for pileup
-        rows = []
-        try:
-            with pysam.AlignmentFile(bam_path, "rb") as bf:
-                for read in bf:
-                    if read.is_unmapped:
-                        continue
-                    # Skip reads that don't span the midpoint of the reference.
-                    # 5' concatemers end before the midpoint; 3' concatemers
-                    # start after it.  Only full-length reads cross it from
-                    # both sides and cover the variable region.
-                    if min_overlap_pos and (
-                        read.reference_end is None
-                        or read.reference_start is None
-                        or read.reference_end <= min_overlap_pos
-                        or read.reference_start >= min_overlap_pos
-                    ):
-                        continue
-                    row = [("-", True)] * ref_len  # default: gap
-                    pairs = read.get_aligned_pairs(with_seq=True)
-                    for qpos, rpos, rbase in pairs:
-                        if rpos is None or rpos >= ref_len:
-                            continue
-                        if qpos is None:
-                            row[rpos] = ("-", True)  # deletion
-                        else:
-                            qbase = read.query_sequence[qpos]
-                            is_match = qbase.upper() == ref_seq[rpos].upper()
-                            row[rpos] = (qbase, is_match)
-                    rows.append(row)
-        except Exception as exc:
-            logger.warning("BAM parsing failed: %s", exc)
+        rows = _rows_from_aligned_bam(bam_path, ref_seq, ref_len, min_overlap_pos)
 
         if not rows:
             logger.warning(
@@ -1202,16 +1262,7 @@ def _build_pileup_grid(
             )
 
     # Cluster reads by mismatch pattern so subpopulations are adjacent.
-    # Key: "." for match/gap, actual base for mismatch → identical patterns
-    # sort together.  Fewest mismatches first (dots sort before letters).
-    if rows:
-        def _mismatch_key(row):
-            return "".join(
-                "." if is_match or base == "-" else base.upper()
-                for base, is_match in row
-            )
-        rows.sort(key=_mismatch_key)
-
+    _sort_rows_by_mismatch(rows)
     return rows
 
 
@@ -1221,11 +1272,23 @@ def _render_pileup_html(well_pos: str, candidate: dict,
     """Render the pileup HTML page for one well.
 
     Uses an HTML5 canvas matrix: each read is a row of colored cells.
-    Green = match, per-base color = mismatch, light gray = gap.
+    Match = neutral gray, per-base color = mismatch, dark red = deletion,
+    white = position not covered by the read.
+
+    Insertions relative to the reference are not shown (the matrix has one
+    column per reference position); the page states this in its legend.
     """
     import html as _html
     import json as _json
 
+    def _js(obj) -> str:
+        """JSON for embedding in an inline <script>.
+
+        ``json.dumps`` does not escape ``</``, so a payload containing
+        ``</script>`` would otherwise close the block early.
+        """
+        return _json.dumps(obj).replace("</", "<\\/")
+
     flanks_js = "null"
     if flank_lengths and (flank_lengths[0] or flank_lengths[1]):
         flanks_js = f"[{flank_lengths[0]},{flank_lengths[1]}]"
@@ -1244,113 +1307,151 @@ def _render_pileup_html(well_pos: str, candidate: dict,
         else:
             status_class = "status-other"
 
-        # Compute per-read identity from pileup data
+        ref_seq = g["ref_seq"]
+        ref_len = len(ref_seq)
+        pileup_rows = g["pileup_rows"]
+
+        # Read counts.  The size of the group and the number of reads that
+        # produced an aligned row are different numbers: the midpoint filter in
+        # the row builder drops concatemers.  Reporting only the latter beside a
+        # fraction derived from the former makes the two contradict each other.
+        n_aligned = g.get("n_aligned", len(pileup_rows))
+        n_group = g.get("n_reads", n_aligned)
+
+        # Identity pooled over reads.  Positions a read does not cover are
+        # excluded, but deletions inside its aligned span count against it —
+        # they are real disagreements with the reference, not missing data.
         identity_str = ""
-        if g["pileup_rows"]:
-            total_bases = 0
-            total_matches = 0
-            for row in g["pileup_rows"]:
-                aligned = [(b, m) for b, m in row if b != "-"]
-                total_bases += len(aligned)
-                total_matches += sum(1 for _, m in aligned if m)
-            if total_bases > 0:
-                identity = total_matches / total_bases
-                identity_str = f" &middot; Read identity: {identity:.1%}"
-
-        ref_len = len(g["ref_seq"])
-
-        header = (
-            f'<div class="group-header">'
-            f'<span class="ref-name">{_html.escape(g["ref_id"])}{star}</span>'
-            f'<span class="group-meta">'
-            f'{g["n_reads"]} reads ({g["frac"]:.0%}) &middot; '
-            f'{ref_len} bp &middot; '
-            f'Consensus: <span class="{status_class}">'
-            f'{_html.escape(g["status"])}</span>'
-            f'{identity_str}'
-            f'</span></div>'
-        )
-
-        # Encode pileup data compactly for JS:
-        # '.' = match, base letter = mismatch, '-' = gap
-        rows_encoded = []
-        for row in g["pileup_rows"]:
-            chars = []
-            for base_char, is_match in row:
-                if base_char == "-":
-                    chars.append("-")
-                elif is_match:
-                    chars.append(".")
-                else:
-                    chars.append(base_char.upper())
-            rows_encoded.append("".join(chars))
-
-        # Build consensus from pileup: majority base at each position.
-        # Also track positions where >10% of reads disagree with the
-        # reference — these are flagged in the ruler as problem positions.
+        total_bases = 0
+        total_matches = 0
+        for row in pileup_rows:
+            covered = [(s, m) for s, m in row if s != PILEUP_NOCOV]
+            total_bases += len(covered)
+            total_matches += sum(1 for _, m in covered if m)
+        if total_bases > 0:
+            identity = total_matches / total_bases
+            identity_str = f" &middot; Identity (all reads): {identity:.1%}"
+
+        # Consensus per reference position: the majority symbol, where a
+        # deletion competes on equal terms (a column deleted in most reads is
+        # deleted in the consensus).  A column is flagged as a problem position
+        # only when covered deeply enough for the fraction to mean anything.
         from collections import Counter
         consensus_encoded = []
         flagged_cols = []
-        ref_seq = g["ref_seq"]
-        _MISMATCH_THRESHOLD = 0.10
         for col_idx in range(ref_len):
             counts = Counter()
-            for row in g["pileup_rows"]:
-                base, _ = row[col_idx]
-                if base != "-":
-                    counts[base.upper()] += 1
+            for row in pileup_rows:
+                sym, _m = row[col_idx]
+                if sym != PILEUP_NOCOV:
+                    counts[sym.upper()] += 1
             total = sum(counts.values())
-            if counts:
-                cons_base = counts.most_common(1)[0][0]
-                if cons_base == ref_seq[col_idx].upper():
-                    consensus_encoded.append(".")
-                else:
-                    consensus_encoded.append(cons_base)
-                # Flag if >10% of reads differ from reference
-                ref_base = ref_seq[col_idx].upper()
-                ref_count = counts.get(ref_base, 0)
-                if total and (total - ref_count) / total > _MISMATCH_THRESHOLD:
-                    flagged_cols.append(col_idx)
-            else:
-                consensus_encoded.append("-")
+            if not counts:
+                consensus_encoded.append(PILEUP_NOCOV)
+                continue
+            cons_sym = counts.most_common(1)[0][0]
+            ref_base = ref_seq[col_idx].upper()
+            consensus_encoded.append("." if cons_sym == ref_base else cons_sym)
+            # Any non-reference symbol counts as disagreement, deletions included.
+            ref_count = counts.get(ref_base, 0)
+            if (
+                total >= PILEUP_MIN_FLAG_DEPTH
+                and (total - ref_count) / total > PILEUP_MISMATCH_THRESHOLD
+            ):
+                flagged_cols.append(col_idx)
         consensus_str = "".join(consensus_encoded)
 
-        # Reconstruct actual consensus DNA and translate the insert region
-        consensus_dna = "".join(
-            ref_seq[i] if c == "." else (c if c != "-" else "N")
-            for i, c in enumerate(consensus_encoded)
-        )
+        # Reconstruct the consensus insert and translate it.  Slice in
+        # reference coordinates first, then drop deleted positions, so that a
+        # deletion actually shortens the sequence and shifts the reading frame
+        # instead of being padded with N and silently preserving it.
         _ins_start = flank_lengths[0] if flank_lengths else 0
         _ins_end = ref_len - (flank_lengths[1] if flank_lengths else 0)
-        insert_dna = consensus_dna[_ins_start:_ins_end]
+        insert_syms = consensus_encoded[_ins_start:_ins_end]
+        insert_dna = "".join(
+            ref_seq[_ins_start + i] if c == "." else ("N" if c == PILEUP_NOCOV else c)
+            for i, c in enumerate(insert_syms)
+            if c != PILEUP_DEL
+        )
         ref_insert_dna = ref_seq[_ins_start:_ins_end]
 
+        # A net indel that is not a multiple of three shifts the frame, so the
+        # translation below diverges from that point on.  That is the correct
+        # display; label it rather than leaving the user to infer it.
+        _net_indel = len(insert_dna) - len(ref_insert_dna)
+        frameshift = bool(_net_indel % 3)
+
         cons_protein = ""
         ref_protein = ""
         try:
             from Bio.Seq import Seq as _BioSeq
-            _trim = len(insert_dna) - len(insert_dna) % 3
-            if _trim:
-                cons_protein = str(_BioSeq(insert_dna[:_trim]).translate())
-                ref_protein = str(_BioSeq(ref_insert_dna[:_trim]).translate())
+            _cons_trim = len(insert_dna) - len(insert_dna) % 3
+            _ref_trim = len(ref_insert_dna) - len(ref_insert_dna) % 3
+            if _cons_trim:
+                cons_protein = str(_BioSeq(insert_dna[:_cons_trim]).translate())
+            if _ref_trim:
+                ref_protein = str(_BioSeq(ref_insert_dna[:_ref_trim]).translate())
         except Exception:
             pass
 
-        ref_seq_js = _json.dumps(ref_seq)
-        rows_js = _json.dumps(rows_encoded)
-        cons_js = _json.dumps(consensus_str)
-        flagged_js = _json.dumps(flagged_cols)
+        indel_str = ""
+        if _net_indel:
+            _kind = "frameshift" if frameshift else "in-frame indel"
+            indel_str = (
+                f' &middot; <span class="indel-flag">'
+                f'{_net_indel:+d} bp {_kind}</span>'
+            )
+
+        reads_str = (
+            f'{n_group} reads ({g["frac"]:.0%})'
+            if n_aligned == n_group
+            else f'{n_aligned} of {n_group} reads aligned ({g["frac"]:.0%})'
+        )
+
+        header = (
+            f'<div class="group-header">'
+            f'<span class="ref-name">{_html.escape(g["ref_id"])}{star}</span>'
+            f'<span class="group-meta">'
+            f'{reads_str} &middot; '
+            f'{ref_len} bp &middot; '
+            f'Consensus: <span class="{status_class}">'
+            f'{_html.escape(g["status"])}</span>'
+            f'{identity_str}'
+            f'{indel_str}'
+            f'</span></div>'
+        )
+
+        # Encode pileup data compactly for JS: '.' = match, base letter =
+        # mismatch, '*' = deletion, '-' = position not covered by the read.
+        rows_encoded = []
+        for row in pileup_rows:
+            chars = []
+            for sym, is_match in row:
+                if sym in (PILEUP_NOCOV, PILEUP_DEL):
+                    chars.append(sym)
+                elif is_match:
+                    chars.append(".")
+                else:
+                    chars.append(sym.upper())
+            rows_encoded.append("".join(chars))
+
+        ref_seq_js = _js(ref_seq)
+        rows_js = _js(rows_encoded)
+        cons_js = _js(consensus_str)
+        flagged_js = _js(flagged_cols)
         n_rows = len(rows_encoded)
         n_cols = ref_len
 
         # Translation data for canvas rendering
-        ref_protein_js = _json.dumps(ref_protein) if ref_protein else "null"
-        cons_protein_js = _json.dumps(cons_protein) if cons_protein else "null"
+        ref_protein_js = _js(ref_protein) if ref_protein else "null"
+        cons_protein_js = _js(cons_protein) if cons_protein else "null"
 
         if n_rows == 0:
             pileup_block = (
                 f'<div class="pileup-empty">'
-                f'No aligned reads available ({g["n_reads"]} reads unaligned)'
+                f'No aligned reads to display &mdash; none of the '
+                f'{n_group} read(s) in this group aligned across the '
+                f'variable region'
                 f'</div>'
             )
         else:
@@ -1462,6 +1563,15 @@ def _render_pileup_html(well_pos: str, candidate: dict,
 .status-other {{
     color: #ef4444;
 }}
+.indel-flag {{
+    color: #991b1b;
+    font-weight: 600;
+}}
+.legend-note {{
+    color: var(--muted);
+    font-size: 0.75rem;
+    font-style: italic;
+}}
 .protein-seq {{
     margin-top: 0.5rem;
     font-family: 'Courier New', Courier, monospace;
@@ -1604,30 +1714,33 @@ def _render_pileup_html(well_pos: str, candidate: dict,
   var canvasW = Math.max(totalW, aaEndPx);
 
   var pileupH = refH + gap + consH + gap + nRows * cellH + aaGap + (hasAA ? aaH * 2 + 2 : 0);
-  canvas.width = canvasW * dpr;
-  canvas.height = pileupH * dpr;
+  // Browsers cap total canvas area (iOS Safari lowest, ~16.8M device px);
+  // exceeding it yields a silently blank canvas.  Scale the backing store down
+  // instead of losing the whole pileup — CSS size is unchanged either way.
+  var MAX_CANVAS_AREA = 16777216;
+  var effDpr = dpr;
+  if (canvasW * pileupH * effDpr * effDpr > MAX_CANVAS_AREA) {{
+    effDpr = Math.max(1, Math.sqrt(MAX_CANVAS_AREA / (canvasW * pileupH)));
+  }}
+  canvas.width = Math.floor(canvasW * effDpr);
+  canvas.height = Math.floor(pileupH * effDpr);
   canvas.style.width = canvasW + 'px';
   canvas.style.height = pileupH + 'px';
   var ctx = canvas.getContext('2d');
-  ctx.scale(dpr, dpr);
+  ctx.scale(effDpr, effDpr);
   var isDark = document.documentElement.getAttribute('data-theme') === 'dark';
   var matchColor = isDark ? '#4a5568' : '#c8ccd0';
   var vectorMatchColor = isDark ? '#3a4455' : '#dfe2e6';
   var gapColor = isDark ? '#ffffff' : '#ffffff';
+  var delColor = isDark ? '#7f1d1d' : '#991b1b';
   var refColor = isDark ? '#e0e0e0' : '#1e293b';
   var consMatchColor = matchColor;
   var baseColors = isDark
     ? {{'A':'#ff6b6b','T':'#339af0','C':'#ffa94d','G':'#ffd43b'}}
     : {{'A':'#e03131','T':'#1971c2','C':'#e8590c','G':'#e67700'}};
-  // Flagged columns: positions where >10% of reads disagree with reference.
-  // Falls back to consensus-derived mismatches when flaggedCols not provided.
+  // Flagged columns: positions where enough reads disagree with the reference,
+  // subject to a minimum depth.  Always supplied by the caller.
   var mismatchCols = flaggedCols || [];
-  if (!flaggedCols) {{
-    for (var _mi = 0; _mi < cons.length; _mi++) {{
-      var _ch = cons[_mi];
-      if (_ch !== '.' && _ch !== '-') mismatchCols.push(_mi);
-    }}
-  }}
   var triRowH = mismatchCols.length > 0 ? 13 : 0;
   function isVector(col) {{
     return flanks && (col < flanks[0] || col >= nCols - flanks[1]);
@@ -1638,12 +1751,12 @@ def _render_pileup_html(well_pos: str, candidate: dict,
   // --- Ruler ---
   var rulerH = (flanks ? 24 : 14) + triRowH;
   if (rulerCanvas) {{
-    rulerCanvas.width = canvasW * dpr;
-    rulerCanvas.height = rulerH * dpr;
+    rulerCanvas.width = Math.floor(canvasW * effDpr);
+    rulerCanvas.height = Math.floor(rulerH * effDpr);
     rulerCanvas.style.width = canvasW + 'px';
     rulerCanvas.style.height = rulerH + 'px';
     var rc = rulerCanvas.getContext('2d');
-    rc.scale(dpr, dpr);
+    rc.scale(effDpr, effDpr);
     var tickColor = isDark ? '#64748b' : '#94a3b8';
     var labelColor = isDark ? '#e0e0e0' : '#1e293b';
     var boundaryColor = isDark ? '#f59e0b' : '#d97706';
@@ -1775,6 +1888,8 @@ def _render_pileup_html(well_pos: str, candidate: dict,
       ctx.fillStyle = isVector(i) ? vectorMatchColor : consMatchColor;
     }} else if (ch === '-') {{
       ctx.fillStyle = gapColor;
+    }} else if (ch === '*') {{
+      ctx.fillStyle = delColor;
     }} else {{
       ctx.fillStyle = baseColors[ch] || '#94a3b8';
     }}
@@ -1790,6 +1905,8 @@ def _render_pileup_html(well_pos: str, candidate: dict,
         ctx.fillStyle = pickMatch(c);
       }} else if (ch === '-') {{
         ctx.fillStyle = gapColor;
+      }} else if (ch === '*') {{
+        ctx.fillStyle = delColor;
       }} else {{
         ctx.fillStyle = baseColors[ch] || '#94a3b8';
       }}
@@ -1818,7 +1935,8 @@ def _render_pileup_html(well_pos: str, candidate: dict,
     for (var ai = 0; ai < refAA.length; ai++) {{
       var ax = insStart * cellW + ai * aaCodonW;
       var rAA = refAA[ai];
-      var cAA = consAA[ai];
+      // consAA can be shorter than refAA once a deletion shortens the insert.
+      var cAA = ai < consAA.length ? consAA[ai] : '–';
       var match = rAA === cAA;
 
       // Ref AA row
@@ -1927,7 +2045,9 @@ def _render_pileup_html(well_pos: str, candidate: dict,
     }} else if (yp < consY + consH) {{
       var ch = cons[col];
       var base = ch === '.' ? refSeq[col] : ch;
-      var note = ch === '.' ? ' (match)' : ch === '-' ? '' : ' (mismatch)';
+      var note = ch === '.' ? ' (match)' : ch === '-' ? ' (no coverage)'
+        : ch === '*' ? '' : ' (mismatch)';
+      if (ch === '*') base = 'deletion';
       tooltip.textContent = rl + 'Consensus pos ' + (col + 1) + ': ' + base + note;
     }} else if (hasAA && yp >= aaY && yp < aaY + aaH * 2 + 2) {{
       var insStart = flanks[0];
@@ -1935,8 +2055,9 @@ def _render_pileup_html(well_pos: str, candidate: dict,
       if (aaIdx >= 0 && aaIdx < refAA.length) {{
         var isRefRow = yp < aaY + aaH;
         var which = isRefRow ? 'Ref' : 'Cons';
-        var aa = isRefRow ? refAA[aaIdx] : consAA[aaIdx];
-        var other = isRefRow ? consAA[aaIdx] : refAA[aaIdx];
+        var cAAt = aaIdx < consAA.length ? consAA[aaIdx] : '(past end)';
+        var aa = isRefRow ? refAA[aaIdx] : cAAt;
+        var other = isRefRow ? cAAt : refAA[aaIdx];
         var note = aa === other ? ' (match)' : ' \u2260 ' + (isRefRow ? 'Cons' : 'Ref') + ': ' + other;
         tooltip.textContent = which + ' AA ' + (aaIdx + 1) + ': ' + aa + note;
       }} else {{
@@ -1946,7 +2067,10 @@ def _render_pileup_html(well_pos: str, candidate: dict,
       var row_idx = Math.floor((yp - readsY) / cellH);
       if (row_idx >= 0 && row_idx < nRows) {{
         var ch = rows[row_idx][col];
-        var label = ch === '.' ? refSeq[col] + ' (match)' : ch === '-' ? 'gap' : ch + ' (mismatch)';
+        var label = ch === '.' ? refSeq[col] + ' (match)'
+          : ch === '-' ? 'not covered by this read'
+          : ch === '*' ? 'deletion'
+          : ch + ' (mismatch)';
         tooltip.textContent = rl + 'Read ' + (row_idx + 1) + ', pos ' + (col + 1) + ': ' + label;
       }} else {{
         tooltip.style.display = 'none'; return;
@@ -1977,9 +2101,10 @@ def _render_pileup_html(well_pos: str, candidate: dict,
     <span class="legend-item"><span class="legend-swatch" style="background:#1971c2;"></span> T</span>
     <span class="legend-item"><span class="legend-swatch" style="background:#e8590c;"></span> C</span>
     <span class="legend-item"><span class="legend-swatch" style="background:#e67700;"></span> G</span>
-    <span class="legend-item"><span class="legend-swatch" style="background:#ffffff;border:1px solid #d1d5db;"></span> Gap</span>
+    <span class="legend-item"><span class="legend-swatch" style="background:#991b1b;"></span> Deletion</span>
+    <span class="legend-item"><span class="legend-swatch" style="background:#ffffff;border:1px solid #d1d5db;"></span> Not covered</span>
     <span class="legend-item"><span class="legend-swatch" style="background:#1e293b;"></span> Reference</span>
-
+    <span class="legend-note">Insertions are not shown &mdash; the matrix has one column per reference position.</span>
 </div>
 {body}
 <script id="usortm-theme-sync">
@@ -2038,15 +2163,16 @@ def _generate_one_pick_pileup(
         minimap2_path, samtools_path, ref_index=ref_index,
     )
 
-    # Use the number of reads that actually aligned to the variable region
-    # (after flank filtering) as the displayed count, not the raw read count
-    # which includes concatemer reads that only cover the 5' flank.
+    # Reads that aligned across the variable region, versus every read routed
+    # to this well (which includes concatemers covering only the 5' flank).
+    # Both are reported now, so neither has to stand in for the other.
     n_variable_reads = len(pileup_rows)
+    n_well_reads = len(well_reads) if well_reads is not None else n_variable_reads
 
     candidate_info = {
         "plate": source_plate,
         "well": source_well,
-        "total_reads": n_variable_reads,
+        "total_reads": n_well_reads,
         "top_frac": consensus_fraction,
         "recoverable_variants": [],
         "groups": [{"variant": variant, "frac": consensus_fraction, "status": ""}],
@@ -2056,7 +2182,8 @@ def _generate_one_pick_pileup(
     _is_recoverable = cons_check in ("Perfect Match", "Silent Mutation")
     group_sections = [{
         "ref_id": variant,
-        "n_reads": n_variable_reads,
+        "n_reads": n_well_reads,
+        "n_aligned": n_variable_reads,
         "frac": consensus_fraction,
         "status": _display_status,
         "is_recoverable": _is_recoverable,
diff --git a/tests/test_pileup_html.py b/tests/test_pileup_html.py
new file mode 100644
index 0000000..6eebc36
--- /dev/null
+++ b/tests/test_pileup_html.py
@@ -0,0 +1,312 @@
+"""Tests for per-well pileup HTML generation.
+
+These cover the metrics layered on top of the alignment — read identity,
+consensus translation, read counts, and problem-position flagging — which are
+what a user actually reads off the page when judging whether a well is
+recoverable.  They exercise :func:`_render_pileup_html` directly with synthetic
+rows, so no minimap2/samtools/BAM is required.
+"""
+
+import os
+import re
+import tempfile
+
+import pytest
+
+from usortm.demux.streakout import (
+    PILEUP_DEL,
+    PILEUP_MIN_FLAG_DEPTH,
+    PILEUP_NOCOV,
+    _render_pileup_html,
+    _run_pileup_alignment,
+)
+
+# 66 bp: ATG, then 30 bp of A, 30 bp of T, then GGG.
+REF = "ATG" + "AAA" * 10 + "TTT" * 10 + "GGG"
+
+
+def make_row(spec, ref=REF):
+    """Build a pileup row from a spec string.
+
+    Characters: a base = observed base, ``*`` = deletion within the aligned
+    span, ``-`` = position not covered by the read.
+    """
+    assert len(spec) == len(ref), "spec must be reference-length"
+    row = []
+    for i, c in enumerate(spec):
+        if c == PILEUP_DEL:
+            row.append((PILEUP_DEL, False))
+        elif c == PILEUP_NOCOV:
+            row.append((PILEUP_NOCOV, False))
+        else:
+            row.append((c, c.upper() == ref[i].upper()))
+    return row
+
+
+def make_section(**overrides):
+    section = {
+        "ref_id": "variant_1",
+        "n_reads": 1,
+        "n_aligned": 1,
+        "frac": 1.0,
+        "status": "Mutation",
+        "is_recoverable": False,
+        "ref_seq": REF,
+        "pileup_rows": [],
+    }
+    section.update(overrides)
+    return section
+
+
+CANDIDATE = {
+    "plate": "1",
+    "well": "A1",
+    "recoverable_variants": [],
+    "total_reads": 40,
+    "top_frac": 0.75,
+}
+
+
+def render(sections, **kwargs):
+    return _render_pileup_html("1A1", CANDIDATE, sections, **kwargs)
+
+
+def identity_of(html):
+    m = re.search(r"Identity \(all reads\): ([\d.]+)%", html)
+    return float(m.group(1)) if m else None
+
+
+def flagged_of(html):
+    return re.search(r"var flaggedCols=(\[[^\]]*\])", html).group(1)
+
+
+def group_meta_of(html):
+    m = re.search(r'<span class="group-meta">(.*?)</span></div>', html, re.S)
+    return re.sub(r"<[^>]+>", "", m.group(1))
+
+
+def js_string(html, var):
+    m = re.search(r'var %s=("(?:[^"\\]|\\.)*")' % var, html)
+    return m.group(1).strip('"') if m else None
+
+
+class TestReadIdentity:
+    """Deletions are disagreements with the reference, not missing data."""
+
+    def test_deletion_is_not_scored_as_a_perfect_read(self):
+        row = make_row(REF[:20] + PILEUP_DEL * 30 + REF[50:])
+        html = render([make_section(pileup_rows=[row])])
+        assert identity_of(html) < 100.0
+
+    def test_deletion_scores_worse_than_a_few_substitutions(self):
+        deleted = make_row(REF[:20] + PILEUP_DEL * 30 + REF[50:])
+        subbed = list(REF)
+        for i in (21, 25, 31, 37, 43):
+            subbed[i] = "C" if REF[i] != "C" else "A"
+        substituted = make_row("".join(subbed))
+
+        id_del = identity_of(render([make_section(pileup_rows=[deleted])]))
+        id_sub = identity_of(render([make_section(pileup_rows=[substituted])]))
+        assert id_del < id_sub
+
+    def test_uncovered_positions_do_not_lower_identity(self):
+        """A short read that matches everywhere it reaches is still 100%."""
+        row = make_row(PILEUP_NOCOV * 20 + REF[20:46] + PILEUP_NOCOV * 20)
+        assert identity_of(render([make_section(pileup_rows=[row])])) == 100.0
+
+    def test_perfect_full_length_read_is_100_percent(self):
+        assert identity_of(render([make_section(pileup_rows=[make_row(REF)])])) == 100.0
+
+
+class TestConsensusTranslation:
+    """A frameshift must shift the frame rather than being padded with N."""
+
+    def test_single_base_deletion_shortens_the_consensus_protein(self):
+        row = make_row(REF[:10] + PILEUP_DEL + REF[11:])
+        html = render([make_section(pileup_rows=[row])], flank_lengths=(3, 3))
+        assert len(js_string(html, "consAA")) < len(js_string(html, "refAA"))
+
+    def test_single_base_deletion_diverges_downstream(self):
+        row = make_row(REF[:10] + PILEUP_DEL + REF[11:])
+        html = render([make_section(pileup_rows=[row])], flank_lengths=(3, 3))
+        ref_aa, cons_aa = js_string(html, "refAA"), js_string(html, "consAA")
+        assert ref_aa[3:] != cons_aa[3:]
+
+    def test_frameshift_is_labelled(self):
+        row = make_row(REF[:10] + PILEUP_DEL + REF[11:])
+        html = render([make_section(pileup_rows=[row])], flank_lengths=(3, 3))
+        assert "-1 bp frameshift" in html
+
+    def test_in_frame_deletion_is_not_called_a_frameshift(self):
+        row = make_row(REF[:12] + PILEUP_DEL * 3 + REF[15:])
+        html = render([make_section(pileup_rows=[row])], flank_lengths=(3, 3))
+        assert "-3 bp in-frame indel" in html
+        assert "frameshift" not in html
+
+    def test_clean_read_reports_no_indel(self):
+        html = render([make_section(pileup_rows=[make_row(REF)])],
+                      flank_lengths=(3, 3))
+        assert "bp frameshift" not in html
+        assert "in-frame indel" not in html
+
+    def test_substitution_only_keeps_reference_length_protein(self):
+        subbed = list(REF)
+        subbed[12] = "C" if REF[12] != "C" else "A"
+        html = render([make_section(pileup_rows=[make_row("".join(subbed))])],
+                      flank_lengths=(3, 3))
+        assert len(js_string(html, "consAA")) == len(js_string(html, "refAA"))
+
+
+class TestReadCounts:
+    """Group size and aligned-row count are different numbers; report both."""
+
+    def test_empty_group_reports_its_real_size(self):
+        html = render([make_section(n_reads=17, n_aligned=0, pileup_rows=[])])
+        empty = re.search(r'<div class="pileup-empty">(.*?)</div>', html, re.S).group(1)
+        assert "17" in empty
+
+    def test_header_shows_aligned_of_group_when_they_differ(self):
+        rows = [make_row(REF)] * 12
+        html = render([make_section(n_reads=30, n_aligned=12, frac=0.45,
+                                   pileup_rows=rows)])
+        assert "12 of 30 reads aligned" in group_meta_of(html)
+
+    def test_header_is_unqualified_when_all_reads_aligned(self):
+        rows = [make_row(REF)] * 5
+        html = render([make_section(n_reads=5, n_aligned=5, pileup_rows=rows)])
+        meta = group_meta_of(html)
+        assert "5 reads" in meta and "of" not in meta.split("&middot;")[0]
+
+
+class TestFlaggedColumns:
+    """Problem positions need a depth floor to mean anything."""
+
+    def test_single_read_column_is_not_flagged(self):
+        row = make_row(PILEUP_NOCOV * 30 + "C" + PILEUP_NOCOV * 35)
+        assert flagged_of(render([make_section(pileup_rows=[row])])) == "[]"
+
+    def test_column_at_min_depth_is_flagged(self):
+        rows = [make_row(PILEUP_NOCOV * 30 + "C" + PILEUP_NOCOV * 35)
+                for _ in range(PILEUP_MIN_FLAG_DEPTH)]
+        html = render([make_section(n_reads=len(rows), n_aligned=len(rows),
+                                   pileup_rows=rows)])
+        assert flagged_of(html) == "[30]"
+
+    def test_deletions_count_toward_flagging(self):
+        """A column deleted in most reads is a problem position."""
+        rows = [make_row(REF[:30] + PILEUP_DEL + REF[31:]) for _ in range(4)]
+        html = render([make_section(n_reads=4, n_aligned=4, pileup_rows=rows)])
+        assert "30" in flagged_of(html)
+
+    def test_clean_deep_pileup_flags_nothing(self):
+        rows = [make_row(REF) for _ in range(10)]
+        html = render([make_section(n_reads=10, n_aligned=10, pileup_rows=rows)])
+        assert flagged_of(html) == "[]"
+
+
+class TestPayloadEncoding:
+    def test_deletion_and_no_coverage_are_distinct_symbols(self):
+        row = make_row(REF[:20] + PILEUP_DEL * 5 + PILEUP_NOCOV * 5 + REF[30:])
+        html = render([make_section(pileup_rows=[row])])
+        encoded = re.search(r'var rows=\["([^"]*)"\]', html).group(1)
+        assert PILEUP_DEL in encoded and PILEUP_NOCOV in encoded
+
+    def test_majority_deletion_appears_in_consensus(self):
+        rows = [make_row(REF[:30] + PILEUP_DEL + REF[31:]) for _ in range(5)]
+        html = render([make_section(n_reads=5, n_aligned=5, pileup_rows=rows)])
+        assert js_string(html, "cons")[30] == PILEUP_DEL
+
+    def test_script_close_cannot_escape_the_inline_block(self):
+        html = render([make_section(pileup_rows=[make_row(REF)])])
+        for var in ("ref", "cons"):
+            assert "</" not in js_string(html, var)
+
+    def test_page_documents_the_insertion_limitation(self):
+        html = render([make_section(pileup_rows=[make_row(REF)])])
+        assert "Insertions are not shown" in html
+
+    def test_canvas_area_is_clamped(self):
+        html = render([make_section(pileup_rows=[make_row(REF)])])
+        assert "MAX_CANVAS_AREA" in html
+
+
+class TestRendersAtAll:
+    @pytest.mark.parametrize("flanks", [None, (3, 3)])
+    def test_produces_a_complete_document(self, flanks):
+        rows = [make_row(REF), make_row(REF[:40] + PILEUP_DEL * 2 + REF[42:])]
+        html = render([make_section(n_reads=2, n_aligned=2, pileup_rows=rows)],
+                      flank_lengths=flanks)
+        assert html.startswith("<!DOCTYPE html>")
+        assert html.rstrip().endswith("</html>")
+        assert "drawPileup(" in html
+
+    def test_multiple_groups_each_get_a_section(self):
+        sections = [
+            make_section(ref_id="var_A", frac=0.7, pileup_rows=[make_row(REF)]),
+            make_section(ref_id="var_B", frac=0.3, pileup_rows=[make_row(REF)]),
+        ]
+        html = render(sections)
+        assert "var_A" in html and "var_B" in html
+        # One call site per group (the shared function definition also matches
+        # a bare "drawPileup(", so count invocations by their first argument).
+        assert html.count('drawPileup("pileup-') == 2
+
+    def test_zero_length_reference_does_not_raise(self):
+        html = render([make_section(ref_seq="", n_reads=3, n_aligned=0,
+                                   pileup_rows=[])])
+        assert html.startswith("<!DOCTYPE html>")
+
+
+class TestAlignmentFailureHandling:
+    """A broken toolchain must fail fast, not block on the minimap2 pipe.
+
+    The parent has to drop its copy of minimap2's stdout once samtools owns it,
+    or minimap2 never receives EPIPE when samtools dies and ``wait()`` blocks
+    forever.  These assert the failure is reported instead.
+    """
+
+    @pytest.fixture
+    def alignment_inputs(self):
+        with tempfile.TemporaryDirectory() as td:
+            fq = os.path.join(td, "reads.fastq")
+            fa = os.path.join(td, "ref.fasta")
+            with open(fq, "w") as fh:
+                fh.write("@r1\nACGTACGTACGT\n+\nIIIIIIIIIIII\n")
+            with open(fa, "w") as fh:
+                fh.write(">ref\nACGTACGTACGTACGTACGT\n")
+            yield fq, fa, os.path.join(td, "out.bam")
+
+    def test_missing_minimap2_returns_false(self, alignment_inputs):
+        fq, fa, out = alignment_inputs
+        assert _run_pileup_alignment(
+            fq, fa, out, "/nonexistent/minimap2", "/nonexistent/samtools"
+        ) is False
+
+    def test_missing_samtools_returns_false(self, alignment_inputs):
+        """The case that could previously hang: samtools exits, minimap2 waits."""
+        fq, fa, out = alignment_inputs
+        mm2 = pytest.importorskip("usortm.demux.deps").find_minimap2()
+        if not mm2:
+            pytest.skip("minimap2 not installed")
+        assert _run_pileup_alignment(
+            fq, fa, out, mm2, "/nonexistent/samtools"
+        ) is False
+
+    def test_missing_reference_returns_false(self, alignment_inputs):
+        fq, _fa, out = alignment_inputs
+        deps = pytest.importorskip("usortm.demux.deps")
+        mm2, st = deps.find_minimap2(), deps.find_samtools()
+        if not (mm2 and st):
+            pytest.skip("minimap2/samtools not installed")
+        assert _run_pileup_alignment(
+            fq, out + ".missing.fasta", out, mm2, st
+        ) is False
+
+    def test_successful_alignment_returns_true(self, alignment_inputs):
+        fq, fa, out = alignment_inputs
+        deps = pytest.importorskip("usortm.demux.deps")
+        mm2, st = deps.find_minimap2(), deps.find_samtools()
+        if not (mm2 and st):
+            pytest.skip("minimap2/samtools not installed")
+        assert _run_pileup_alignment(fq, fa, out, mm2, st) is True
+        assert os.path.exists(out)