From 096999dc40f6996f13cf41647b391831cfd60a4a Mon Sep 17 00:00:00 2001
From: "Mitchell R. Vollger" <mvollger@gmail.com>
Date: Fri, 14 Aug 2026 11:53:42 -0700
Subject: [PATCH 1/4] style: format shell blocks with snakefmt

snakefmt's shell formatting (shfmt) is active in the current pixi env
and reformats every shell block. Two pre-existing '((cat $(...)'
subshells are disambiguated to '( (cat ...' because shfmt correctly
rejects the former as arithmetic. No functional changes.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
---
 workflow/rules/apply-model.smk     |  32 ++++----
 workflow/rules/coverages.smk       |  46 ++++++------
 workflow/rules/decorated-reads.smk |  20 ++---
 workflow/rules/fire-peaks.smk      | 115 +++++++++++++++--------------
 workflow/rules/levio.smk           |   8 +-
 workflow/rules/stats.smk           |  28 +++----
 workflow/rules/track-hub.smk       |  37 +++++-----
 7 files changed, 144 insertions(+), 142 deletions(-)

diff --git a/workflow/rules/apply-model.smk b/workflow/rules/apply-model.smk
index 8a2a6bf7d..af721bfd2 100644
--- a/workflow/rules/apply-model.smk
+++ b/workflow/rules/apply-model.smk
@@ -8,6 +8,10 @@ rule fire:
     output:
         cram="results/{sm}/{sm}-fire-{v}-filtered.cram",
         crai="results/{sm}/{sm}-fire-{v}-filtered.cram.crai",
+    benchmark:
+        "results/{sm}/additional-outputs-{v}/benchmarks/{sm}-fire-bam.txt"
+    conda:
+        DEFAULT_ENV
     threads: 32
     resources:
         mem_mb=32 * 1024,
@@ -17,10 +21,6 @@ rule fire:
         min_ave_msp_size=config.get("min_ave_msp_size", 10),
         use_ont=USE_ONT,
         flag=FILTER_FLAG,
-    benchmark:
-        "results/{sm}/additional-outputs-{v}/benchmarks/{sm}-fire-bam.txt"
-    conda:
-        DEFAULT_ENV
     shell:
         """
         samtools view -@ {threads} -u -F {params.flag} {input.bam} \
@@ -37,7 +37,7 @@ rule fire:
                 --input-fmt-option required_fields=0x1bff \
                 --write-index -o {output.cram}
 
-        # check if the cram file has zero reads 
+        # check if the cram file has zero reads
         reads_in_header=$(samtools view {output.cram} | head | wc -l || true)
         if [ $reads_in_header -eq 0 ]; then
             printf "\nNo reads passed filters exiting...\n\nPlease review https://fiberseq.github.io/quick-start.html to make sure the input BAM has been correctly processed.\n\n"
@@ -51,9 +51,9 @@ rule fire_sites_chrom:
         cram=rules.fire.output.cram,
     output:
         bed=temp("temp/{sm}/chrom/{v}-{chrom}.sorted.bed.gz"),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     resources:
         mem_mb=16 * 1024,
     params:
@@ -62,13 +62,13 @@ rule fire_sites_chrom:
         """
         samtools view -@ {threads} -u {input.cram} {wildcards.chrom} \
             | {FT_EXE} fire -t {threads} --extract - \
-                | LC_ALL=C sort --parallel={threads} \
-                    -k1,1 -k2,2n -k3,3n -k4,4 \
-                | bioawk -tc hdr '$10<={params.min_fdr}' \
-                | (grep '\\S' || true) \
-                | (grep -v '^#' || true) \
-                | bgzip -@ {threads} \
-            > {output.bed}
+            | LC_ALL=C sort --parallel={threads} \
+                -k1,1 -k2,2n -k3,3n -k4,4 \
+            | bioawk -tc hdr '$10<={params.min_fdr}' \
+            | (grep '\\S' || true) \
+            | (grep -v '^#' || true) \
+            | bgzip -@ {threads} \
+                >{output.bed}
         """
 
 
@@ -79,12 +79,12 @@ rule fire_sites:
         ),
     output:
         bed="results/{sm}/additional-outputs-{v}/fire-peaks/{sm}-{v}-fire-elements.bed.gz",
-    threads: 1
     conda:
         DEFAULT_ENV
+    threads: 1
     shell:
         """
-        cat {input.beds} > {output.bed}
+        cat {input.beds} >{output.bed}
         """
 
 
@@ -93,9 +93,9 @@ rule fire_sites_index:
         bed=rules.fire_sites.output.bed,
     output:
         tbi=rules.fire_sites.output.bed + ".tbi",
-    threads: 1
     conda:
         DEFAULT_ENV
+    threads: 1
     shell:
         """
         tabix -p bed {input.bed}
diff --git a/workflow/rules/coverages.smk b/workflow/rules/coverages.smk
index e7eb5ea89..ec3b6f955 100644
--- a/workflow/rules/coverages.smk
+++ b/workflow/rules/coverages.smk
@@ -10,18 +10,18 @@ rule genome_bedgraph:
     output:
         bg=temp("temp/{sm}/coverage/{sm}-{v}.bed.gz"),
         tbi=temp("temp/{sm}/coverage/{sm}-{v}.bed.gz.tbi"),
-    threads: 16
     shadow:
         "minimal"
     conda:
         DEFAULT_ENV
+    threads: 16
     shell:
-        """ 
+        """
         mosdepth -F 4 -f {input.ref} -t {threads} tmp {input.cram}
         bgzip -cd tmp.per-base.bed.gz \
-            | LC_ALL=C sort --parallel={threads} -k1,1 -k2,2n -k3,3n -k4,4  \
+            | LC_ALL=C sort --parallel={threads} -k1,1 -k2,2n -k3,3n -k4,4 \
             | bgzip -@ {threads} \
-        > {output.bg}
+                >{output.bg}
         tabix -f -p bed {output.bg}
         """
 
@@ -33,13 +33,13 @@ rule coverage:
         cov="results/{sm}/additional-outputs-{v}/coverage/{sm}-{v}-median-coverage.txt",
         minimum="results/{sm}/additional-outputs-{v}/coverage/{sm}-{v}-minimum-coverage.txt",
         maximum="results/{sm}/additional-outputs-{v}/coverage/{sm}-{v}-maximum-coverage.txt",
+    benchmark:
+        "results/{sm}/additional-outputs-{v}/benchmarks/coverage/{sm}.txt"
     conda:
         "../envs/python.yaml"
     threads: 1
     resources:
         mem_mb=64 * 1024,
-    benchmark:
-        "results/{sm}/additional-outputs-{v}/benchmarks/coverage/{sm}.txt"
     params:
         coverage_within_n_sd=COVERAGE_WITHIN_N_SD,
         min_coverage=MIN_COVERAGE,
@@ -57,18 +57,18 @@ rule fiber_locations_chromosome:
         crai=rules.fire.output.crai,
     output:
         bed=temp("temp/{sm}/coverage/{v}-{chrom}.fiber-locations.bed.gz"),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
         # get fiber locations
         (samtools view -@ {threads} -u {input.cram} {wildcards.chrom} \
             | {FT_EXE} extract -t {threads} -s --all - \
-            | hck -F '#ct' -F st -F en -F fiber -F strand -F HP ) \
+            | hck -F '#ct' -F st -F en -F fiber -F strand -F HP) \
             | (grep -v "^#" || true) \
             | bgzip -@ {threads} \
-        > {output.bed}
+                >{output.bed}
         """
 
 
@@ -91,24 +91,24 @@ rule fiber_locations:
         filtered_tbi=temp(
             "temp/{sm}/coverage/filtered-for-coverage/{v}-fiber-locations.bed.gz.tbi"
         ),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     params:
         max_frac_overlap=0.2,
     shell:
         """
-        cat {input.fibers} > {output.bed}
+        cat {input.fibers} >{output.bed}
         tabix -f -p bed {output.bed}
-        
+
         # get filtered fiber locations
         MIN=$(cat {input.minimum})
         MAX=$(cat {input.maximum})
         bedtools intersect -header -sorted -v -f {params.max_frac_overlap} \
             -a {output.bed} \
             -b <(bgzip -cd {input.bg} | awk -v MAX="$MAX" -v MIN="$MIN" '$4 <= MIN || $4 >= MAX') \
-        | bgzip -@ {threads} \
-        > {output.filtered}
+            | bgzip -@ {threads} \
+                >{output.filtered}
         tabix -f -p bed {output.filtered}
         """
 
@@ -122,9 +122,9 @@ rule exclude_from_shuffle:
         fai=ancient(FAI),
     output:
         bed="results/{sm}/additional-outputs-{v}/coverage/exclude-from-shuffles.bed.gz",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     params:
         exclude=lambda wc: " ".join(EXCLUDES) if EXCLUDES else "",
     shell:
@@ -139,7 +139,7 @@ rule exclude_from_shuffle:
             | bedtools sort \
             | bedtools merge \
             | bgzip -@ {threads} \
-        > {output.bed}
+                >{output.bed}
         """
 
 
@@ -154,12 +154,12 @@ rule unreliable_coverage_regions:
         bed_tbi="results/{sm}/additional-outputs-{v}/coverage/unreliable-coverage-regions.bed.gz.tbi",
         tmp=temp("temp/{sm}/additional-outputs-{v}/unreliable-coverage-regions.bed"),
         bb="results/{sm}/trackHub-{v}/bb/unreliable-coverage-regions.bb",
+    conda:
+        DEFAULT_ENV
     threads: 4
     params:
         min_len=MIN_UNRELIABLE_COVERAGE_LEN,
         bed3_as=workflow.source_path("../templates/bed3.as"),
-    conda:
-        DEFAULT_ENV
     shell:
         """
         MIN=$(cat {input.minimum})
@@ -170,15 +170,15 @@ rule unreliable_coverage_regions:
             | bedtools merge -i - \
             | awk '$3-$2 >= {params.min_len}' \
             | bgzip -@ {threads} \
-        > {output.bed}
+                >{output.bed}
 
         # bigbed
         # for some reason bigtools gives a too many files open error when reading from stdin
-        bedtools merge -i {output.bed} > {output.tmp}
+        bedtools merge -i {output.bed} >{output.tmp}
         bigtools bedtobigbed \
-                -s start -a {params.bed3_as} \
-                {output.tmp} {input.fai} {output.bb}
+            -s start -a {params.bed3_as} \
+            {output.tmp} {input.fai} {output.bb}
 
-        # index 
+        # index
         tabix -f -p bed {output.bed}
         """
diff --git a/workflow/rules/decorated-reads.smk b/workflow/rules/decorated-reads.smk
index b39eb5f6f..6a92ee80d 100644
--- a/workflow/rules/decorated-reads.smk
+++ b/workflow/rules/decorated-reads.smk
@@ -11,20 +11,20 @@ rule decorate_fibers_chromosome:
     output:
         bed=temp("temp/{sm}/decorate/{v}-{chrom}.bed.gz"),
         decorated=temp("temp/{sm}/decorate/{v}-{chrom}.dec.bed.gz"),
-    threads: 4
-    resources:
-        mem_mb=get_mem_mb,
     # the following steps can take a while so this helps the pipeline start this earlier.
     priority: 10
     conda:
         DEFAULT_ENV
+    threads: 4
+    resources:
+        mem_mb=get_mem_mb,
     shell:
         """
         samtools view -@ {threads} -u {input.cram} {wildcards.chrom} \
             | {FT_EXE} track-decorators -t {threads} --bed12 {output.bed} \
             | sort -k1,1 -k2,2n -k3,3n -k4,4 \
             | bgzip -@ {threads} \
-        > {output.decorated}
+                >{output.decorated}
         """
 
 
@@ -41,12 +41,12 @@ rule decorate_fibers_1:
         bb="results/{sm}/trackHub-{v}/bb/fire-fibers.bb",
     benchmark:
         "results/{sm}/additional-outputs-{v}/benchmarks/decorate_fibers_1/{sm}.txt"
-    threads: 8
-    resources:
-        runtime=240,
     priority: 10
     conda:
         DEFAULT_ENV
+    threads: 8
+    resources:
+        runtime=240,
     params:
         bed_as=workflow.source_path("../templates/bed12_filter.as"),
         nzooms=NZOOMS,
@@ -79,12 +79,12 @@ rule decorate_fibers_2:
         #bed=temp("temp/{sm}/trackHub-{v}/bb/fire-fiber-decorators.bed.gz"),
     benchmark:
         "results/{sm}/additional-outputs-{v}/benchmarks/decorate_fibers_2/{sm}.txt"
-    threads: 8
-    resources:
-        runtime=60 * 16,
     priority: 10
     conda:
         DEFAULT_ENV
+    threads: 8
+    resources:
+        runtime=60 * 16,
     params:
         dec_as=workflow.source_path("../templates/decoration.as"),
         nzooms=NZOOMS,
diff --git a/workflow/rules/fire-peaks.smk b/workflow/rules/fire-peaks.smk
index e1be2144a..fbccff0e5 100644
--- a/workflow/rules/fire-peaks.smk
+++ b/workflow/rules/fire-peaks.smk
@@ -6,9 +6,9 @@ rule filtered_and_shuffled_fiber_locations_chromosome:
         fai=ancient(FAI),
     output:
         shuffled=temp("temp/{sm}/shuffle/{v}-{chrom}.fiber-locations-shuffled.bed.gz"),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
         tabix {input.filtered} {wildcards.chrom} \
@@ -17,9 +17,9 @@ rule filtered_and_shuffled_fiber_locations_chromosome:
                 -excl {input.exclude} \
                 -i - \
                 -g {input.fai} \
-            |  sort -k1,1 -k2,2n -k3,3n -k4,4 \
+            | sort -k1,1 -k2,2n -k3,3n -k4,4 \
             | bgzip -@ {threads} \
-        > {output.shuffled}
+                >{output.shuffled}
         """
 
 
@@ -29,16 +29,16 @@ rule shuffled_pileup_chromosome:
         shuffled=rules.filtered_and_shuffled_fiber_locations_chromosome.output.shuffled,
     output:
         bed=temp("temp/{sm}/shuffle/{v}-{chrom}.pileup.bed.gz"),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
         {FT_EXE} pileup {input.cram} --rgn {wildcards.chrom} -t {threads} \
             --fiber-coverage --shuffle {input.shuffled} \
             --no-msp --no-nuc \
             | bgzip -@ {threads} \
-        > {output.bed}    
+                >{output.bed}
         """
 
 
@@ -52,12 +52,12 @@ rule shuffled_pileup:
     output:
         bed=temp("temp/{sm}/shuffle/{v}-shuffled-pileup.bed.gz"),
         tbi=temp("temp/{sm}/shuffle/{v}-shuffled-pileup.bed.gz.tbi"),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
-        cat {input.beds} > {output.bed}
+        cat {input.beds} >{output.bed}
         tabix -p bed {output.bed}
         """
 
@@ -74,11 +74,11 @@ rule fdr_table:
         tbl="results/{sm}/additional-outputs-{v}/fire-peaks/{sm}-{v}-fire-score-to-fdr.tbl",
     conda:
         "../envs/python.yaml"
-    params:
-        script=workflow.source_path("../scripts/fdr-table.py"),
     threads: 16
     resources:
         mem_mb=get_mem_mb_xl,
+    params:
+        script=workflow.source_path("../scripts/fdr-table.py"),
     shell:
         """
         export POLARS_MAX_THREADS={threads}
@@ -98,16 +98,16 @@ rule pileup_chromosome:
         bam=rules.fire.output.cram,
     output:
         bed=temp("temp/{sm}/{v}-{chrom}.pileup.bed.gz"),
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
         {FT_EXE} pileup -t {threads} \
             --haps --fiber-coverage \
             {input.bam} --rgn {wildcards.chrom} \
             | bgzip -@ {threads} \
-            > {output.bed}
+                >{output.bed}
         """
 
 
@@ -117,13 +117,13 @@ rule fdr_track_chromosome:
         fdr_tbl=rules.fdr_table.output.tbl,
     output:
         bed=temp("temp/{sm}/fire-peaks/{v}-{chrom}-FDR.track.bed"),
-    threads: 4
     conda:
         "../envs/python.yaml"
-    params:
-        script=workflow.source_path("../scripts/fdr-table.py"),
+    threads: 4
     resources:
         mem_mb=get_mem_mb_xl,
+    params:
+        script=workflow.source_path("../scripts/fdr-table.py"),
     shell:
         """
         export POLARS_MAX_THREADS={threads}
@@ -144,24 +144,24 @@ rule pileup:
         fofn=temp("temp/{sm}/fire/fire-{v}-pileup.fofn"),
         bed="results/{sm}/{sm}-fire-{v}-pileup.bed.gz",
         tbi="results/{sm}/{sm}-fire-{v}-pileup.bed.gz.tbi",
-    threads: 8
     conda:
         DEFAULT_ENV
+    threads: 8
     shell:
         """
         printf '\nMaking FOFN\n'
-        echo {input.beds} > {output.fofn}
-        
+        echo {input.beds} >{output.fofn}
+
         printf '\nMake header\n'
-        ((cat $(cat {output.fofn}) | grep "^#" | head -n 1) || true) \
+        ( (cat $(cat {output.fofn}) | grep "^#" | head -n 1) || true) \
             | bgzip -@ {threads} \
-            > {output.bed}
+                >{output.bed}
 
         printf '\nConcatenating\n'
         cat $(cat {output.fofn}) \
             | grep -v "^#" \
             | bgzip -@ {threads} \
-        >> {output.bed}
+                >>{output.bed}
 
         printf '\nIndexing\n'
         tabix -f -p bed {output.bed}
@@ -176,9 +176,9 @@ rule helper_fdr_peaks_by_fire_elements:
         fire_tbi=rules.fire_sites_index.output.tbi,
     output:
         bed=temp("temp/{sm}/fire-peaks/{v}-{chrom}-fire-peaks.bed.gz"),
-    threads: 2
     conda:
         DEFAULT_ENV
+    threads: 2
     params:
         max_peak_fdr=MAX_PEAK_FDR,
         min_per_acc_peak=MIN_PER_ACC_PEAK,
@@ -186,36 +186,37 @@ rule helper_fdr_peaks_by_fire_elements:
         """
         HEADER=$(bgzip -cd {input.bed} | head -n 1 || true)
         NC=$(echo $HEADER | awk '{{print NF}}' || true)
-        FIRE_CT=$((NC+1))
-        FIRE_ST=$((NC+2))
-        FIRE_EN=$((NC+3))
-        FIRE_SIZE=$((NC+4))
-        FIRE_ID=$((NC+5))
+        FIRE_CT=$((NC + 1))
+        FIRE_ST=$((NC + 2))
+        FIRE_EN=$((NC + 3))
+        FIRE_SIZE=$((NC + 4))
+        FIRE_ID=$((NC + 5))
 
         OUT_HEADER=$(printf "$HEADER\\tpeak_chrom\\tpeak_start\\tpeak_end\\tFIRE_IDs\\tFIRE_size_mean\\tFIRE_size_ssd\\tFIRE_start_ssd\\tFIRE_end_ssd")
         echo $OUT_HEADER
-        
-        ( \
-            printf "$OUT_HEADER\\n"; \
+
+        (
+            printf "$OUT_HEADER\\n"
             tabix -h {input.bed} {wildcards.chrom} \
                 | bioawk -tc hdr '(NR==1)||($is_local_max=="true")' \
                 | csvtk filter -tT -C '$' -f "FDR<={params.max_peak_fdr}" \
                 | csvtk filter -tT -C '$' -f "fire_coverage>1" \
                 | bioawk -tc hdr '(NR==1)||(NF>0 && $fire_coverage/$coverage>={params.min_per_acc_peak})' \
                 | bedtools intersect -wa -wb -sorted -a - \
-                    -b <(tabix {input.fire} {wildcards.chrom} \
+                    -b <(
+                        tabix {input.fire} {wildcards.chrom} \
                             | cut -f 1-3 \
-                            | awk -v OFMT="%f" '{{print $0"\t"$3-$2"\t"NR}}' \
-                        ) \
+                            | awk -v OFMT="%f" '{{print $0"\t"$3-$2"\t"NR}}'
+                    ) \
                 | bedtools groupby -g 1-$NC \
                     -o first,median,median,collapse,mean,sstdev,sstdev,sstdev \
-                    -c $FIRE_CT,$FIRE_ST,$FIRE_EN,$FIRE_ID,$FIRE_SIZE,$FIRE_SIZE,$FIRE_ST,$FIRE_EN \
+                    -c $FIRE_CT,$FIRE_ST,$FIRE_EN,$FIRE_ID,$FIRE_SIZE,$FIRE_SIZE,$FIRE_ST,$FIRE_EN
         ) \
             | hck -f 1,$FIRE_ST,$FIRE_EN,2-$NC,$FIRE_SIZE- \
             | csvtk round -tT -C '$' -n 0 -f 2,3 \
             | bedtools sort -header -i - \
             | bgzip -@ {threads} \
-            > {output.bed}
+                >{output.bed}
         """
 
 
@@ -226,9 +227,9 @@ rule fdr_peaks_by_fire_elements_chromosome:
         maximum=rules.coverage.output.maximum,
     output:
         bed=temp("temp/{sm}/fire-peaks/{v}-grouped-{chrom}-fire-peaks.bed.gz"),
-    threads: 4
     conda:
         "../envs/python.yaml"
+    threads: 4
     params:
         script=workflow.source_path("../scripts/merge_fire_peaks.py"),
         min_frac_accessible=MIN_FRAC_ACCESSIBLE,
@@ -241,7 +242,7 @@ rule fdr_peaks_by_fire_elements_chromosome:
                 --min-cov $(cat {input.minimum}) \
                 --min-frac-accessible {params.min_frac_accessible} \
             | bgzip -@ {threads} \
-        > {output.bed}
+                >{output.bed}
         """
 
 
@@ -256,22 +257,22 @@ rule fire_peaks:
         fofn=temp("temp/{sm}/fire-peaks/{sm}-fire-{v}-peaks.fofn"),
         bed="results/{sm}/{sm}-fire-{v}-peaks.bed.gz",
         tbi="results/{sm}/{sm}-fire-{v}-peaks.bed.gz.tbi",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
         printf "\nMaking FOFN\n"
-        echo {input.beds} > {output.fofn}
+        echo {input.beds} >{output.fofn}
 
-        printf "\nMaking header\n"        
-        ((cat $(cat {output.fofn}) | bgzip -cd | grep "^#" | head -n 1) || true) \
-            | bgzip -@ {threads} > {output.bed}
+        printf "\nMaking header\n"
+        ( (cat $(cat {output.fofn}) | bgzip -cd | grep "^#" | head -n 1) || true) \
+            | bgzip -@ {threads} >{output.bed}
 
         printf "\nConcatenating\n"
         cat $(cat {output.fofn}) | bgzip -cd -@ {threads} | grep -v "^#" \
-            | bgzip -@ {threads} >> {output.bed}
-        
+            | bgzip -@ {threads} >>{output.bed}
+
         printf "\nIndexing\n"
         tabix -f -p bed {output.bed}
         """
@@ -296,22 +297,22 @@ rule wide_fire_peaks:
         bed3_as=workflow.source_path("../templates/bed3.as"),
     shell:
         """
-        ( \
-            bgzip -cd {input.bed}; \
+        (
+            bgzip -cd {input.bed}
             bioawk -tc hdr 'NR==1 || $FDR<={params.max_peak_fdr}' {input.track} \
-                | bioawk -tc hdr 'NR==1 || (NF>0 && $coverage>0 && $fire_coverage/$coverage>={params.min_frac_acc})' \
+                | bioawk -tc hdr 'NR==1 || (NF>0 && $coverage>0 && $fire_coverage/$coverage>={params.min_frac_acc})'
         ) \
             | cut -f 1-3 \
             | bedtools sort \
             | bedtools merge -d {params.nuc_size} \
             | bgzip -@ {threads} \
-        > {output.bed}
-        
+                >{output.bed}
+
         bgzip -cd -@ 16 {output.bed} \
             | bigtools bedtobigbed \
                 -s start -a {params.bed3_as} \
-                - {input.fai} {output.bb}        
-        
+                - {input.fai} {output.bb}
+
         tabix -p bed {output.bed}
         """
 
@@ -325,9 +326,9 @@ rule one_percent_fire_peaks:
         tbi="results/{sm}/additional-outputs-{v}/fire-peaks/one-percent-FDR/{sm}-fire-{v}-01-fire-peaks.bed.gz.tbi",
         wide="results/{sm}/additional-outputs-{v}/fire-peaks/one-percent-FDR/{sm}-fire-{v}-01-fire-wide-peaks.bed.gz",
         wtbi="results/{sm}/additional-outputs-{v}/fire-peaks/one-percent-FDR/{sm}-fire-{v}-01-fire-wide-peaks.bed.gz.tbi",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     params:
         nuc_size=config.get("nucleosome_size", 147),
     shell:
@@ -335,18 +336,18 @@ rule one_percent_fire_peaks:
         bgzip -cd {input.bed} \
             | csvtk filter -tT -C '$' -f "FDR<=0.01" \
             | bgzip -@ {threads} \
-            > {output.bed}
+                >{output.bed}
         tabix -f -p bed {output.bed}
 
-        ( \
-            bgzip -cd {output.bed}; \
-            bioawk -tc hdr '$FDR<=0.01' {input.track} \
+        (
+            bgzip -cd {output.bed}
+            bioawk -tc hdr '$FDR<=0.01' {input.track}
         ) \
             | cut -f 1-3 \
             | bedtools sort \
             | bedtools merge -d {params.nuc_size} \
             | bgzip -@ {threads} \
-        > {output.wide}
+                >{output.wide}
         tabix -f -p bed {output.wide}
         """
 
@@ -359,8 +360,8 @@ rule peaks_vs_percent:
             "results/{sm}/additional-outputs-{v}/figures/{sm}-fire-{v}-peaks-vs-percent.pdf",
             category="Peak calls",
         ),
-    threads: 4
     conda:
         "../envs/R.yaml"
+    threads: 4
     script:
         "../scripts/peaks-vs-percent.R"
diff --git a/workflow/rules/levio.smk b/workflow/rules/levio.smk
index e6d4b1b16..409303abc 100644
--- a/workflow/rules/levio.smk
+++ b/workflow/rules/levio.smk
@@ -41,12 +41,12 @@ rule leviosam2:
         lifted=temp("temp/{sm}/leviosam2/{sm}-{chrom}-committed.bam"),
         deferred=temp("temp/{sm}/leviosam2/{sm}-{chrom}-deferred.bam"),
         unliftable=temp("temp/{sm}/leviosam2/{sm}-{chrom}-unliftable.bam"),
+    conda:
+        DEFAULT_ENV
     threads: MAX_THREADS
     resources:
         mem_mb=MAX_THREADS * 4 * 1024,
         runtime=16 * 60,
-    conda:
-        DEFAULT_ENV
     params:
         # maximum number of CIGAR opts to change, also the max gap size that can be spanned
         G=config.get("levio_G", 100_000),
@@ -83,12 +83,12 @@ rule leviosam2_sorted:
         ref=REF,
     output:
         bam=temp("temp/{sm}/leviosam2/{sm}-{chrom}-sorted.bam"),
+    conda:
+        DEFAULT_ENV
     threads: SORT_THREADS
     resources:
         mem_mb=SORT_THREADS * 4 * 1024,
         runtime=16 * 60,
-    conda:
-        DEFAULT_ENV
     shell:
         """
         samtools sort {input.lifted} \
diff --git a/workflow/rules/stats.smk b/workflow/rules/stats.smk
index 366eca87e..f462c8592 100644
--- a/workflow/rules/stats.smk
+++ b/workflow/rules/stats.smk
@@ -9,19 +9,21 @@ rule clustering_vs_null:
         tmp=temp("temp/{sm}/tmp.pre.calls.bed"),
         null=temp("temp/{sm}/null.calls.bed"),
         bed="results/{sm}/clustering-vs-null.bed.gz",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     shell:
         """
-        bgzip -cd -@{threads} {input.bed} | cut -f 1-3 > {output.tmp}
-        bedtools shuffle -chrom -i {output.tmp} -g {input.fai} > {output.null}
+        bgzip -cd -@{threads} {input.bed} | cut -f 1-3 >{output.tmp}
+        bedtools shuffle -chrom -i {output.tmp} -g {input.fai} >{output.null}
 
-        ( bedtools genomecov -bg -i {output.tmp} -g {input.fai} | sed 's/$/\\tReal/g' ; \
-          bedtools genomecov -bg -i {output.null} -g {input.fai} | sed 's/$/\\tNull/g' ) \
+        (
+            bedtools genomecov -bg -i {output.tmp} -g {input.fai} | sed 's/$/\\tReal/g'
+            bedtools genomecov -bg -i {output.null} -g {input.fai} | sed 's/$/\\tNull/g'
+        ) \
             | bedtools sort \
             | bgzip -@ {threads} \
-        > {output.bed}
+                >{output.bed}
         """
 
 
@@ -33,20 +35,20 @@ rule fires_in_peaks:
     output:
         tmp=temp("temp/{sm}/tmp.FIREs-{v}-in-peaks.bed"),
         txt="results/{sm}/additional-outputs-{v}/fire-peaks/{sm}-{v}-fires-in-peaks.txt",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     params:
         script=workflow.source_path("../scripts/percent-in-clusters.sh"),
     shell:
         """
-        bedtools intersect -sorted -a {input.fire} -b {input.exclude} -v > {output.tmp}
+        bedtools intersect -sorted -a {input.fire} -b {input.exclude} -v >{output.tmp}
 
-        echo "Total # of FIREs within normal coverage regions" >> {output.txt}
-        wc -l {output.tmp} >> {output.txt}
+        echo "Total # of FIREs within normal coverage regions" >>{output.txt}
+        wc -l {output.tmp} >>{output.txt}
 
-        echo "# of FIREs within peaks" >> {output.txt}
-        bedtools intersect -sorted -u -a {input.fire} -b {input.peaks} | wc -l >> {output.txt} 
+        echo "# of FIREs within peaks" >>{output.txt}
+        bedtools intersect -sorted -u -a {input.fire} -b {input.peaks} | wc -l >>{output.txt}
         """
 
 
@@ -78,8 +80,8 @@ rule hap_differences:
         ),
         bed="results/{sm}/{sm}-fire-{v}-hap-differences.bed.gz",
         bed9=temp("temp/{sm}/hap1-vs-hap2/FIRE-{v}.hap.differences.bed9"),
-    threads: 4
     conda:
         "../envs/R.yaml"
+    threads: 4
     script:
         "../scripts/hap-diffs.R"
diff --git a/workflow/rules/track-hub.smk b/workflow/rules/track-hub.smk
index 422fde256..1be0fcb55 100644
--- a/workflow/rules/track-hub.smk
+++ b/workflow/rules/track-hub.smk
@@ -5,9 +5,9 @@ rule percent_accessible:
     output:
         tmp=temp("temp/{sm}/{hp}/{v}-percent.accessible.bed"),
         bw="results/{sm}/trackHub-{v}/bw/{hp}.percent.accessible.bw",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     resources:
         mem_mb=get_mem_mb,
     params:
@@ -19,18 +19,17 @@ rule percent_accessible:
         bgzip -cd {input.bed} \
             | bioawk -tc hdr '$coverage{params.suffix}>0' \
             | bioawk -tc hdr \
-              'NR>1{{print $1,$2,$3,100*$fire_coverage{params.suffix}/$coverage{params.suffix}}}' \
-        > {output.tmp}
+                'NR>1{{print $1,$2,$3,100*$fire_coverage{params.suffix}/$coverage{params.suffix}}}' \
+                >{output.tmp}
 
         # add fake if file is empty
         if [[ -s {output.tmp} ]]; then
             echo "File is not empty"
         else
             echo "File is empty"
-            printf "{params.chrom}\t0\t1\t0\\n" > {output.tmp}
+            printf "{params.chrom}\t0\t1\t0\\n" >{output.tmp}
         fi
 
-
         bigtools bedgraphtobigwig \
             --nzooms {params.nzooms} -s start \
             {output.tmp} {input.fai} {output.bw}
@@ -65,9 +64,9 @@ rule fdr_track_to_bw:
         fai=ancient(FAI),
     output:
         bw="results/{sm}/trackHub-{v}/bw/{col}.bw",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     params:
         nzooms=NZOOMS,
     shell:
@@ -86,9 +85,9 @@ rule fire_peaks_bb:
         fai=ancient(FAI),
     output:
         bb="results/{sm}/trackHub-{v}/bb/fire-peaks.bb",
-    threads: 4
     conda:
         DEFAULT_ENV
+    threads: 4
     params:
         bedfmt=workflow.source_path("../templates/fire_peak.as"),
     shell:
@@ -109,19 +108,19 @@ rule hap_differences_track:
         fai=ancient(FAI),
     output:
         bb="results/{sm}/trackHub-{v}/bb/hap_differences.bb",
+    conda:
+        DEFAULT_ENV
     threads: 4
     resources:
         mem_mb=get_mem_mb,
-    conda:
-        DEFAULT_ENV
     params:
         chrom=get_chroms()[0],
         bed9_as=workflow.source_path("../templates/bed9.as"),
     shell:
         """
-        ( \
-            printf "{params.chrom}\t0\t1\tfake\t100\t+\t0\t1\t230,230,230\\n"; \
-            bedtools sort -i {input.bed9} \
+        (
+            printf "{params.chrom}\t0\t1\tfake\t100\t+\t0\t1\t230,230,230\\n"
+            bedtools sort -i {input.bed9}
         ) \
             | bigtools bedtobigbed \
                 -s start -a {params.bed9_as} \
@@ -135,11 +134,11 @@ rule trackhub:
     output:
         hub="results/{sm}/trackHub-{v}/hub.txt",
         description="results/{sm}/trackHub-{v}/fire-description.html",
-    resources:
-        load=get_load,
-    threads: 4
     conda:
         "../envs/python.yaml"
+    threads: 4
+    resources:
+        load=get_load,
     params:
         ref=REF_NAME,
         script=workflow.source_path("../scripts/trackhub.py"),
@@ -147,9 +146,9 @@ rule trackhub:
     shell:
         """
         python {params.script} -v 2 \
-          --trackhub-dir results/{wildcards.sm}/trackHub-{wildcards.v} \
-          --reference {params.ref} \
-          --sample {wildcards.sm} \
-          --average-coverage $(cat {input.cov}) 
+            --trackhub-dir results/{wildcards.sm}/trackHub-{wildcards.v} \
+            --reference {params.ref} \
+            --sample {wildcards.sm} \
+            --average-coverage $(cat {input.cov})
         cp {params.description} {output.description}
         """

From cb299846312f8d843261f263ec55d4c217bc3c46 Mon Sep 17 00:00:00 2001
From: "Mitchell R. Vollger" <mvollger@gmail.com>
Date: Fri, 14 Aug 2026 14:08:01 -0700
Subject: [PATCH 2/4] ci: run checks on pull requests against any base branch

Stacked pull requests have a non-main base and previously got no checks.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
---
 .github/workflows/main.yml | 1 -
 1 file changed, 1 deletion(-)

diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index 037f5d4e8..f5ba8c199 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -4,7 +4,6 @@ on:
   push:
     branches: [main]
   pull_request:
-    branches: [main]
 
 jobs:
   # pixi test

From ddddc8067931bb5d4374b0fbb2c10408e036665f Mon Sep 17 00:00:00 2001
From: "Mitchell R. Vollger" <mvollger@gmail.com>
Date: Fri, 14 Aug 2026 14:08:32 -0700
Subject: [PATCH 3/4] refactor: remove the unused leviosam2 chain-lift rules

Nothing consumes the leviosam2 outputs, chain mode has been broken since
the fire rule stopped producing a bam (its input referenced a bam output
that no longer exists), and the rules carried further latent bugs (index
written to a hardcoded shared path instead of the declared output).
Delete the rules and the chain/levio_exe config options.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
---
 workflow/Snakefile       |  10 ----
 workflow/rules/levio.smk | 102 ---------------------------------------
 2 files changed, 112 deletions(-)
 delete mode 100644 workflow/rules/levio.smk

diff --git a/workflow/Snakefile b/workflow/Snakefile
index 967e69371..4cb4ca372 100644
--- a/workflow/Snakefile
+++ b/workflow/Snakefile
@@ -74,11 +74,6 @@ bw_types = ["log_FDR"]  # "score", "FDR",
 el_types = ["fire", "linker", "nucleosome"]
 
 
-# DSA options
-DSA_CHAIN = config.get("chain", None)
-DSA = DSA_CHAIN is not None
-LEVIO_EXE = config.get("levio_exe", "leviosam2")
-
 # developer options
 FT_EXE = config.get("ft_exe", "ft")
 if FT_EXE != "ft":
@@ -104,11 +99,6 @@ include: "rules/decorated-reads.smk"
 include: "rules/track-hub.smk"
 
 
-if DSA:
-
-    include: "rules/levio.smk"
-
-
 wildcard_constraints:
     chrom="|".join(get_chroms()),
     call="|".join(["msp", "m6a"]),
diff --git a/workflow/rules/levio.smk b/workflow/rules/levio.smk
deleted file mode 100644
index 409303abc..000000000
--- a/workflow/rules/levio.smk
+++ /dev/null
@@ -1,102 +0,0 @@
-#
-# Index the chain file for leviosam2.
-#
-# Input here is a chain file that defines the alignment between the DSA and the reference genome
-# at a contig level (>100 kbp of alignment).
-#
-# The output is a special leviosam2 index file that is used to lift over the alignments from the DSA to the reference genome.
-#
-rule leviosam2_index:
-    input:
-        chain=DSA_CHAIN,
-        fai=FAI,
-    output:
-        index=temp("temp/{sm}/leviosam2-index/index.clft"),
-    conda:
-        DEFAULT_ENV
-    threads: 1
-    resources:
-        mem_mb=64 * 1024,
-        runtime=16 * 60,
-    shell:
-        """
-        {LEVIO_EXE} index \
-            -p results/leviosam2-index/index \
-            -c {input.chain} \
-            -F {input.fai}
-        """
-
-
-#
-# Lift over the alignments from the DSA to the reference genome using the chain file / leviosam2 index.
-#
-# This is not a realignment, but a lift over of the reads from the DSA to the reference genome.
-#
-rule leviosam2:
-    input:
-        bam=rules.fire.output.bam,
-        levio_index=rules.leviosam2_index.output.index,
-        ref=REF,
-    output:
-        lifted=temp("temp/{sm}/leviosam2/{sm}-{chrom}-committed.bam"),
-        deferred=temp("temp/{sm}/leviosam2/{sm}-{chrom}-deferred.bam"),
-        unliftable=temp("temp/{sm}/leviosam2/{sm}-{chrom}-unliftable.bam"),
-    conda:
-        DEFAULT_ENV
-    threads: MAX_THREADS
-    resources:
-        mem_mb=MAX_THREADS * 4 * 1024,
-        runtime=16 * 60,
-    params:
-        # maximum number of CIGAR opts to change, also the max gap size that can be spanned
-        G=config.get("levio_G", 100_000),
-        # Using -S clipped_frac 0.05 means when a read has >5% clipped bases, it is deferred. A lower value is more stringent (by deferring more reads).
-        # aln_score is the minumum score before the alignment is lifted over
-        S=config.get(
-            "levio_S",
-            f"-S mapq:0 -S hdist:{100_000} -S isize:{100_000} -S clipped_frac:0.95 -S aln_score:100",
-        ),
-        # number of reads per thread
-        T=config.get("levio_T", 4 * 256),
-    shell:
-        """
-        PRE="temp/{wildcards.sm}/leviosam2/{wildcards.sm}-{wildcards.chrom}"
-        {LEVIO_EXE} lift -t {threads} -a {input.cram} \
-            -T {params.T} -G {params.G} {params.S} \
-            -C {input.levio_index} -p $PRE -f {input.ref} -m -O bam
-        """
-
-
-#
-# This step sorted the leviosam2 output and fixes some tags in the CRAM file.
-#
-# Specifically, the MAPQ is reset to 60 for all reads that were previously aligned to the DSA.
-# And the XS tag is set to zero for all reads that were aligned to the DSA.
-# This is a hueristic that we may need to return to in the future.
-#
-# Other tags and fields like CIGAR, bitflags, and MD are correctly updated by
-# leviosam2 during liftover.
-#
-rule leviosam2_sorted:
-    input:
-        lifted=rules.leviosam2.output.lifted,
-        ref=REF,
-    output:
-        bam=temp("temp/{sm}/leviosam2/{sm}-{chrom}-sorted.bam"),
-    conda:
-        DEFAULT_ENV
-    threads: SORT_THREADS
-    resources:
-        mem_mb=SORT_THREADS * 4 * 1024,
-        runtime=16 * 60,
-    shell:
-        """
-        samtools sort {input.lifted} \
-            -@ {threads} -m 3G \
-            -o {output.bam}
-        """
-
-
-# params:
-# reset_mapq=workflow.source_path("../scripts/reset-mapq.py"),
-# python {params.reset_mapq} -t {threads} {input.lifted} \

From 32b9e5d745250caf0bb2e1d9968c43bc0cc7438b Mon Sep 17 00:00:00 2001
From: "Mitchell R. Vollger" <mvollger@gmail.com>
Date: Fri, 14 Aug 2026 14:08:52 -0700
Subject: [PATCH 4/4] feat: per-sample reference support via manifest
 ref/ref_name columns

The manifest can now carry ref and ref_name columns so each sample can
use its own reference; the two-column manifest with config-level values
keeps working, and a '.' cell falls back to the config value. Chromosome
names and lengths come only from each sample's BAM header, in header
order (the sort order of the data), and a generated per-sample genome
file replaces the fai for all bedtools and bigtools sizes arguments
(with a parse-time staleness guard, so a leftover genome file can never
disagree with its bam header). All reference and chromosome lookups are
functions of the sm wildcard. Includes the fixes from the multi-agent
code review (manifest cells parsed as literal text, config-specific
errors for null ref values, sample-attributed test assertions).
Supersedes #44.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
---
 .github/workflows/main.yml               |   2 +
 CONTRIBUTING.md                          |  13 +-
 config/README.md                         |  18 +-
 pixi.lock                                | 259 ++++++++++++++---------
 pixi.toml                                |  14 +-
 tests/README.md                          |  44 ++++
 tests/config/err-config-ref-only.yaml    |   2 +
 tests/config/err-dup-sample.tbl          |   3 +
 tests/config/err-dup-sample.yaml         |   3 +
 tests/config/err-keepchrs.yaml           |   4 +
 tests/config/err-missing-bam.tbl         |   2 +
 tests/config/err-missing-bam.yaml        |   1 +
 tests/config/err-missing-fai.tbl         |   2 +
 tests/config/err-missing-fai.yaml        |   1 +
 tests/config/err-missing-refpath.tbl     |   2 +
 tests/config/err-missing-refpath.yaml    |   1 +
 tests/config/err-nan-cell.tbl            |   3 +
 tests/config/err-nan-cell.yaml           |   1 +
 tests/config/err-no-ref-anywhere.yaml    |   1 +
 tests/config/err-ref-col-only.tbl        |   2 +
 tests/config/err-ref-col-only.yaml       |   1 +
 tests/config/err-sentinel-no-config.tbl  |   2 +
 tests/config/err-sentinel-no-config.yaml |   1 +
 tests/config/extra-exclude.bed           |   1 +
 tests/config/multi.tbl                   |   5 +
 tests/config/multi.yaml                  |   5 +
 tests/config/no-excludes.yaml            |   1 +
 tests/config/no-fai.fa                   |   2 +
 tests/config/override-info.yaml          |   3 +
 tests/config/sentinel.tbl                |   2 +
 tests/config/sentinel.yaml               |   3 +
 tests/config/single-4col.tbl             |   2 +
 tests/config/single-4col.yaml            |   1 +
 tests/config/two-sample.tbl              |   3 +
 tests/config/two-sample.yaml             |   3 +
 tests/dry-run-checks.sh                  | 153 +++++++++++++
 tests/generate-test-bams.py              |  69 ++++++
 tests/verify-multi.sh                    |  52 +++++
 tests/verify-results.sh                  |  28 +++
 workflow/Snakefile                       |  17 +-
 workflow/envs/runner.yaml                |   5 +-
 workflow/rules/apply-model.smk           |   4 +-
 workflow/rules/common.smk                | 253 +++++++++++++++++-----
 workflow/rules/coverages.smk             |  40 +++-
 workflow/rules/decorated-reads.smk       |  14 +-
 workflow/rules/fire-peaks.smk            |  14 +-
 workflow/rules/stats.smk                 |   8 +-
 workflow/rules/track-hub.smk             |  26 +--
 48 files changed, 883 insertions(+), 213 deletions(-)
 create mode 100644 tests/README.md
 create mode 100644 tests/config/err-config-ref-only.yaml
 create mode 100644 tests/config/err-dup-sample.tbl
 create mode 100644 tests/config/err-dup-sample.yaml
 create mode 100644 tests/config/err-keepchrs.yaml
 create mode 100644 tests/config/err-missing-bam.tbl
 create mode 100644 tests/config/err-missing-bam.yaml
 create mode 100644 tests/config/err-missing-fai.tbl
 create mode 100644 tests/config/err-missing-fai.yaml
 create mode 100644 tests/config/err-missing-refpath.tbl
 create mode 100644 tests/config/err-missing-refpath.yaml
 create mode 100644 tests/config/err-nan-cell.tbl
 create mode 100644 tests/config/err-nan-cell.yaml
 create mode 100644 tests/config/err-no-ref-anywhere.yaml
 create mode 100644 tests/config/err-ref-col-only.tbl
 create mode 100644 tests/config/err-ref-col-only.yaml
 create mode 100644 tests/config/err-sentinel-no-config.tbl
 create mode 100644 tests/config/err-sentinel-no-config.yaml
 create mode 100644 tests/config/extra-exclude.bed
 create mode 100644 tests/config/multi.tbl
 create mode 100644 tests/config/multi.yaml
 create mode 100644 tests/config/no-excludes.yaml
 create mode 100644 tests/config/no-fai.fa
 create mode 100644 tests/config/override-info.yaml
 create mode 100644 tests/config/sentinel.tbl
 create mode 100644 tests/config/sentinel.yaml
 create mode 100644 tests/config/single-4col.tbl
 create mode 100644 tests/config/single-4col.yaml
 create mode 100644 tests/config/two-sample.tbl
 create mode 100644 tests/config/two-sample.yaml
 create mode 100644 tests/dry-run-checks.sh
 create mode 100644 tests/generate-test-bams.py
 create mode 100644 tests/verify-multi.sh
 create mode 100644 tests/verify-results.sh

diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml
index f5ba8c199..abf9d5d36 100644
--- a/.github/workflows/main.yml
+++ b/.github/workflows/main.yml
@@ -15,4 +15,6 @@ jobs:
         with:
           pixi-version: v0.70.2
           cache: true
+      - run: pixi run test-dry
       - run: pixi run test
+      - run: pixi run test-verify
diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md
index 964441543..5008f9972 100644
--- a/CONTRIBUTING.md
+++ b/CONTRIBUTING.md
@@ -1,13 +1,22 @@
 # Contributing to FIRE
 
-PRs are welcome! Please run the test suite and the formatters before
-submitting:
+PRs are welcome! Please run the validation checks, the test suite, and the
+formatters before submitting:
 
 ```bash
+pixi run test-dry
 pixi run test
+pixi run test-verify
 pixi run fmt
 ```
 
+If your change touches manifest or reference handling, also run the
+multi-sample test:
+
+```bash
+pixi run test-multi
+```
+
 ## Conventional commits and releases
 
 Releases are automated with
diff --git a/config/README.md b/config/README.md
index 01fe3a763..3a7de2e3d 100644
--- a/config/README.md
+++ b/config/README.md
@@ -2,15 +2,27 @@
 see `config.yaml` for an example.
 
 ## Required input options
-Reference name, this is the name that will be used as the genome name in the UCSC track hub, so be sure to use a valid UCSC genome name when possible. 
+Reference name, this is the name that will be used as the genome name in the UCSC track hub, so be sure to use a valid UCSC genome name when possible. A reference name of `hg38` or `GRCh38` also turns on the default `excludes` files (see below).
 ```
 ref_name: hg38
 ```
-Reference `fasta` file:
+Reference `fasta` file (a `.fai` index must exist beside it):
 ```
 ref: /path/to/hg38.fa
 ```
-Manifest of input sample(s), must have two white-space separated columns: sample name (`sample`) and input bam file path (`bam`). See `config.tbl` for an example. The `bam` file must be indexed and aligned to the reference genome in the `ref` option.
+Manifest of input sample(s), white-space separated with a header row. See `config.tbl` for an example. The two-column form gives a sample name (`sample`) and an input bam path (`bam`); every sample then uses the `ref` and `ref_name` from this config file:
+```
+sample	bam
+sample1	/path/to/sample1.bam
+```
+The manifest can instead carry per-sample references with two more columns, `ref` and `ref_name`. Add both columns together. Every cell must be non-empty. A cell containing `.` uses the config value for that row. Filled cells override the config values:
+```
+sample	bam	ref	ref_name
+sample1	/path/to/sample1.bam	/path/to/hg38.fa	hg38
+sample2	/path/to/sample2.bam	/path/to/chm13.fa	GCA_009914755.4
+sample3	/path/to/sample3.bam	.	.
+```
+Each `bam` file must be indexed and aligned to its reference genome. FIRE reads the chromosome names and lengths from the bam header, in header order, not from the fasta, so the fasta can contain extra contigs that the bam does not use. Because of this, FIRE opens every manifest bam when it starts, for every command including dry-runs — keep the input bams readable for the lifetime of the results.
 ```
 manifest: config/config.tbl
 ```
diff --git a/pixi.lock b/pixi.lock
index 060ac2f1c..d74bcf4a6 100644
--- a/pixi.lock
+++ b/pixi.lock
@@ -18,7 +18,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/noarch/boltons-25.0.0-pyhd8ed1ab_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/brotli-python-1.2.0-py312h67db365_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/bzip2-1.0.8-hda65f42_8.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/c-ares-1.34.5-hb9d3cd8_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/c-ares-1.34.8-h280c20c_1.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/ca-certificates-2025.11.12-hbd8a1cb_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/certifi-2025.11.12-pyhd8ed1ab_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/cffi-2.0.0-py312h460c074_1.conda
@@ -53,7 +53,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/noarch/hpack-4.1.0-pyhd8ed1ab_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/humanfriendly-10.0-pyh707e725_8.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/hyperframe-6.1.0-pyhd8ed1ab_0.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/icu-75.1-he02047a_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/icu-78.3-py310h44b86e0_2.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/idna-3.11-pyhd8ed1ab_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/immutables-0.21-py312h4c3975b_2.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/iniconfig-2.3.0-pyhd8ed1ab_0.conda
@@ -64,12 +64,13 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/noarch/jsonschema-specifications-2025.9.1-pyhcf101f3_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/jupyter_core-5.9.1-pyhc90fa1f_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/keyutils-1.6.3-hb9d3cd8_0.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/krb5-1.21.3-h659f571_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/krb5-1.22.2-hbde042b_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/ld_impl_linux-64-2.44-h1aa0949_5.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libarchive-3.8.2-gpl_h7be2006_100.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libblas-3.9.0-38_h4a7cf45_openblas.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libcblas-3.9.0-38_h0358290_openblas.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libcurl-8.17.0-h4e3cde8_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libcurl-8.21.0-heca4667_4.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libdeflate-1.25-hd45a770_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libedit-3.1.20250104-pl5321h7949ede_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libev-4.33-hd590300_2.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libexpat-2.7.1-hecca717_0.conda
@@ -82,22 +83,22 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libiconv-1.18-h3b78370_2.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/liblapack-3.9.0-38_h47877c9_openblas.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/liblapacke-3.9.0-38_h6ae95b6_openblas.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/liblzma-5.8.1-hb9d3cd8_2.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/liblzma-5.8.3-hb03c661_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libmamba-2.3.3-h09abcb8_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libmambapy-2.3.3-py312h9b917ee_1.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libnghttp2-1.67.0-had1ee68_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libnghttp2-1.68.1-h877daf1_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libnsl-2.0.1-hb9d3cd8_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libopenblas-0.3.30-pthreads_h94d23a6_3.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libpsl-0.23.1-hf670292_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libsolv-0.7.35-h9463b59_0.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libsqlite-3.51.0-hee844dc_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libsqlite-3.53.4-hf4e2dac_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libssh2-1.11.1-hcf80075_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libstdcxx-15.2.0-h8f9b012_7.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libstdcxx-ng-15.2.0-h4852527_7.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libuuid-2.41.2-he9a06e4_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/libxcrypt-4.4.36-hd590300_1.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-16-2.15.1-ha9997c6_0.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-2.15.1-h26afc86_0.conda
-      - conda: https://conda.anaconda.org/conda-forge/linux-64/libzlib-1.3.1-hb9d3cd8_2.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-16-2.15.3-hca6bf5a_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-2.15.3-h49c6c72_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/linux-64/libzlib-1.3.2-h25fd6f3_3.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/logmuse-0.2.8-pyhd8ed1ab_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/lz4-c-1.10.0-h5888daf_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/lzo-2.10-h280c20c_1002.conda
@@ -128,6 +129,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/pydantic-core-2.41.5-py312h868fb18_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pygments-2.19.2-pyhd8ed1ab_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pyparsing-3.2.5-pyhcf101f3_0.conda
+      - conda: https://conda.anaconda.org/bioconda/linux-64/pysam-0.24.0-py312hf5ad864_1.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pysocks-1.7.1-pyha55dd90_7.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pytest-9.0.0-pyhcf101f3_0.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-64/python-3.12.12-hd63d673_1_cpython.conda
@@ -271,6 +273,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libcblas-3.11.0-9_hb0561ab_openblas.conda
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libcurl-8.21.0-hf618e03_4.conda
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libcxx-22.1.8-h55c6f16_0.conda
+      - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libdeflate-1.25-he7e0567_1.conda
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libedit-3.1.20250104-pl5321h26f1114_1.conda
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libev-4.33-h1a92334_3.conda
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libexpat-2.8.1-hf6b4638_1.conda
@@ -329,6 +332,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/noarch/pygments-2.20.0-pyhd8ed1ab_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pyparsing-3.3.2-pyhcf101f3_0.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pyproject_hooks-1.2.0-pyhd8ed1ab_1.conda
+      - conda: https://conda.anaconda.org/bioconda/osx-arm64/pysam-0.24.0-py312h12d0683_1.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/pysocks-1.7.1-pyha55dd90_7.conda
       - conda: https://conda.anaconda.org/conda-forge/osx-arm64/python-3.12.13-hd1323d7_1_cpython.conda
       - conda: https://conda.anaconda.org/conda-forge/noarch/python-build-1.5.0-pyhc364b38_0.conda
@@ -667,16 +671,18 @@ packages:
   license_family: BSD
   size: 124965
   timestamp: 1785906749812
-- conda: https://conda.anaconda.org/conda-forge/linux-64/c-ares-1.34.5-hb9d3cd8_0.conda
-  sha256: f8003bef369f57396593ccd03d08a8e21966157269426f71e943f96e4b579aeb
-  md5: f7f0d6cc2dc986d42ac2689ec88192be
+- conda: https://conda.anaconda.org/conda-forge/linux-64/c-ares-1.34.8-h280c20c_1.conda
+  sha256: 5139b6afbfaca91c47104ba9a6a40f81211e1c9200e96b73ce3a56f0ca1902f4
+  md5: 2cef891b791040aab83e218c7d137679
   depends:
   - __glibc >=2.17,<3.0.a0
-  - libgcc >=13
+  - libgcc >=14
+  constrains:
+  - c-ares-static <0a0
   license: MIT
   license_family: MIT
-  size: 206884
-  timestamp: 1744127994291
+  size: 226755
+  timestamp: 1786116641939
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/c-ares-1.34.8-h1a92334_1.conda
   sha256: 104b41473845649101ba8ebf8221c7431256465d34ce380b10e9a90558ed33ac
   md5: 7d9390a4d4b43f91652823d870f68065
@@ -1621,17 +1627,17 @@ packages:
   license_family: MIT
   size: 17397
   timestamp: 1737618427549
-- conda: https://conda.anaconda.org/conda-forge/linux-64/icu-75.1-he02047a_0.conda
-  sha256: 71e750d509f5fa3421087ba88ef9a7b9be11c53174af3aa4d06aff4c18b38e8e
-  md5: 8b189310083baabfb622af68fd9d3ae3
+- conda: https://conda.anaconda.org/conda-forge/linux-64/icu-78.3-py310h44b86e0_2.conda
+  sha256: 9f07834f0c546ab14d885ce0366285f61f44e326c0edd1fc63b8294e113ae432
+  md5: 72a381cbad04f24b1c2a43ef707f45b4
   depends:
   - __glibc >=2.17,<3.0.a0
-  - libgcc-ng >=12
-  - libstdcxx-ng >=12
+  - libstdcxx >=14
+  - libgcc >=14
   license: MIT
   license_family: MIT
-  size: 12129203
-  timestamp: 1720853576813
+  size: 14459115
+  timestamp: 1786545741408
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/icu-78.3-py310h579977c_2.conda
   sha256: 6cdb5dee54c72e56ab189fb3ad33cb28533553d42590e7e831160248f4416a43
   md5: a5efc0b42bb8b42e97d0a29ae3e3c187
@@ -1819,20 +1825,21 @@ packages:
   license: LGPL-2.1-or-later
   size: 134088
   timestamp: 1754905959823
-- conda: https://conda.anaconda.org/conda-forge/linux-64/krb5-1.21.3-h659f571_0.conda
-  sha256: 99df692f7a8a5c27cd14b5fb1374ee55e756631b9c3d659ed3ee60830249b238
-  md5: 3f43953b7d3fb3aaa1d0d0723d91e368
+- conda: https://conda.anaconda.org/conda-forge/linux-64/krb5-1.22.2-hbde042b_1.conda
+  sha256: 9b07046870772f28740e3f6149f09ff222843733087a33c5540b169c6289652d
+  md5: 54157a1c8c0bb70f62dd0b17fba7e7f2
   depends:
-  - keyutils >=1.6.1,<2.0a0
-  - libedit >=3.1.20191231,<3.2.0a0
-  - libedit >=3.1.20191231,<4.0a0
-  - libgcc-ng >=12
-  - libstdcxx-ng >=12
-  - openssl >=3.3.1,<4.0a0
+  - __glibc >=2.17,<3.0.a0
+  - keyutils >=1.6.3,<2.0a0
+  - libedit >=3.1.20250104,<3.2.0a0
+  - libedit >=3.1.20250104,<4.0a0
+  - libgcc >=14
+  - libstdcxx >=14
+  - openssl >=3.5.7,<4.0a0
   license: MIT
   license_family: MIT
-  size: 1370023
-  timestamp: 1719463201255
+  size: 1388990
+  timestamp: 1781859420533
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/krb5-1.22.2-hfd3d5f3_1.conda
   sha256: c740e4a2e7247776a9883158fdab50ae0732c8f67f96d8f1db8ad9da5e0b5222
   md5: 8780f41b013d19219faef9c82260744b
@@ -1957,22 +1964,23 @@ packages:
   license_family: BSD
   size: 18110
   timestamp: 1786058893756
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libcurl-8.17.0-h4e3cde8_0.conda
-  sha256: 100e29ca864c32af15a5cc354f502d07b2600218740fdf2439fa7d66b50b3529
-  md5: 01e149d4a53185622dc2e788281961f2
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libcurl-8.21.0-heca4667_4.conda
+  sha256: aff8ef75636d0825ce34911e0bef2f26816e7afd090a9dcb4b9bacab75cbf584
+  md5: 3f2fd5617cfacac49c85f6dc63842ea1
   depends:
   - __glibc >=2.17,<3.0.a0
-  - krb5 >=1.21.3,<1.22.0a0
+  - krb5 >=1.22.2,<1.23.0a0
   - libgcc >=14
-  - libnghttp2 >=1.67.0,<2.0a0
+  - libnghttp2 >=1.68.1,<2.0a0
+  - libpsl >=0.23.0,<0.24.0a0
   - libssh2 >=1.11.1,<2.0a0
-  - libzlib >=1.3.1,<2.0a0
-  - openssl >=3.5.4,<4.0a0
+  - libzlib >=1.3.2,<2.0a0
+  - openssl >=3.5.7,<4.0a0
   - zstd >=1.5.7,<1.6.0a0
   license: curl
   license_family: MIT
-  size: 460366
-  timestamp: 1762333743748
+  size: 480565
+  timestamp: 1785500108494
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libcurl-8.21.0-hf618e03_4.conda
   sha256: d25be36712d7f854a5f93513b9fbf1b0ee3219975b7fef4a3ee071b021b10167
   md5: 66cf9c5003ee81ecdb0dc9f9df17bbe5
@@ -1998,6 +2006,25 @@ packages:
   license_family: Apache
   size: 569349
   timestamp: 1781670209146
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libdeflate-1.25-hd45a770_1.conda
+  sha256: 82e134c8a08b1eed9a2ed8ab578b89aa1730dcde3dea8dd87645ed0637878e54
+  md5: 40f9b31aa9cf007789867df0decd0492
+  depends:
+  - __glibc >=2.17,<3.0.a0
+  - libgcc >=14
+  license: MIT
+  license_family: MIT
+  size: 73710
+  timestamp: 1785908694612
+- conda: https://conda.anaconda.org/conda-forge/osx-arm64/libdeflate-1.25-he7e0567_1.conda
+  sha256: d896f4aa4ce4c590c2838678cb1917356fdb461d2a189991c0280c818c362172
+  md5: 78650d671cb56909bb3e5c13bce310f9
+  depends:
+  - __osx >=11.0
+  license: MIT
+  license_family: MIT
+  size: 55727
+  timestamp: 1785909153744
 - conda: https://conda.anaconda.org/conda-forge/linux-64/libedit-3.1.20250104-pl5321h7949ede_0.conda
   sha256: d789471216e7aba3c184cd054ed61ce3f6dac6f87a50ec69291b9297f8c18724
   md5: c277e0a4d549b03ac1e9d6cbbe3d017b
@@ -2242,17 +2269,17 @@ packages:
   license_family: BSD
   size: 18174
   timestamp: 1786058909286
-- conda: https://conda.anaconda.org/conda-forge/linux-64/liblzma-5.8.1-hb9d3cd8_2.conda
-  sha256: f2591c0069447bbe28d4d696b7fcb0c5bd0b4ac582769b89addbcf26fb3430d8
-  md5: 1a580f7796c7bf6393fddb8bbbde58dc
+- conda: https://conda.anaconda.org/conda-forge/linux-64/liblzma-5.8.3-hb03c661_1.conda
+  sha256: 9787df8c22a59c9a70d3e5a10db9ad663485e75e9ccc3f09bd092cb7b95e0dab
+  md5: 1390b7c5ac0b1d8e447bc5efa6d3c8c2
   depends:
   - __glibc >=2.17,<3.0.a0
-  - libgcc >=13
+  - libgcc >=14
   constrains:
-  - xz 5.8.1.*
+  - xz 5.8.3.*
   license: 0BSD
-  size: 112894
-  timestamp: 1749230047870
+  size: 112995
+  timestamp: 1786348617826
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/liblzma-5.8.3-h8088a28_1.conda
   sha256: 23d0630046a3e8b164d8f80f2b74ed2605af2e7050ab9913018056402fae4311
   md5: 8ab10323068b107661a4b9a4af84f3b5
@@ -2376,22 +2403,22 @@ packages:
   license: BSL-1.0
   size: 39302
   timestamp: 1786189396450
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libnghttp2-1.67.0-had1ee68_0.conda
-  sha256: a4a7dab8db4dc81c736e9a9b42bdfd97b087816e029e221380511960ac46c690
-  md5: b499ce4b026493a13774bcf0f4c33849
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libnghttp2-1.68.1-h877daf1_0.conda
+  sha256: 663444d77a42f2265f54fb8b48c5450bfff4388d9c0f8253dd7855f0d993153f
+  md5: 2a45e7f8af083626f009645a6481f12d
   depends:
   - __glibc >=2.17,<3.0.a0
-  - c-ares >=1.34.5,<2.0a0
+  - c-ares >=1.34.6,<2.0a0
   - libev >=4.33,<4.34.0a0
   - libev >=4.33,<5.0a0
   - libgcc >=14
   - libstdcxx >=14
   - libzlib >=1.3.1,<2.0a0
-  - openssl >=3.5.2,<4.0a0
+  - openssl >=3.5.5,<4.0a0
   license: MIT
   license_family: MIT
-  size: 666600
-  timestamp: 1756834976695
+  size: 663344
+  timestamp: 1773854035739
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libnghttp2-1.68.1-h8f3e76b_0.conda
   sha256: 2bc7bc3978066f2c274ebcbf711850cc9ab92e023e433b9631958a098d11e10a
   md5: 6ea18834adbc3b33df9bd9fb45eaf95b
@@ -2445,6 +2472,18 @@ packages:
   license_family: BSD
   size: 4318474
   timestamp: 1784288246205
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libpsl-0.23.1-hf670292_0.conda
+  sha256: 2e6405feb59e3f40a5a4641dfa73afda158046893aa73d478e23415e18997ece
+  md5: c8217f5bdd5b018087bdea081912cda5
+  depends:
+  - libgcc >=14
+  - libstdcxx >=14
+  - __glibc >=2.17,<3.0.a0
+  - icu >=78.3,<79.0a0
+  license: MIT
+  license_family: MIT
+  size: 72505
+  timestamp: 1786443494718
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libpsl-0.23.1-h7a62e17_0.conda
   sha256: 0140cb0d059ac531e476b85a543668b69b7dffec16f0d26269ab6c9b70546921
   md5: 75a4cdf128d016141db61732ea462276
@@ -2479,17 +2518,17 @@ packages:
   license_family: BSD
   size: 430365
   timestamp: 1780057267477
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libsqlite-3.51.0-hee844dc_0.conda
-  sha256: 4c992dcd0e34b68f843e75406f7f303b1b97c248d18f3c7c330bdc0bc26ae0b3
-  md5: 729a572a3ebb8c43933b30edcc628ceb
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libsqlite-3.53.4-hf4e2dac_0.conda
+  sha256: 72023efc207fe681e26b65fc9d668062cf0b4f0eacf3431e6eb099b95c1f2efd
+  md5: df088a279cd5e6fd2790b4c196434da1
   depends:
   - __glibc >=2.17,<3.0.a0
-  - icu >=75.1,<76.0a0
+  - icu >=78.3,<79.0a0
   - libgcc >=14
-  - libzlib >=1.3.1,<2.0a0
+  - libzlib >=1.3.2,<2.0a0
   license: blessing
-  size: 945576
-  timestamp: 1762299687230
+  size: 964200
+  timestamp: 1785016112246
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libsqlite-3.53.4-h1ae2325_0.conda
   sha256: 745662565e103f290e9dc4263bbd88285082f8cf699854fe2d5f1e35a4a0d326
   md5: 0e3477c0c3e718dcf2eb74ccc8f68570
@@ -2534,15 +2573,6 @@ packages:
   license_family: GPL
   size: 3898269
   timestamp: 1759968103436
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libstdcxx-ng-15.2.0-h4852527_7.conda
-  sha256: 024fd46ac3ea8032a5ec3ea7b91c4c235701a8bf0e6520fe5e6539992a6bd05f
-  md5: f627678cf829bd70bccf141a19c3ad3e
-  depends:
-  - libstdcxx 15.2.0 h8f9b012_7
-  license: GPL-3.0-only WITH GCC-exception-3.1
-  license_family: GPL
-  size: 29343
-  timestamp: 1759968157195
 - conda: https://conda.anaconda.org/conda-forge/linux-64/libuuid-2.41.2-he9a06e4_0.conda
   sha256: e5ec6d2ad7eef538ddcb9ea62ad4346fde70a4736342c4ad87bd713641eb9808
   md5: 80c07c68d2f6870250959dcc95b209d1
@@ -2561,21 +2591,21 @@ packages:
   license: LGPL-2.1-or-later
   size: 100393
   timestamp: 1702724383534
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-2.15.1-h26afc86_0.conda
-  sha256: ec0735ae56c3549149eebd7dc22c0bed91fd50c02eaa77ff418613ddda190aa8
-  md5: e512be7dc1f84966d50959e900ca121f
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-2.15.3-h49c6c72_0.conda
+  sha256: 3bc5551720c58591f6ea1146f7d1539c734ed1c40e7b9f5cb8cb7e900c509aba
+  md5: 995d8c8bad2a3cc8db14675a153dec2b
   depends:
   - __glibc >=2.17,<3.0.a0
-  - icu >=75.1,<76.0a0
+  - icu >=78.3,<79.0a0
   - libgcc >=14
   - libiconv >=1.18,<2.0a0
-  - liblzma >=5.8.1,<6.0a0
-  - libxml2-16 2.15.1 ha9997c6_0
-  - libzlib >=1.3.1,<2.0a0
+  - liblzma >=5.8.3,<6.0a0
+  - libxml2-16 2.15.3 hca6bf5a_0
+  - libzlib >=1.3.2,<2.0a0
   license: MIT
   license_family: MIT
-  size: 45283
-  timestamp: 1761015644057
+  size: 46810
+  timestamp: 1776376751152
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libxml2-2.15.3-h5654f7c_0.conda
   sha256: 2fe1d8de0854342ae9cabe408b476935f82f5636e153b3b497456264dc8ff3a1
   md5: 8e037d73747d6fe34e12d7bcac10cf21
@@ -2590,22 +2620,22 @@ packages:
   license_family: MIT
   size: 41102
   timestamp: 1776377119495
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-16-2.15.1-ha9997c6_0.conda
-  sha256: 71436e72a286ef8b57d6f4287626ff91991eb03c7bdbe835280521791efd1434
-  md5: e7733bc6785ec009e47a224a71917e84
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libxml2-16-2.15.3-hca6bf5a_0.conda
+  sha256: 3d44f737c5ae52d5af32682cc1530df433f401f8e58a7533926536244127572a
+  md5: e79d2c2f24b027aa8d5ab1b1ba3061e7
   depends:
   - __glibc >=2.17,<3.0.a0
-  - icu >=75.1,<76.0a0
+  - icu >=78.3,<79.0a0
   - libgcc >=14
   - libiconv >=1.18,<2.0a0
-  - liblzma >=5.8.1,<6.0a0
-  - libzlib >=1.3.1,<2.0a0
+  - liblzma >=5.8.3,<6.0a0
+  - libzlib >=1.3.2,<2.0a0
   constrains:
-  - libxml2 2.15.1
+  - libxml2 2.15.3
   license: MIT
   license_family: MIT
-  size: 556302
-  timestamp: 1761015637262
+  size: 559775
+  timestamp: 1776376739004
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libxml2-16-2.15.3-h5ef1a60_0.conda
   sha256: ff75b84cdb9e8d123db2fa694a8ac2c2059516b6cbc98ac21fb68e235d0fd354
   md5: 19edaa53885fc8205614b03da2482282
@@ -2621,18 +2651,17 @@ packages:
   license_family: MIT
   size: 466360
   timestamp: 1776377102261
-- conda: https://conda.anaconda.org/conda-forge/linux-64/libzlib-1.3.1-hb9d3cd8_2.conda
-  sha256: d4bfe88d7cb447768e31650f06257995601f89076080e76df55e3112d4e47dc4
-  md5: edb0dca6bc32e4f4789199455a1dbeb8
+- conda: https://conda.anaconda.org/conda-forge/linux-64/libzlib-1.3.2-h25fd6f3_3.conda
+  sha256: eb8a0db0aa570124f7d2a93d7c7f596e3390df5e047818d873baad32985fc736
+  md5: 0de0122d9570a8ab637c6b73db268389
   depends:
   - __glibc >=2.17,<3.0.a0
-  - libgcc >=13
   constrains:
-  - zlib 1.3.1 *_2
+  - zlib 1.3.2 *_3
   license: Zlib
   license_family: Other
-  size: 60963
-  timestamp: 1727963148474
+  size: 63713
+  timestamp: 1785362952714
 - conda: https://conda.anaconda.org/conda-forge/osx-arm64/libzlib-1.3.2-h8088a28_3.conda
   sha256: a18fa5d5bac452401459f966cf0d872224e8080c4ff93c77e168d43ab42ef9d7
   md5: f39288f0ea63ae962e1a2e4f355a0d75
@@ -3400,6 +3429,38 @@ packages:
   license_family: MIT
   size: 15528
   timestamp: 1733710122949
+- conda: https://conda.anaconda.org/bioconda/linux-64/pysam-0.24.0-py312hf5ad864_1.conda
+  sha256: cea331ad0156a060b4ee9f534daee29da451439679e5b4d2d1558ddaa4102f61
+  md5: f94309496d6273c25f0f0a499b75b2fb
+  depends:
+  - __glibc >=2.17,<3.0.a0
+  - bzip2 >=1.0.8,<2.0a0
+  - libcurl >=8.20.0,<9.0a0
+  - libdeflate >=1.25,<1.26.0a0
+  - libgcc >=14
+  - liblzma >=5.8.3,<6.0a0
+  - libzlib >=1.3.2,<2.0a0
+  - openssl >=3.5.7,<4.0a0
+  - python >=3.12,<3.13.0a0
+  - python_abi 3.12.* *_cp312
+  license: MIT
+  size: 3908681
+  timestamp: 1781873141976
+- conda: https://conda.anaconda.org/bioconda/osx-arm64/pysam-0.24.0-py312h12d0683_1.conda
+  sha256: bab067cb5ad9c4603783533b55c0ae8c999e2ae23183a4d150c3c71c9036840f
+  md5: 6959148824e8c370a1578914aa2f7d0c
+  depends:
+  - __osx >=11.0
+  - bzip2 >=1.0.8,<2.0a0
+  - libcurl >=8.20.0,<9.0a0
+  - libdeflate >=1.25,<1.26.0a0
+  - liblzma >=5.8.3,<6.0a0
+  - libzlib >=1.3.2,<2.0a0
+  - python >=3.12,<3.13.0a0
+  - python_abi 3.12.* *_cp312
+  license: MIT
+  size: 3660664
+  timestamp: 1781872373639
 - conda: https://conda.anaconda.org/conda-forge/noarch/pysocks-1.7.1-pyha55dd90_7.conda
   sha256: ba3b032fa52709ce0d9fd388f63d330a026754587a2f461117cac9ab73d8d0d8
   md5: 461219d1a5bd61342293efa2c0c90eac
diff --git a/pixi.toml b/pixi.toml
index ff7499ff9..9e1d148a5 100644
--- a/pixi.toml
+++ b/pixi.toml
@@ -8,9 +8,16 @@ version = "0.2.0"
 
 [tasks]
 fmt = "ruff format . && taplo format pixi.toml && snakefmt workflow/"
-test-data = { cmd = '''bash -c 'if [ -f "$INIT_CWD/fire-test-data/test.cram" ]; then echo "test data already present, skipping download"; else mkdir -p "$INIT_CWD/fire-test-data" && rclone sync ":s3,env_auth=false,provider=Other,endpoint=s3.kopah.orci.washington.edu:stergachis/public/FIRE/test-data" "$INIT_CWD/fire-test-data/"; fi' ''' }
-test-clean = { cmd = '''bash -c 'cd "$INIT_CWD/fire-test-data" && rm -rf results temp .snakemake' ''' }
-test = { cmd = '''bash -c 'cd "$INIT_CWD/fire-test-data" && rm -rf results temp && snakemake -s "$PIXI_PROJECT_ROOT/workflow/Snakefile" --configfile test.yaml -k' ''', depends-on = [
+test-data = { cmd = '''bash -c 'if [ -f "$PIXI_PROJECT_ROOT/fire-test-data/test.cram" ]; then echo "test data already present, skipping download"; else mkdir -p "$PIXI_PROJECT_ROOT/fire-test-data" && rclone sync ":s3,env_auth=false,provider=Other,endpoint=s3.kopah.orci.washington.edu:stergachis/public/FIRE/test-data" "$PIXI_PROJECT_ROOT/fire-test-data/"; fi' ''' }
+test-clean = { cmd = '''bash -c 'cd "$PIXI_PROJECT_ROOT/fire-test-data" && rm -rf results temp .snakemake' ''' }
+test = { cmd = '''bash -c 'cd "$PIXI_PROJECT_ROOT/fire-test-data" && rm -rf results temp && snakemake -s "$PIXI_PROJECT_ROOT/workflow/Snakefile" --configfile test.yaml -k' ''', depends-on = [
+  "test-data",
+], clean-env = true }
+test-dry = { cmd = '''bash -c 'cd "$PIXI_PROJECT_ROOT/fire-test-data" && bash "$PIXI_PROJECT_ROOT/tests/dry-run-checks.sh"' ''', depends-on = [
+  "test-data",
+], clean-env = true }
+test-verify = { cmd = '''bash -c 'cd "$PIXI_PROJECT_ROOT/fire-test-data" && bash "$PIXI_PROJECT_ROOT/tests/verify-results.sh"' ''' }
+test-multi = { cmd = '''bash -c 'cd "$PIXI_PROJECT_ROOT/fire-test-data" && rm -rf results temp && python "$PIXI_PROJECT_ROOT/tests/generate-test-bams.py" && snakemake -s "$PIXI_PROJECT_ROOT/workflow/Snakefile" --configfile "$PIXI_PROJECT_ROOT/tests/config/multi.yaml" -k && bash "$PIXI_PROJECT_ROOT/tests/verify-multi.sh"' ''', depends-on = [
   "test-data",
 ], clean-env = true }
 fire = { cmd = [
@@ -34,6 +41,7 @@ slurm = { cmd = [
 
 [dependencies]
 conda = "*"
+pysam = "*"
 snakemake = "==8.21"
 snakemake-executor-plugin-slurm = ">=0.11.2"
 snakefmt = "*"
diff --git a/tests/README.md b/tests/README.md
new file mode 100644
index 000000000..36e475cf6
--- /dev/null
+++ b/tests/README.md
@@ -0,0 +1,44 @@
+# FIRE test fixtures and scripts
+
+Everything in this directory is test material. Nothing here is an example
+configuration for real use — see `config/` for that. All configs reference
+the downloaded `fire-test-data/` directory and run with that directory as
+the working directory.
+
+## Entry points (pixi tasks)
+
+| Task | What it does |
+|---|---|
+| `pixi run test-dry` | Runs `dry-run-checks.sh`: the full validation matrix below as fast `snakemake -n` checks. Runs first in CI. |
+| `pixi run test` | Full single-sample regression run (unchanged, uses `fire-test-data/test.yaml`). |
+| `pixi run test-verify` | Runs `verify-results.sh`: golden output counts for the regression run. Runs in CI after `test`. |
+| `pixi run test-multi` | Local only. Generates a chr20-only bam and a reversed-header bam (`generate-test-bams.py`), runs the four-sample `config/multi.yaml`, then checks it with `verify-multi.sh`. |
+
+## Fixtures in `config/`
+
+Positive cases (the dry run must pass):
+
+| Fixture | Exercises |
+|---|---|
+| `single-4col.yaml` + `.tbl` | Four-column manifest, no config-level reference. |
+| `two-sample.yaml` + `.tbl` | Two samples, one with `ref_name: hg38` (auto-excludes branch) and one without (empty-excludes branch). Sets `excludes:` to guard the list-copy bug in `get_excludes`. |
+| `sentinel.yaml` + `.tbl` | A `.` cell in `ref`/`ref_name` falls back to the config value. |
+| `override-info.yaml` | Manifest columns and config values both set; asserts the override INFO line. |
+| `multi.yaml` + `.tbl` | Four samples for `test-multi`: `generated/test-chr20.bam` (bam contigs are a subset of the fasta), `generated/test-rev.bam` (header order chr21,chr20 — checks that outputs keep bam header order), and a sentinel row. |
+| `no-excludes.yaml` | Two samples with no `excludes:` key; the rendered shell must show the truly-empty excludes branch. |
+| `extra-exclude.bed` | Small exclude file used by `two-sample.yaml` and `multi.yaml`. |
+
+Error cases (`err-*`, the dry run must fail with a specific message):
+
+| Fixture | Validator branch |
+|---|---|
+| `err-no-ref-anywhere.yaml` | No reference in the manifest or the config. |
+| `err-ref-col-only.yaml` + `.tbl` | Manifest has `ref` without `ref_name`. |
+| `err-config-ref-only.yaml` | Config has `ref` without `ref_name`. |
+| `err-nan-cell.yaml` + `.tbl` | Short manifest row (missing cells). |
+| `err-sentinel-no-config.yaml` + `.tbl` | `.` cell with no config value to fall back to. |
+| `err-dup-sample.yaml` + `.tbl` | Duplicate sample names. |
+| `err-missing-refpath.yaml` + `.tbl` | Reference fasta does not exist. |
+| `err-missing-fai.yaml` + `.tbl` + `no-fai.fa` | Fasta exists but has no `.fai`. |
+| `err-missing-bam.yaml` + `.tbl` | Input bam does not exist. |
+| `err-keepchrs.yaml` | `keep_chromosomes` filters out every chromosome. |
diff --git a/tests/config/err-config-ref-only.yaml b/tests/config/err-config-ref-only.yaml
new file mode 100644
index 000000000..9e84fb19f
--- /dev/null
+++ b/tests/config/err-config-ref-only.yaml
@@ -0,0 +1,2 @@
+manifest: test.tbl
+ref: test.fa.gz
diff --git a/tests/config/err-dup-sample.tbl b/tests/config/err-dup-sample.tbl
new file mode 100644
index 000000000..50cc1fa84
--- /dev/null
+++ b/tests/config/err-dup-sample.tbl
@@ -0,0 +1,3 @@
+sample	bam
+test	test.cram
+test	test.cram
diff --git a/tests/config/err-dup-sample.yaml b/tests/config/err-dup-sample.yaml
new file mode 100644
index 000000000..487ec70bd
--- /dev/null
+++ b/tests/config/err-dup-sample.yaml
@@ -0,0 +1,3 @@
+manifest: ../tests/config/err-dup-sample.tbl
+ref: test.fa.gz
+ref_name: hg38
diff --git a/tests/config/err-keepchrs.yaml b/tests/config/err-keepchrs.yaml
new file mode 100644
index 000000000..97f254ab3
--- /dev/null
+++ b/tests/config/err-keepchrs.yaml
@@ -0,0 +1,4 @@
+ref: test.fa.gz
+ref_name: hg38
+manifest: test.tbl
+keep_chromosomes: "chrX"
diff --git a/tests/config/err-missing-bam.tbl b/tests/config/err-missing-bam.tbl
new file mode 100644
index 000000000..ef0be9158
--- /dev/null
+++ b/tests/config/err-missing-bam.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref	ref_name
+test	nope.cram	test.fa.gz	hg38
diff --git a/tests/config/err-missing-bam.yaml b/tests/config/err-missing-bam.yaml
new file mode 100644
index 000000000..912850091
--- /dev/null
+++ b/tests/config/err-missing-bam.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/err-missing-bam.tbl
diff --git a/tests/config/err-missing-fai.tbl b/tests/config/err-missing-fai.tbl
new file mode 100644
index 000000000..08caed378
--- /dev/null
+++ b/tests/config/err-missing-fai.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref	ref_name
+test	test.cram	../tests/config/no-fai.fa	hg38
diff --git a/tests/config/err-missing-fai.yaml b/tests/config/err-missing-fai.yaml
new file mode 100644
index 000000000..54d2f4923
--- /dev/null
+++ b/tests/config/err-missing-fai.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/err-missing-fai.tbl
diff --git a/tests/config/err-missing-refpath.tbl b/tests/config/err-missing-refpath.tbl
new file mode 100644
index 000000000..83d42911e
--- /dev/null
+++ b/tests/config/err-missing-refpath.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref	ref_name
+test	test.cram	nope.fa	hg38
diff --git a/tests/config/err-missing-refpath.yaml b/tests/config/err-missing-refpath.yaml
new file mode 100644
index 000000000..4bbcfa376
--- /dev/null
+++ b/tests/config/err-missing-refpath.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/err-missing-refpath.tbl
diff --git a/tests/config/err-nan-cell.tbl b/tests/config/err-nan-cell.tbl
new file mode 100644
index 000000000..7f2924577
--- /dev/null
+++ b/tests/config/err-nan-cell.tbl
@@ -0,0 +1,3 @@
+sample	bam	ref	ref_name
+test	test.cram	test.fa.gz	hg38
+test2	test.cram	test.fa.gz
diff --git a/tests/config/err-nan-cell.yaml b/tests/config/err-nan-cell.yaml
new file mode 100644
index 000000000..906d3274a
--- /dev/null
+++ b/tests/config/err-nan-cell.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/err-nan-cell.tbl
diff --git a/tests/config/err-no-ref-anywhere.yaml b/tests/config/err-no-ref-anywhere.yaml
new file mode 100644
index 000000000..0cb6d126b
--- /dev/null
+++ b/tests/config/err-no-ref-anywhere.yaml
@@ -0,0 +1 @@
+manifest: test.tbl
diff --git a/tests/config/err-ref-col-only.tbl b/tests/config/err-ref-col-only.tbl
new file mode 100644
index 000000000..dec20c0a8
--- /dev/null
+++ b/tests/config/err-ref-col-only.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref
+test	test.cram	test.fa.gz
diff --git a/tests/config/err-ref-col-only.yaml b/tests/config/err-ref-col-only.yaml
new file mode 100644
index 000000000..50f5a8801
--- /dev/null
+++ b/tests/config/err-ref-col-only.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/err-ref-col-only.tbl
diff --git a/tests/config/err-sentinel-no-config.tbl b/tests/config/err-sentinel-no-config.tbl
new file mode 100644
index 000000000..b5a0612b2
--- /dev/null
+++ b/tests/config/err-sentinel-no-config.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref	ref_name
+test	test.cram	.	.
diff --git a/tests/config/err-sentinel-no-config.yaml b/tests/config/err-sentinel-no-config.yaml
new file mode 100644
index 000000000..d8cb5ece0
--- /dev/null
+++ b/tests/config/err-sentinel-no-config.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/err-sentinel-no-config.tbl
diff --git a/tests/config/extra-exclude.bed b/tests/config/extra-exclude.bed
new file mode 100644
index 000000000..44655c5ff
--- /dev/null
+++ b/tests/config/extra-exclude.bed
@@ -0,0 +1 @@
+chr20	0	10000
diff --git a/tests/config/multi.tbl b/tests/config/multi.tbl
new file mode 100644
index 000000000..cd6f818fa
--- /dev/null
+++ b/tests/config/multi.tbl
@@ -0,0 +1,5 @@
+sample	bam	ref	ref_name
+test	test.cram	test.fa.gz	hg38
+test2	test.cram	test.fa.gz	testref
+test20	generated/test-chr20.bam	test.fa.gz	.
+test-rev	generated/test-rev.bam	test.fa.gz	hg38
diff --git a/tests/config/multi.yaml b/tests/config/multi.yaml
new file mode 100644
index 000000000..24ff66189
--- /dev/null
+++ b/tests/config/multi.yaml
@@ -0,0 +1,5 @@
+manifest: ../tests/config/multi.tbl
+ref: test.fa.gz
+ref_name: testref2
+excludes:
+  - ../tests/config/extra-exclude.bed
diff --git a/tests/config/no-excludes.yaml b/tests/config/no-excludes.yaml
new file mode 100644
index 000000000..cc1265e71
--- /dev/null
+++ b/tests/config/no-excludes.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/two-sample.tbl
diff --git a/tests/config/no-fai.fa b/tests/config/no-fai.fa
new file mode 100644
index 000000000..f5fc23851
--- /dev/null
+++ b/tests/config/no-fai.fa
@@ -0,0 +1,2 @@
+>tiny
+ACGT
diff --git a/tests/config/override-info.yaml b/tests/config/override-info.yaml
new file mode 100644
index 000000000..e77abe1d5
--- /dev/null
+++ b/tests/config/override-info.yaml
@@ -0,0 +1,3 @@
+manifest: ../tests/config/single-4col.tbl
+ref: test.fa.gz
+ref_name: hg38
diff --git a/tests/config/sentinel.tbl b/tests/config/sentinel.tbl
new file mode 100644
index 000000000..b5a0612b2
--- /dev/null
+++ b/tests/config/sentinel.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref	ref_name
+test	test.cram	.	.
diff --git a/tests/config/sentinel.yaml b/tests/config/sentinel.yaml
new file mode 100644
index 000000000..13f516a20
--- /dev/null
+++ b/tests/config/sentinel.yaml
@@ -0,0 +1,3 @@
+manifest: ../tests/config/sentinel.tbl
+ref: test.fa.gz
+ref_name: hg38
diff --git a/tests/config/single-4col.tbl b/tests/config/single-4col.tbl
new file mode 100644
index 000000000..306a70b75
--- /dev/null
+++ b/tests/config/single-4col.tbl
@@ -0,0 +1,2 @@
+sample	bam	ref	ref_name
+test	test.cram	test.fa.gz	hg38
diff --git a/tests/config/single-4col.yaml b/tests/config/single-4col.yaml
new file mode 100644
index 000000000..0c7cde333
--- /dev/null
+++ b/tests/config/single-4col.yaml
@@ -0,0 +1 @@
+manifest: ../tests/config/single-4col.tbl
diff --git a/tests/config/two-sample.tbl b/tests/config/two-sample.tbl
new file mode 100644
index 000000000..040266d51
--- /dev/null
+++ b/tests/config/two-sample.tbl
@@ -0,0 +1,3 @@
+sample	bam	ref	ref_name
+test	test.cram	test.fa.gz	hg38
+test2	test.cram	test.fa.gz	testref
diff --git a/tests/config/two-sample.yaml b/tests/config/two-sample.yaml
new file mode 100644
index 000000000..31cc2a867
--- /dev/null
+++ b/tests/config/two-sample.yaml
@@ -0,0 +1,3 @@
+manifest: ../tests/config/two-sample.tbl
+excludes:
+  - ../tests/config/extra-exclude.bed
diff --git a/tests/dry-run-checks.sh b/tests/dry-run-checks.sh
new file mode 100644
index 000000000..d97ad04cc
--- /dev/null
+++ b/tests/dry-run-checks.sh
@@ -0,0 +1,153 @@
+#!/usr/bin/env bash
+# Dry-run validation matrix for the FIRE manifest and reference handling.
+# Run with the working directory set to fire-test-data (see the test-dry pixi task).
+set -uo pipefail
+
+SNAKEFILE="$PIXI_PROJECT_ROOT/workflow/Snakefile"
+CFG="$PIXI_PROJECT_ROOT/tests/config"
+FAILURES=0
+
+snk() {
+    snakemake -s "$SNAKEFILE" -n -q rules --configfile "$@"
+}
+
+expect_pass() {
+    local config=$1
+    if ! out=$(snk "$config" 2>&1); then
+        echo "FAIL (expected pass): $config"
+        echo "$out" | tail -5
+        FAILURES=$((FAILURES + 1))
+    else
+        echo "ok (pass): $config"
+    fi
+}
+
+expect_fail() {
+    local config=$1 substring=$2
+    if out=$(snk "$config" 2>&1); then
+        echo "FAIL (expected failure): $config"
+        FAILURES=$((FAILURES + 1))
+    elif ! grep -qF "$substring" <<<"$out"; then
+        echo "FAIL (wrong message): $config, wanted: $substring"
+        echo "$out" | tail -5
+        FAILURES=$((FAILURES + 1))
+    else
+        echo "ok (fail): $config"
+    fi
+}
+
+expect_stderr() {
+    local config=$1 substring=$2
+    if ! out=$(snk "$config" 2>&1); then
+        echo "FAIL (expected pass): $config"
+        FAILURES=$((FAILURES + 1))
+    elif ! grep -qF "$substring" <<<"$out"; then
+        echo "FAIL (missing stderr line): $config, wanted: $substring"
+        FAILURES=$((FAILURES + 1))
+    else
+        echo "ok (stderr): $config"
+    fi
+}
+
+# positive cases
+expect_pass test.yaml
+expect_pass "$CFG/single-4col.yaml"
+expect_pass "$CFG/two-sample.yaml"
+expect_pass "$CFG/sentinel.yaml"
+expect_stderr "$CFG/override-info.yaml" "manifest ref/ref_name columns override config-level values"
+
+# sample-attributed rendered-shell checks for exclude_from_shuffle;
+# -R forces the rule so cached results do not hide the rendering. The
+# rendered shell contains the sample-scoped output path, so paragraphs
+# (awk RS='') can be attributed to one sample.
+shell_block() {
+    # config, sample -> the rendered exclude_from_shuffle shell for sample
+    snakemake -s "$SNAKEFILE" -n -p -R exclude_from_shuffle --configfile "$1" 2>&1 \
+        | awk -v RS='' -v sm="results/$2/" '/bedtools genomecov/ && index($0, sm)'
+}
+
+check_block() {
+    local label=$1 block=$2 must=$3 must_not=$4
+    if [ -z "$block" ]; then
+        echo "FAIL: no rendered exclude_from_shuffle shell for $label"
+        FAILURES=$((FAILURES + 1))
+        return
+    fi
+    if [ -n "$must" ] && ! grep -qF "$must" <<<"$block"; then
+        echo "FAIL: $label rendered shell is missing: $must"
+        FAILURES=$((FAILURES + 1))
+        return
+    fi
+    if [ -n "$must_not" ] && grep -qF "$must_not" <<<"$block"; then
+        echo "FAIL: $label rendered shell must not contain: $must_not"
+        FAILURES=$((FAILURES + 1))
+        return
+    fi
+    echo "ok (rendered shell): $label"
+}
+
+# with config excludes set: test gets hg38 blacklists + the config exclude,
+# test2 gets only the config exclude (a leak of hg38 paths into test2 fails)
+block_test=$(shell_block "$CFG/two-sample.yaml" test)
+block_test2=$(shell_block "$CFG/two-sample.yaml" test2)
+check_block "test (hg38 excludes)" "$block_test" "hg38.blacklist" ""
+check_block "test (config exclude)" "$block_test" "extra-exclude.bed" ""
+check_block "test2 (no hg38 leak)" "$block_test2" "extra-exclude.bed" "hg38.blacklist"
+# the shell template references {params.exclude} twice (the [ -n ] guard
+# and gunzip), so each exclude path renders exactly twice; more means a
+# get_excludes mutation bug duplicated the list
+gap_count=$(grep -oF "hg38.gap.bed.gz" <<<"$block_test" | wc -l | tr -d ' ')
+if [ "$gap_count" -ne 2 ]; then
+    echo "FAIL: hg38.gap.bed.gz appears $gap_count times for test, expected 2"
+    FAILURES=$((FAILURES + 1))
+else
+    echo "ok (rendered shell): no exclude duplication for test"
+fi
+
+# without config excludes: test2 exercises the truly-empty excludes branch
+# (the [ -n ] guard must render with an empty parameter)
+block_empty=$(shell_block "$CFG/no-excludes.yaml" test2)
+check_block "test2 (empty excludes)" "$block_empty" "" "gunzip -cf ."
+if grep -qF 'if [ -n "" ]' <<<"$block_empty"; then
+    echo "ok (rendered shell): empty excludes guard renders for test2"
+else
+    echo "FAIL: empty excludes guard did not render empty for test2"
+    FAILURES=$((FAILURES + 1))
+fi
+
+# a leftover genome file that no longer matches the bam header must be
+# removed at parse time (the staleness guard for the ancient() bam input)
+mkdir -p temp/test
+printf 'chrBOGUS\t12345\n' >temp/test/test.genome
+out=$(snk test.yaml 2>&1) || {
+    echo "FAIL: stale-genome dry-run did not build"
+    FAILURES=$((FAILURES + 1))
+}
+if ! grep -qF "removing stale genome file" <<<"$out"; then
+    echo "FAIL: stale genome file was not reported"
+    FAILURES=$((FAILURES + 1))
+elif [ -f temp/test/test.genome ]; then
+    echo "FAIL: stale genome file was not removed"
+    FAILURES=$((FAILURES + 1))
+else
+    echo "ok (stale guard): mismatched genome file removed at parse time"
+fi
+
+# error cases, one per validator branch; substrings are specific enough
+# that a wrong error cannot satisfy the assertion
+expect_fail "$CFG/err-no-ref-anywhere.yaml" "no reference specified"
+expect_fail "$CFG/err-ref-col-only.yaml" "manifest columns 'ref' and 'ref_name' must be provided together"
+expect_fail "$CFG/err-config-ref-only.yaml" "config options 'ref' and 'ref_name' must be provided together"
+expect_fail "$CFG/err-nan-cell.yaml" "missing or malformed manifest fields"
+expect_fail "$CFG/err-sentinel-no-config.yaml" "is not set in config.yaml"
+expect_fail "$CFG/err-dup-sample.yaml" "duplicate sample names"
+expect_fail "$CFG/err-missing-refpath.yaml" "reference file"
+expect_fail "$CFG/err-missing-fai.yaml" "reference index file"
+expect_fail "$CFG/err-missing-bam.yaml" "cannot read input bam"
+expect_fail "$CFG/err-keepchrs.yaml" "no chromosomes left"
+
+if [ "$FAILURES" -gt 0 ]; then
+    echo "dry-run-checks: $FAILURES failure(s)"
+    exit 1
+fi
+echo "dry-run-checks: all checks passed"
diff --git a/tests/generate-test-bams.py b/tests/generate-test-bams.py
new file mode 100644
index 000000000..d734c2f67
--- /dev/null
+++ b/tests/generate-test-bams.py
@@ -0,0 +1,69 @@
+"""Make the generated test BAMs from the test CRAM.
+
+Two BAMs, each skipped when already present:
+
+- generated/test-chr20.bam: header lists only chr20, so the sample
+  exercises the BAM-contigs-are-a-subset-of-the-fasta case.
+- generated/test-rev.bam: header lists chr21 before chr20, so the sample
+  exercises header-order preservation (header order differs from
+  lexicographic order).
+
+Run with the working directory set to fire-test-data (see the test-multi
+pixi task).
+"""
+
+from pathlib import Path
+
+import pysam
+
+CRAM = "test.cram"
+REF = "test.fa.gz"
+OUT_DIR = Path("generated")
+
+
+def remap(read: pysam.AlignedSegment, out: pysam.AlignmentFile) -> pysam.AlignedSegment:
+    """Point the read's tids at `out`'s header by contig name."""
+    if read.reference_name is not None:
+        read.reference_id = out.get_tid(read.reference_name)
+    if read.next_reference_id >= 0 and read.next_reference_name is not None:
+        read.next_reference_id = out.get_tid(read.next_reference_name)
+    return read
+
+
+def make_chr20() -> None:
+    out = OUT_DIR / "test-chr20.bam"
+    if out.exists():
+        print(f"{out} already present, skipping generation")
+        return
+    with pysam.AlignmentFile(CRAM, "rc", reference_filename=REF) as cram:
+        header = cram.header.to_dict()
+        header["SQ"] = [sq for sq in header["SQ"] if sq["SN"] == "chr20"]
+        with pysam.AlignmentFile(out, "wb", header=header) as bam:
+            for read in cram.fetch("chr20"):
+                bam.write(remap(read, bam))
+    pysam.index(str(out))
+    print(f"wrote {out}")
+
+
+def make_reversed() -> None:
+    out = OUT_DIR / "test-rev.bam"
+    if out.exists():
+        print(f"{out} already present, skipping generation")
+        return
+    unsorted = OUT_DIR / "test-rev.unsorted.bam"
+    with pysam.AlignmentFile(CRAM, "rc", reference_filename=REF) as cram:
+        header = cram.header.to_dict()
+        header["SQ"] = list(reversed(header["SQ"]))
+        with pysam.AlignmentFile(unsorted, "wb", header=header) as bam:
+            for read in cram.fetch(until_eof=True):
+                bam.write(remap(read, bam))
+    pysam.sort("-o", str(out), str(unsorted))
+    unsorted.unlink()
+    pysam.index(str(out))
+    print(f"wrote {out}")
+
+
+if __name__ == "__main__":
+    OUT_DIR.mkdir(exist_ok=True)
+    make_chr20()
+    make_reversed()
diff --git a/tests/verify-multi.sh b/tests/verify-multi.sh
new file mode 100644
index 000000000..e4222e707
--- /dev/null
+++ b/tests/verify-multi.sh
@@ -0,0 +1,52 @@
+#!/usr/bin/env bash
+# Checks for the multi-sample run (four samples, per-sample references).
+# Run with the working directory set to fire-test-data after `pixi run test-multi`.
+set -euo pipefail
+
+V="v$(echo "$PIXI_PROJECT_VERSION" | cut -d. -f1-2)"
+FAILURES=0
+
+assert_genome() {
+    local sm=$1 name=$2
+    local genomes="results/$sm/trackHub-$V/genomes.txt"
+    if ! grep -qx "genome $name" "$genomes"; then
+        echo "FAIL: $genomes does not contain 'genome $name'"
+        FAILURES=$((FAILURES + 1))
+    else
+        echo "ok: $sm uses genome $name"
+    fi
+}
+
+assert_genome test hg38
+assert_genome test2 testref
+assert_genome test20 testref2
+assert_genome test-rev hg38
+
+# the chr20-only sample must have no chr21 peaks; grep -c reads the whole
+# stream, so no SIGPIPE can mask a hit under pipefail
+peaks20="results/test20/test20-fire-$V-peaks.bed.gz"
+chr21_count=$(gunzip -c "$peaks20" | cut -f 1 | { grep -cx "chr21" || true; })
+if [ "$chr21_count" -gt 0 ]; then
+    echo "FAIL: $peaks20 contains $chr21_count chr21 records"
+    FAILURES=$((FAILURES + 1))
+else
+    echo "ok: test20 output is chr20 only"
+fi
+
+# the reversed-header sample must keep bam header order (chr21 first) in
+# concatenated per-chromosome outputs
+pileup_rev="results/test-rev/test-rev-fire-$V-pileup.bed.gz"
+# || true absorbs the SIGPIPE that head sends up the pipeline
+first_chrom=$(gunzip -c "$pileup_rev" | grep -v "^#" | head -n 1 | cut -f 1 || true)
+if [ "$first_chrom" != "chr21" ]; then
+    echo "FAIL: $pileup_rev starts with $first_chrom, expected chr21 (header order)"
+    FAILURES=$((FAILURES + 1))
+else
+    echo "ok: test-rev output preserves bam header order"
+fi
+
+if [ "$FAILURES" -gt 0 ]; then
+    echo "verify-multi: $FAILURES failure(s)"
+    exit 1
+fi
+echo "verify-multi: all checks passed"
diff --git a/tests/verify-results.sh b/tests/verify-results.sh
new file mode 100644
index 000000000..1566c2d29
--- /dev/null
+++ b/tests/verify-results.sh
@@ -0,0 +1,28 @@
+#!/usr/bin/env bash
+# Golden-count checks for the single-sample regression test.
+# Run with the working directory set to fire-test-data after `pixi run test`.
+set -euo pipefail
+
+V="v$(echo "$PIXI_PROJECT_VERSION" | cut -d. -f1-2)"
+PEAKS="results/test/test-fire-$V-peaks.bed.gz"
+ELEMENTS="results/test/additional-outputs-$V/fire-peaks/test-$V-fire-elements.bed.gz"
+GENOMES="results/test/trackHub-$V/genomes.txt"
+
+peaks=$(gunzip -c "$PEAKS" | wc -l | tr -d ' ')
+elements=$(gunzip -c "$ELEMENTS" | wc -l | tr -d ' ')
+echo "peaks file lines: $peaks (expected 1181)"
+echo "fire elements: $elements (expected 46223)"
+
+FAILURES=0
+[ "$peaks" -eq 1181 ] || FAILURES=$((FAILURES + 1))
+[ "$elements" -eq 46223 ] || FAILURES=$((FAILURES + 1))
+if ! grep -qx "genome hg38" "$GENOMES"; then
+    echo "FAIL: $GENOMES does not contain 'genome hg38'"
+    FAILURES=$((FAILURES + 1))
+fi
+
+if [ "$FAILURES" -gt 0 ]; then
+    echo "verify-results: $FAILURES failure(s)"
+    exit 1
+fi
+echo "verify-results: all checks passed"
diff --git a/workflow/Snakefile b/workflow/Snakefile
index 4cb4ca372..1cf3304eb 100644
--- a/workflow/Snakefile
+++ b/workflow/Snakefile
@@ -1,5 +1,6 @@
-# pandas is a dependency of snakemake so a runner environment is not needed
-# if I ever add this back because of new depedancies this must be the first line of the Snakefile
+# pandas and pysam come from the pixi runner environment (pixi.toml [dependencies]);
+# envs/runner.yaml mirrors them for users who run snakemake without pixi.
+# If a conda runner env is ever re-enabled, this must be the first line of the Snakefile:
 # conda: "envs/runner.yaml"
 
 # TODO get the container working
@@ -32,11 +33,7 @@ if config.get("full-version", False):
 MAX_THREADS = config.get("max_threads", 4)
 SORT_THREADS = config.get("sort_threads", 8)
 
-# reference genome and reference regions
-REF = get_ref()
-FAI = get_fai()
-REF_NAME = config["ref_name"]
-EXCLUDES = get_excludes()
+# chromosome filtering; references are per-sample (see the manifest columns)
 KEEP_CHRS = config.get("keep_chromosomes", ".*")
 
 # coverage requirements
@@ -44,9 +41,9 @@ MIN_COVERAGE = config.get("min_coverage", 4)
 COVERAGE_WITHIN_N_SD = config.get("coverage_within_n_sd", 5)
 
 # sample, haplotype, and chromosome wildcard building
-FAI_DF = get_fai_df()
 DEFAULT_ENV = config.get("env", "../envs/env.yaml")
 MANIFEST = get_manifest()
+remove_stale_genome_files()
 MIN_FIRE_FDR = config.get("min_fire_fdr", 0.10)
 
 # FDR / peak calling thresholds
@@ -100,9 +97,9 @@ include: "rules/track-hub.smk"
 
 
 wildcard_constraints:
-    chrom="|".join(get_chroms()),
+    chrom="|".join(re.escape(chrom) for chrom in all_chroms()),
     call="|".join(["msp", "m6a"]),
-    sm="|".join(MANIFEST.index),
+    sm="|".join(re.escape(sm) for sm in MANIFEST.index),
     types="|".join(types),
     fdr=r"\d+",
     hp="|".join(haps),
diff --git a/workflow/envs/runner.yaml b/workflow/envs/runner.yaml
index dddf584ec..b13f49efe 100644
--- a/workflow/envs/runner.yaml
+++ b/workflow/envs/runner.yaml
@@ -5,7 +5,4 @@ channels:
 dependencies:
   - numpy
   - pandas==2.2.3
-#- tqdm
-#- pip
-#- pip:
-#- pysam==0.22.1
+  - pysam
diff --git a/workflow/rules/apply-model.smk b/workflow/rules/apply-model.smk
index af721bfd2..18ec9ad50 100644
--- a/workflow/rules/apply-model.smk
+++ b/workflow/rules/apply-model.smk
@@ -4,7 +4,7 @@
 rule fire:
     input:
         bam=ancient(get_input_bam),
-        ref=ancient(REF),
+        ref=lambda wc: ancient(get_ref(wc)),
     output:
         cram="results/{sm}/{sm}-fire-{v}-filtered.cram",
         crai="results/{sm}/{sm}-fire-{v}-filtered.cram.crai",
@@ -75,7 +75,7 @@ rule fire_sites_chrom:
 rule fire_sites:
     input:
         beds=expand(
-            rules.fire_sites_chrom.output.bed, chrom=get_chroms(), allow_missing=True
+            rules.fire_sites_chrom.output.bed, chrom=get_chroms, allow_missing=True
         ),
     output:
         bed="results/{sm}/additional-outputs-{v}/fire-peaks/{sm}-{v}-fire-elements.bed.gz",
diff --git a/workflow/rules/common.smk b/workflow/rules/common.smk
index 2dec22b79..00211edf4 100644
--- a/workflow/rules/common.smk
+++ b/workflow/rules/common.smk
@@ -1,29 +1,95 @@
 import re
 import logging
 import sys
+from functools import lru_cache
+
+import pysam
+
+# marks a manifest cell that falls back to the config value
+MANIFEST_NA = "."
+
+
+@lru_cache(maxsize=None)
+def _bam_contigs(bam):
+    """Read (name, length) pairs from the header of a BAM/CRAM, once per file."""
+    verbosity = pysam.set_verbosity(0)
+    try:
+        with pysam.AlignmentFile(bam, check_sq=False, require_index=False) as f:
+            contigs = tuple(zip(f.references, f.lengths))
+    except (OSError, ValueError) as e:
+        raise ValueError(f"FIRE: cannot read input bam {bam}: {e}") from e
+    finally:
+        pysam.set_verbosity(verbosity)
+    if not contigs:
+        raise ValueError(
+            f"FIRE: input bam {bam} has no reference sequences (@SQ) in its header; "
+            "FIRE requires an aligned bam"
+        )
+    return contigs
+
 
-FIRST_REPORT = True
+@lru_cache(maxsize=None)
+def _sample_chroms(sm):
+    """Filtered chromosome names for one sample, in BAM header order.
 
+    Header order is the sort order of the data (mosdepth output, coordinate
+    sorted reads), so it must be preserved for the sorted bedtools
+    operations downstream.
+    """
+    min_contig_length = config.get("min_contig_length", 0)
+    try:
+        contigs = _bam_contigs(MANIFEST.loc[sm, "bam"])
+    except ValueError as e:
+        raise ValueError(f"{e} (sample '{sm}')") from e
+    skipped = [name for name, length in contigs if length < min_contig_length]
+    if skipped:
+        print(
+            f"WARNING: {sm}: skipping contigs with length < {min_contig_length:,}: {skipped}",
+            file=sys.stderr,
+        )
+    chroms = tuple(
+        name
+        for name, length in contigs
+        if length >= min_contig_length
+        and "chrUn_" not in name
+        and "_random" not in name
+        and re.fullmatch(KEEP_CHRS, name)
+    )
+    print(f"INFO: {sm}: using {len(chroms)} chromosomes", file=sys.stderr)
+    if not chroms:
+        raise ValueError(
+            f"FIRE: no chromosomes left for sample '{sm}' after filtering. "
+            "Check the keep_chromosomes and min_contig_length options in config.yaml. "
+            f"The bam header contains: {[name for name, _ in contigs]}"
+        )
+    return chroms
 
-def get_ref():
-    if "ref" not in config:
-        raise ValueError("FIRE: ref parameter is missing in config.yaml")
-    ref = config["ref"]
-    if not os.path.exists(ref):
-        raise ValueError(f"FIRE: reference file {ref} does not exist")
-    return os.path.abspath(ref)
 
+def get_ref(wc):
+    return MANIFEST.loc[wc.sm, "ref"]
 
-def get_fai():
-    fai = f"{get_ref()}.fai"
-    if not os.path.exists(fai):
-        raise ValueError(f"FIRE: reference index file {fai} does not exist")
-    return fai
 
+def get_fai(wc):
+    return f"{get_ref(wc)}.fai"
 
-def get_excludes():
-    excludes = config.get("excludes", [])
-    if REF_NAME == "hg38" or REF_NAME == "GRCh38":
+
+def get_ref_name(wc):
+    return MANIFEST.loc[wc.sm, "ref_name"]
+
+
+def get_chroms(wc):
+    return list(_sample_chroms(wc.sm))
+
+
+def all_chroms():
+    # sorted only for a deterministic wildcard-constraint regex; the
+    # alternation order has no effect on matching
+    return sorted({chrom for sm in MANIFEST.index for chrom in _sample_chroms(sm)})
+
+
+def get_excludes(wc):
+    excludes = list(config.get("excludes", []))
+    if get_ref_name(wc) in ["hg38", "GRCh38"]:
         files = [
             "../annotations/hg38.gap.bed.gz",
             "../annotations/hg38.blacklist.ENCFF356LFX.bed.gz",
@@ -33,47 +99,110 @@ def get_excludes():
     return excludes
 
 
-def get_fai_df():
-    fai = get_fai()
-    return pd.read_csv(fai, sep="\t", names=["chr", "length", "x", "y", "z"])
+def _config_ref_value(col):
+    """The config value for ref/ref_name, or None when absent or empty."""
+    value = config.get(col)
+    if value is None or str(value) == "":
+        return None
+    return str(value)
 
 
-def get_chroms():
-    global FIRST_REPORT
-    min_contig_length = config.get("min_contig_length", 0)
-    skipped_contigs = FAI_DF["chr"][FAI_DF["length"] < min_contig_length]
-    if len(skipped_contigs) > 0 and FIRST_REPORT:
+def _fill_manifest_refs(manifest):
+    """Fill and validate the ref and ref_name manifest columns."""
+    empty_in_config = [
+        col
+        for col in ["ref", "ref_name"]
+        if col in config and _config_ref_value(col) is None
+    ]
+    if empty_in_config:
+        raise ValueError(
+            f"FIRE: config options {empty_in_config} are empty in config.yaml; "
+            "set a value or remove the key"
+        )
+    in_manifest = [col for col in ["ref", "ref_name"] if col in manifest.columns]
+    in_config = [
+        col for col in ["ref", "ref_name"] if _config_ref_value(col) is not None
+    ]
+    if len(in_manifest) == 1:
+        raise ValueError(
+            "FIRE: manifest columns 'ref' and 'ref_name' must be provided together "
+            f"(found only '{in_manifest[0]}')"
+        )
+    if len(in_config) == 1:
+        raise ValueError(
+            "FIRE: config options 'ref' and 'ref_name' must be provided together "
+            f"(found only '{in_config[0]}')"
+        )
+    if not in_manifest and not in_config:
+        raise ValueError(
+            "FIRE: no reference specified: add 'ref' and 'ref_name' columns to the "
+            "manifest, or set 'ref' and 'ref_name' in config.yaml"
+        )
+    if in_manifest and in_config:
         print(
-            f"WARNING: Skipping contigs with length < {min_contig_length:,}: {skipped_contigs}",
+            "INFO: manifest ref/ref_name columns override config-level values",
             file=sys.stderr,
         )
+    for col in ["ref", "ref_name"]:
+        if col not in manifest.columns:
+            manifest[col] = _config_ref_value(col)
+            continue
+        sentinel = manifest[col] == MANIFEST_NA
+        if sentinel.any():
+            if _config_ref_value(col) is None:
+                missing = manifest.index[sentinel].tolist()
+                raise ValueError(
+                    f"FIRE: samples {missing} use '{MANIFEST_NA}' for '{col}' in the "
+                    f"manifest, but '{col}' is not set in config.yaml"
+                )
+            manifest.loc[sentinel, col] = _config_ref_value(col)
+    return manifest
 
-    chroms = FAI_DF["chr"][FAI_DF["length"] >= min_contig_length]
-    chroms = sorted([chrom for chrom in chroms if "chrUn_" not in chrom])
-    chroms = [chrom for chrom in chroms if "_random" not in chrom]
-    chroms = [chrom for chrom in chroms if re.fullmatch(KEEP_CHRS, chrom)]
-
-    if FIRST_REPORT:
-        FIRST_REPORT = False
-        print(f"INFO: Using N chromosomes: {len(chroms)}", file=sys.stderr)
 
-    if len(chroms) == 0:
+def get_manifest():
+    manifest_path = config.get("manifest")
+    if manifest_path is None:
+        raise ValueError("FIRE: manifest parameter is missing in config.yaml")
+    if not os.path.exists(manifest_path):
+        raise ValueError(f"FIRE: manifest file {manifest_path} does not exist")
+    try:
+        # dtype=str + keep_default_na=False keep every cell as literal text:
+        # numeric sample names stay strings, a sample named NA stays "NA",
+        # and missing trailing cells parse as "" instead of NaN
+        manifest = pd.read_csv(
+            manifest_path, sep=r"\s+", comment="#", dtype=str, keep_default_na=False
+        )
+    except (pd.errors.ParserError, pd.errors.EmptyDataError) as e:
+        raise ValueError(f"FIRE: cannot parse manifest {manifest_path}: {e}") from e
+    for col in ["sample", "bam"]:
+        if col not in manifest.columns:
+            raise ValueError(
+                f"FIRE: manifest must have 'sample' and 'bam' columns; "
+                f"found: {list(manifest.columns)}"
+            )
+    if len(manifest) == 0:
+        raise ValueError(f"FIRE: manifest {manifest_path} has no samples")
+    dups = manifest["sample"][manifest["sample"].duplicated()].tolist()
+    if dups:
+        raise ValueError(f"FIRE: duplicate sample names in manifest: {dups}")
+    for sm in manifest["sample"]:
+        if not re.fullmatch(r"[A-Za-z0-9_.-]+", sm):
+            raise ValueError(f"FIRE: sample name '{sm}' must match [A-Za-z0-9_.-]+")
+    manifest = manifest.set_index("sample")
+    manifest = _fill_manifest_refs(manifest)
+    malformed = manifest.index[manifest[["bam", "ref", "ref_name"]].eq("").any(axis=1)]
+    if len(malformed) > 0:
         raise ValueError(
-            f"No chromosomes left after filtering. Check your keep_chromosomes parameter in config.yaml. "
-            f"Your fai file contains the following chromosomes: {FAI_DF['chr']}"
+            f"FIRE: samples {malformed.tolist()} have missing or malformed manifest "
+            "fields; every row must fill all manifest columns "
+            f"(use '{MANIFEST_NA}' in ref/ref_name to fall back to the config value)"
         )
-    return chroms
-
-
-def get_manifest():
-    manifest = config.get("manifest")
-    if manifest is None:
-        raise ValueError("manifest parameter is missing in config.yaml")
-    if not os.path.exists(manifest):
-        raise ValueError(f"Manifest file {manifest} does not exist")
-    manifest = pd.read_csv(config["manifest"], sep=r"\s+", comment="#").set_index(
-        "sample"
-    )
+    manifest["ref"] = manifest["ref"].map(os.path.abspath)
+    for ref in manifest["ref"].unique():
+        if not os.path.isfile(ref):
+            raise ValueError(f"FIRE: reference file {ref} does not exist")
+        if not os.path.isfile(f"{ref}.fai"):
+            raise ValueError(f"FIRE: reference index file {ref}.fai does not exist")
     return manifest
 
 
@@ -81,6 +210,34 @@ def get_input_bam(wc):
     return MANIFEST.loc[wc.sm, "bam"]
 
 
+def genome_file_content(sm):
+    """The genome (chrom sizes) file text for one sample: all bam header
+    contigs, in header order."""
+    return "".join(
+        f"{name}\t{length}\n" for name, length in _bam_contigs(MANIFEST.loc[sm, "bam"])
+    )
+
+
+def remove_stale_genome_files():
+    """Delete leftover genome files that no longer match the bam header.
+
+    The genome_file rule has an ancient() input, so a leftover temp file
+    from an interrupted run is not regenerated when the manifest bam is
+    swapped. Every parse reads the live bam headers anyway, so a
+    mismatched file is removed here, and the missing output then forces
+    regeneration before any consumer runs.
+    """
+    for sm in MANIFEST.index:
+        path = f"temp/{sm}/{sm}.genome"
+        if not os.path.exists(path):
+            continue
+        with open(path) as f:
+            existing = f.read()
+        if existing != genome_file_content(sm):
+            print(f"INFO: {sm}: removing stale genome file {path}", file=sys.stderr)
+            os.remove(path)
+
+
 def get_mem_mb(wildcards, attempt):
     if attempt < 3:
         return attempt * 1024 * 32
diff --git a/workflow/rules/coverages.smk b/workflow/rules/coverages.smk
index ec3b6f955..9fc949259 100644
--- a/workflow/rules/coverages.smk
+++ b/workflow/rules/coverages.smk
@@ -1,10 +1,26 @@
+rule genome_file:
+    input:
+        # ancient for the same reason as rule fire: a re-copied bam with a
+        # fresh mtime must not cascade reruns through every genome-file
+        # consumer. Staleness is impossible: remove_stale_genome_files()
+        # deletes any leftover file that no longer matches the bam header
+        # at parse time, which forces regeneration.
+        bam=lambda wc: ancient(get_input_bam(wc)),
+    output:
+        genome=temp("temp/{sm}/{sm}.genome"),
+    localrule: True
+    run:
+        with open(output.genome, "w") as out:
+            out.write(genome_file_content(wildcards.sm))
+
+
 #
 # Coverage calculations
 #
 rule genome_bedgraph:
     input:
-        ref=ancient(REF),
-        fai=ancient(FAI),
+        ref=lambda wc: ancient(get_ref(wc)),
+        fai=lambda wc: ancient(get_fai(wc)),
         cram=rules.fire.output.cram,
         crai=rules.fire.output.crai,
     output:
@@ -17,9 +33,10 @@ rule genome_bedgraph:
     threads: 16
     shell:
         """
+        # mosdepth output is position sorted in bam header order; keep that
+        # order, it is the order of every other per-sample file
         mosdepth -F 4 -f {input.ref} -t {threads} tmp {input.cram}
         bgzip -cd tmp.per-base.bed.gz \
-            | LC_ALL=C sort --parallel={threads} -k1,1 -k2,2n -k3,3n -k4,4 \
             | bgzip -@ {threads} \
                 >{output.bg}
         tabix -f -p bed {output.bg}
@@ -43,7 +60,7 @@ rule coverage:
     params:
         coverage_within_n_sd=COVERAGE_WITHIN_N_SD,
         min_coverage=MIN_COVERAGE,
-        chroms=get_chroms(),
+        chroms=get_chroms,
     script:
         "../scripts/cov.py"
 
@@ -76,12 +93,13 @@ rule fiber_locations:
     input:
         fibers=expand(
             rules.fiber_locations_chromosome.output.bed,
-            chrom=get_chroms(),
+            chrom=get_chroms,
             allow_missing=True,
         ),
         bg=rules.genome_bedgraph.output.bg,
         minimum=rules.coverage.output.minimum,
         maximum=rules.coverage.output.maximum,
+        genome=rules.genome_file.output.genome,
     output:
         bed=temp("temp/{sm}/coverage/{v}-fiber-locations.bed.gz"),
         bed_tbi=temp("temp/{sm}/coverage/{v}-fiber-locations.bed.gz.tbi"),
@@ -104,7 +122,7 @@ rule fiber_locations:
         # get filtered fiber locations
         MIN=$(cat {input.minimum})
         MAX=$(cat {input.maximum})
-        bedtools intersect -header -sorted -v -f {params.max_frac_overlap} \
+        bedtools intersect -header -sorted -g {input.genome} -v -f {params.max_frac_overlap} \
             -a {output.bed} \
             -b <(bgzip -cd {input.bg} | awk -v MAX="$MAX" -v MIN="$MIN" '$4 <= MIN || $4 >= MAX') \
             | bgzip -@ {threads} \
@@ -119,18 +137,18 @@ rule fiber_locations:
 rule exclude_from_shuffle:
     input:
         filtered=rules.fiber_locations.output.filtered,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bed="results/{sm}/additional-outputs-{v}/coverage/exclude-from-shuffles.bed.gz",
     conda:
         DEFAULT_ENV
     threads: 4
     params:
-        exclude=lambda wc: " ".join(EXCLUDES) if EXCLUDES else "",
+        exclude=lambda wc: " ".join(get_excludes(wc)),
     shell:
         """
         (
-            bedtools genomecov -bga -i {input.filtered} -g {input.fai} | awk '$4 == 0'
+            bedtools genomecov -bga -i {input.filtered} -g {input.genome} | awk '$4 == 0'
             if [ -n "{params.exclude}" ]; then
                 gunzip -cf {params.exclude}
             fi
@@ -148,7 +166,7 @@ rule unreliable_coverage_regions:
         bg=rules.genome_bedgraph.output.bg,
         minimum=rules.coverage.output.minimum,
         maximum=rules.coverage.output.maximum,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bed="results/{sm}/additional-outputs-{v}/coverage/unreliable-coverage-regions.bed.gz",
         bed_tbi="results/{sm}/additional-outputs-{v}/coverage/unreliable-coverage-regions.bed.gz.tbi",
@@ -177,7 +195,7 @@ rule unreliable_coverage_regions:
         bedtools merge -i {output.bed} >{output.tmp}
         bigtools bedtobigbed \
             -s start -a {params.bed3_as} \
-            {output.tmp} {input.fai} {output.bb}
+            {output.tmp} {input.genome} {output.bb}
 
         # index
         tabix -f -p bed {output.bed}
diff --git a/workflow/rules/decorated-reads.smk b/workflow/rules/decorated-reads.smk
index 6a92ee80d..d985f2de7 100644
--- a/workflow/rules/decorated-reads.smk
+++ b/workflow/rules/decorated-reads.smk
@@ -32,10 +32,10 @@ rule decorate_fibers_1:
     input:
         bed=expand(
             rules.decorate_fibers_chromosome.output.bed,
-            chrom=get_chroms(),
+            chrom=get_chroms,
             allow_missing=True,
         ),
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         #bed=temp("temp/{sm}/fiber-calls/fire-fibers.bed.gz"),
         bb="results/{sm}/trackHub-{v}/bb/fire-fibers.bb",
@@ -62,7 +62,7 @@ rule decorate_fibers_1:
                 --block-size {params.block_size} --items-per-slot {params.items_per_slot} \
                 --nzooms {params.nzooms} \
                 -s start -a {params.bed_as} \
-                - {input.fai} {output.bb}
+                - {input.genome} {output.bb}
         """
 
 
@@ -70,10 +70,10 @@ rule decorate_fibers_2:
     input:
         decorated=expand(
             rules.decorate_fibers_chromosome.output.decorated,
-            chrom=get_chroms(),
+            chrom=get_chroms,
             allow_missing=True,
         ),
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bb="results/{sm}/trackHub-{v}/bb/fire-fiber-decorators.bb",
         #bed=temp("temp/{sm}/trackHub-{v}/bb/fire-fiber-decorators.bed.gz"),
@@ -104,7 +104,7 @@ rule decorate_fibers_2:
                 --block-size {params.block_size} --items-per-slot {params.items_per_slot} \
                 --nzooms {params.nzooms} \
                 -s start -a {params.dec_as} \
-                - {input.fai} {output.bb}
+                - {input.genome} {output.bb}
         """
 
 
@@ -114,5 +114,5 @@ if False:
         cat {input.decorated} > {output.bed}
         bedToBigBed \
             -allow1bpOverlap -type=bed12+ -as={params.dec_as} \
-            {output.bed} {input.fai} {output.bb}
+            {output.bed} {input.genome} {output.bb}
         """
diff --git a/workflow/rules/fire-peaks.smk b/workflow/rules/fire-peaks.smk
index fbccff0e5..a0e2d2dbd 100644
--- a/workflow/rules/fire-peaks.smk
+++ b/workflow/rules/fire-peaks.smk
@@ -3,7 +3,7 @@ rule filtered_and_shuffled_fiber_locations_chromosome:
         filtered=rules.fiber_locations.output.filtered,
         filtered_tbi=rules.fiber_locations.output.filtered_tbi,
         exclude=rules.exclude_from_shuffle.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         shuffled=temp("temp/{sm}/shuffle/{v}-{chrom}.fiber-locations-shuffled.bed.gz"),
     conda:
@@ -16,7 +16,7 @@ rule filtered_and_shuffled_fiber_locations_chromosome:
             | bedtools shuffle -chrom -seed 42 \
                 -excl {input.exclude} \
                 -i - \
-                -g {input.fai} \
+                -g {input.genome} \
             | sort -k1,1 -k2,2n -k3,3n -k4,4 \
             | bgzip -@ {threads} \
                 >{output.shuffled}
@@ -46,7 +46,7 @@ rule shuffled_pileup:
     input:
         beds=expand(
             rules.shuffled_pileup_chromosome.output.bed,
-            chrom=get_chroms(),
+            chrom=get_chroms,
             allow_missing=True,
         ),
     output:
@@ -137,7 +137,7 @@ rule pileup:
     input:
         beds=expand(
             rules.fdr_track_chromosome.output.bed,
-            chrom=get_chroms(),
+            chrom=get_chroms,
             allow_missing=True,
         ),
     output:
@@ -250,7 +250,7 @@ rule fire_peaks:
     input:
         beds=expand(
             rules.fdr_peaks_by_fire_elements_chromosome.output.bed,
-            chrom=get_chroms(),
+            chrom=get_chroms,
             allow_missing=True,
         ),
     output:
@@ -282,7 +282,7 @@ rule wide_fire_peaks:
     input:
         bed=rules.fire_peaks.output.bed,
         track=rules.pileup.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bed="results/{sm}/additional-outputs-{v}/fire-peaks/{sm}-fire-{v}-wide-peaks.bed.gz",
         tbi="results/{sm}/additional-outputs-{v}/fire-peaks/{sm}-fire-{v}-wide-peaks.bed.gz.tbi",
@@ -311,7 +311,7 @@ rule wide_fire_peaks:
         bgzip -cd -@ 16 {output.bed} \
             | bigtools bedtobigbed \
                 -s start -a {params.bed3_as} \
-                - {input.fai} {output.bb}
+                - {input.genome} {output.bb}
 
         tabix -p bed {output.bed}
         """
diff --git a/workflow/rules/stats.smk b/workflow/rules/stats.smk
index f462c8592..56cdc0764 100644
--- a/workflow/rules/stats.smk
+++ b/workflow/rules/stats.smk
@@ -4,7 +4,7 @@
 rule clustering_vs_null:
     input:
         bed=rules.fire_sites.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         tmp=temp("temp/{sm}/tmp.pre.calls.bed"),
         null=temp("temp/{sm}/null.calls.bed"),
@@ -15,11 +15,11 @@ rule clustering_vs_null:
     shell:
         """
         bgzip -cd -@{threads} {input.bed} | cut -f 1-3 >{output.tmp}
-        bedtools shuffle -chrom -i {output.tmp} -g {input.fai} >{output.null}
+        bedtools shuffle -chrom -i {output.tmp} -g {input.genome} >{output.null}
 
         (
-            bedtools genomecov -bg -i {output.tmp} -g {input.fai} | sed 's/$/\\tReal/g'
-            bedtools genomecov -bg -i {output.null} -g {input.fai} | sed 's/$/\\tNull/g'
+            bedtools genomecov -bg -i {output.tmp} -g {input.genome} | sed 's/$/\\tReal/g'
+            bedtools genomecov -bg -i {output.null} -g {input.genome} | sed 's/$/\\tNull/g'
         ) \
             | bedtools sort \
             | bgzip -@ {threads} \
diff --git a/workflow/rules/track-hub.smk b/workflow/rules/track-hub.smk
index 1be0fcb55..8dd036a6b 100644
--- a/workflow/rules/track-hub.smk
+++ b/workflow/rules/track-hub.smk
@@ -1,7 +1,7 @@
 rule percent_accessible:
     input:
         bed=rules.pileup.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         tmp=temp("temp/{sm}/{hp}/{v}-percent.accessible.bed"),
         bw="results/{sm}/trackHub-{v}/bw/{hp}.percent.accessible.bw",
@@ -13,7 +13,7 @@ rule percent_accessible:
     params:
         suffix=get_hap_col_suffix,
         nzooms=NZOOMS,
-        chrom=get_chroms()[0],
+        chrom=lambda wc: get_chroms(wc)[0],
     shell:
         """
         bgzip -cd {input.bed} \
@@ -32,14 +32,14 @@ rule percent_accessible:
 
         bigtools bedgraphtobigwig \
             --nzooms {params.nzooms} -s start \
-            {output.tmp} {input.fai} {output.bw}
+            {output.tmp} {input.genome} {output.bw}
         """
 
 
 rule element_coverages_bw:
     input:
         bed=rules.pileup.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bw="results/{sm}/trackHub-{v}/bw/{hp}.{el_type}.coverage.bw",
     conda:
@@ -54,14 +54,14 @@ rule element_coverages_bw:
             | grep -v "^#" \
             | bigtools bedgraphtobigwig \
                 -s start --nzooms {params.nzooms} \
-                - {input.fai} {output.bw}
+                - {input.genome} {output.bw}
         """
 
 
 rule fdr_track_to_bw:
     input:
         bed=rules.pileup.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bw="results/{sm}/trackHub-{v}/bw/{col}.bw",
     conda:
@@ -75,14 +75,14 @@ rule fdr_track_to_bw:
             | grep -v "^#" \
             | bigtools bedgraphtobigwig \
                 -s start --nzooms {params.nzooms} \
-                - {input.fai} {output.bw}
+                - {input.genome} {output.bw}
         """
 
 
 rule fire_peaks_bb:
     input:
         bed=rules.fire_peaks.output.bed,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bb="results/{sm}/trackHub-{v}/bb/fire-peaks.bb",
     conda:
@@ -98,14 +98,14 @@ rule fire_peaks_bb:
             | rg -v '^#' \
             | bigtools bedtobigbed \
                 -a {params.bedfmt} -s start \
-                - {input.fai} {output.bb}
+                - {input.genome} {output.bb}
         """
 
 
 rule hap_differences_track:
     input:
         bed9=rules.hap_differences.output.bed9,
-        fai=ancient(FAI),
+        genome=rules.genome_file.output.genome,
     output:
         bb="results/{sm}/trackHub-{v}/bb/hap_differences.bb",
     conda:
@@ -114,7 +114,7 @@ rule hap_differences_track:
     resources:
         mem_mb=get_mem_mb,
     params:
-        chrom=get_chroms()[0],
+        chrom=lambda wc: get_chroms(wc)[0],
         bed9_as=workflow.source_path("../templates/bed9.as"),
     shell:
         """
@@ -124,7 +124,7 @@ rule hap_differences_track:
         ) \
             | bigtools bedtobigbed \
                 -s start -a {params.bed9_as} \
-                - {input.fai} {output.bb}
+                - {input.genome} {output.bb}
         """
 
 
@@ -140,7 +140,7 @@ rule trackhub:
     resources:
         load=get_load,
     params:
-        ref=REF_NAME,
+        ref=get_ref_name,
         script=workflow.source_path("../scripts/trackhub.py"),
         description=workflow.source_path("../templates/fire-description.html"),
     shell: