From 289cfbe97c51b4816dd490c57e79fb97536f4041 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Thu, 30 Jul 2026 15:20:52 +0200 Subject: [PATCH 1/6] perf(score): slide the splice-motif window across the junction scan `detect_splice_motif` and `find_best_junction_position` were the two largest self-time entries in a sampled 2M-read run, together about 45% of the on-CPU samples. The scan moves the junction one base per iteration, so the four bases that decide the motif (the intron's first two and last two) overlap the previous position in two of four places. Carrying them in a `MotifWindow` and sliding turns four genome reads per position into two. Whether the motif branch runs at all is decided once before the loop rather than per iteration, since `del` does not change across the scan. An out-of-range position becomes a sentinel that no motif arm matches, which is what `get_base` returning `None` already meant. Output-neutral: SAM byte-identical to the parent commit on 200k real reads. Interleaved A/B at 16 threads on 2M reads, machine at 87-95% CPU idle: 22.94s median before, 22.29s after, the same direction in all four rounds. About 2.8%. Co-Authored-By: Claude Opus 5 (1M context) --- src/align/score.rs | 113 ++++++++++++++++++++++++++++++++++++++------- 1 file changed, 97 insertions(+), 16 deletions(-) diff --git a/src/align/score.rs b/src/align/score.rs index 7a34f55..0e89881 100644 --- a/src/align/score.rs +++ b/src/align/score.rs @@ -353,6 +353,18 @@ impl AlignmentScorer { let mut best_motif = SpliceMotif::NonCanonical; let mut best_motif_score = self.score_gap_noncan; + // `del` does not change across the scan, so whether the motif branch + // runs is decided once rather than per iteration. + let motif_in_range = + del >= self.align_intron_min as i64 && del <= self.align_intron_max as i64; + + // The four motif bases sit at `donor`, `donor+1`, `donor+del-2` and + // `donor+del-1`, and `donor` moves exactly one base per iteration, so + // consecutive iterations share two of the four. Carrying the window + // instead of re-fetching it turns four genome reads per position into + // two. `window` is `None` until the first motif iteration primes it. + let mut window: Option = None; + loop { let ri = r_a_end_inc as i64 + jr1 as i64; if ri >= 0 && (ri as usize) < read_seq.len() { @@ -378,7 +390,7 @@ impl AlignmentScorer { } // Check splice motif at this junction position - if del >= self.align_intron_min as i64 && del <= self.align_intron_max as i64 { + if motif_in_range { // Donor position in SA space: one past the last donor-exon base let donor_sa = (g_a_end_inc as i64 + jr1 as i64 + 1) as u64; // Convert to forward genome coordinates for motif detection @@ -387,7 +399,14 @@ impl AlignmentScorer { } else { donor_sa }; - let motif = self.detect_splice_motif(donor_fwd, del as u32, genome); + let w = match window.as_mut() { + Some(w) => { + w.slide_to(donor_fwd, del as u64, genome); + &*w + } + None => window.insert(MotifWindow::at(donor_fwd, del as u64, genome)), + }; + let motif = w.motif(); let motif_score = self.score_splice_junction(motif); let score2 = score1 + motif_score; @@ -535,20 +554,82 @@ impl AlignmentScorer { /// `donor_pos` is the 0-based position of the intron's first base on the /// forward strand; `intron_len` is the intron length in bases. pub fn detect_splice_motif(donor_pos: u64, intron_len: u32, genome: &Genome) -> SpliceMotif { - let d1 = genome.get_base(donor_pos); - let d2 = genome.get_base(donor_pos + 1); - let a1 = genome.get_base(donor_pos + intron_len as u64 - 2); - let a2 = genome.get_base(donor_pos + intron_len as u64 - 1); - - // Base encoding: A=0, C=1, G=2, T=3. - match (d1, d2, a1, a2) { - (Some(2), Some(3), Some(0), Some(2)) => SpliceMotif::GtAg, - (Some(2), Some(1), Some(0), Some(2)) => SpliceMotif::GcAg, - (Some(0), Some(3), Some(0), Some(1)) => SpliceMotif::AtAc, - (Some(1), Some(3), Some(0), Some(1)) => SpliceMotif::CtAc, - (Some(1), Some(3), Some(2), Some(1)) => SpliceMotif::CtGc, - (Some(2), Some(3), Some(0), Some(3)) => SpliceMotif::GtAt, - _ => SpliceMotif::NonCanonical, + MotifWindow::at(donor_pos, intron_len as u64, genome).motif() +} + +/// A position off the end of the genome. No motif arm matches it, so it falls +/// through to `NonCanonical` exactly as `get_base` returning `None` did. +const OUT_OF_RANGE: u8 = u8::MAX; + +#[inline] +fn base_or_out_of_range(genome: &Genome, pos: u64) -> u8 { + genome.get_base(pos).unwrap_or(OUT_OF_RANGE) +} + +/// The four bases that decide a splice motif: the intron's first two and last +/// two, on the forward strand. +/// +/// Kept as a struct so the junction scan can slide it. Successive junction +/// positions differ by one base, so three of the four positions overlap the +/// previous window and only two bases have to be read from the genome. +#[derive(Clone, Copy)] +struct MotifWindow { + donor: u64, + d1: u8, + d2: u8, + a1: u8, + a2: u8, +} + +impl MotifWindow { + #[inline] + fn at(donor: u64, intron_len: u64, genome: &Genome) -> Self { + Self { + donor, + d1: base_or_out_of_range(genome, donor), + d2: base_or_out_of_range(genome, donor + 1), + a1: base_or_out_of_range(genome, donor + intron_len - 2), + a2: base_or_out_of_range(genome, donor + intron_len - 1), + } + } + + /// Move the window to `donor`, reusing what overlaps. + /// + /// A step of exactly one base either way shares two of the four: moving + /// right, the old `d2`/`a2` become the new `d1`/`a1`; moving left, the old + /// `d1`/`a1` become the new `d2`/`a2`. Any other step is rare enough that + /// re-reading all four is the simpler answer. + #[inline] + fn slide_to(&mut self, donor: u64, intron_len: u64, genome: &Genome) { + if donor == self.donor + 1 { + self.d1 = self.d2; + self.a1 = self.a2; + self.d2 = base_or_out_of_range(genome, donor + 1); + self.a2 = base_or_out_of_range(genome, donor + intron_len - 1); + self.donor = donor; + } else if donor + 1 == self.donor { + self.d2 = self.d1; + self.a2 = self.a1; + self.d1 = base_or_out_of_range(genome, donor); + self.a1 = base_or_out_of_range(genome, donor + intron_len - 2); + self.donor = donor; + } else if donor != self.donor { + *self = Self::at(donor, intron_len, genome); + } + } + + /// Base encoding: A=0, C=1, G=2, T=3. + #[inline] + fn motif(&self) -> SpliceMotif { + match (self.d1, self.d2, self.a1, self.a2) { + (2, 3, 0, 2) => SpliceMotif::GtAg, + (2, 1, 0, 2) => SpliceMotif::GcAg, + (0, 3, 0, 1) => SpliceMotif::AtAc, + (1, 3, 0, 1) => SpliceMotif::CtAc, + (1, 3, 2, 1) => SpliceMotif::CtGc, + (2, 3, 0, 3) => SpliceMotif::GtAt, + _ => SpliceMotif::NonCanonical, + } } } From 0c0ec815920ab8f95ee97d1c155acb6d6b78ebe5 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Thu, 30 Jul 2026 20:00:04 +0200 Subject: [PATCH 2/6] perf(score): look the splice motif up in a table instead of matching The four-base motif decision compiled to a chain of comparisons whose outcome is data-dependent and close to unpredictable, and the junction scan pays it at every position. A 256-entry table indexed by the four bases packed two bits each replaces it with one load. `examples/jrbench.rs` measures the scan in isolation: 8.5 ns per iteration before, 5.2 ns after, a 38% cut, with identical results. That example exists because the earlier attempts on this path were guesses about whether the loop was memory-bound, and two of them were wrong. It answers the question by construction, holding the iteration count and the instruction mix fixed and changing only how far the acceptor sits from the donor: del ns/iter 64 8.42 1024 8.47 16384 8.52 262144 8.41 4194304 8.62 67108864 8.56 Flat from 64 bases to 67 million, so the acceptor distance costs nothing and the loop is not stalled on that stream. Which is why skipping the acceptor read behind a donor-pair test measured *slower*: it traded a prefetchable access for an unpredictable branch. The branch was the cost all along, and this removes it. `the_motif_table_agrees_with_the_original_match_on_every_input` checks the table against the match it replaces over every combination of `A`/`C`/`G`/`T`, `N`, the chromosome-boundary byte and the out-of-range sentinel, so the inputs that no match arm covered are covered explicitly. Output-neutral: SAM byte-identical on 200k real reads. Co-Authored-By: Claude Opus 5 (1M context) --- examples/jrbench.rs | 97 +++++++++++++++++++++++++++++++++++++++++++++ src/align/score.rs | 82 ++++++++++++++++++++++++++++++++++---- 2 files changed, 171 insertions(+), 8 deletions(-) create mode 100644 examples/jrbench.rs diff --git a/examples/jrbench.rs b/examples/jrbench.rs new file mode 100644 index 0000000..d987c9b --- /dev/null +++ b/examples/jrbench.rs @@ -0,0 +1,97 @@ +//! Is the junction scan stalled on memory, and is it the acceptor stream? +//! +//! The stack sampler gives self time, not cache misses, so this asks the +//! question by construction instead: run the same scan over the same genome, +//! changing only how far the acceptor sits from the donor. Everything else +//! (iteration count, branch pattern, instruction mix) is held fixed, because +//! `del` enters the loop only as an address offset once it is inside the +//! intron-length range. +//! +//! If time per iteration is flat across `del`, the loop is not memory-bound on +//! the acceptor and prefetching it would be wasted work. If it climbs with +//! `del`, that stream is the cost. +//! +//! Run: cargo run --release --example jrbench + +use rustar_aligner::align::score::AlignmentScorer; +use rustar_aligner::genome::{Genome, GenomeSeq}; +use std::time::Instant; + +/// Deterministic pseudo-random bases, so the run is reproducible and the motif +/// hit rate is the same at every `del`. +fn synthetic_genome(n: usize) -> Genome { + let mut state = 0x2545_F491_4F6C_DD1Du64; + let mut seq = Vec::with_capacity(n); + for _ in 0..n { + state ^= state << 13; + state ^= state >> 7; + state ^= state << 17; + seq.push((state >> 33) as u8 & 3); + } + Genome { + transform_blocks: None, + sequence: GenomeSeq::Owned(seq), + n_genome: n as u64, + n_genome_real: n as u64, + n_chr_real: 1, + chr_name: vec!["chr1".to_string()], + chr_length: vec![n as u64], + chr_start: vec![0, n as u64], + } +} + +fn main() { + // Large enough that the donor and acceptor streams cannot both stay in + // cache when they are far apart. + const N: usize = 256 << 20; // 256 MB of bases, one byte each + const READ_LEN: usize = 150; + const SCANS: usize = 20_000; + + let genome = synthetic_genome(N); + let mut scorer = AlignmentScorer::from_params_minimal(); + scorer.align_intron_min = 20; + scorer.align_intron_max = u32::MAX; + + let mut read = vec![0u8; READ_LEN]; + let mut state = 0x9E37_79B9_7F4A_7C15u64; + for b in read.iter_mut() { + state ^= state << 13; + state ^= state >> 7; + state ^= state << 17; + *b = (state >> 33) as u8 & 3; + } + + println!("{:>12} {:>10} {:>12}", "del", "wall", "ns/iter"); + for del in [64u64, 1_024, 16_384, 262_144, 4_194_304, 67_108_864] { + // Spread the scans over the genome so no single window stays resident. + let stride = (N as u64 - del - 4096) / SCANS as u64; + let iters_per_scan = 100usize; // r_gap 0 + next_seed_len 100 + let t0 = Instant::now(); + let mut sink = 0i64; + for k in 0..SCANS { + let g_a_end = 2048 + k as u64 * stride; + let (jr, _motif, score, _l, _r) = scorer.find_best_junction_position( + &read, + 75, + g_a_end, + 0, + del as i64, + &genome, + false, + N as u64, + 75, + iters_per_scan, + ); + sink += jr as i64 + score as i64; + } + let dt = t0.elapsed(); + let iters = (SCANS * iters_per_scan) as f64; + println!( + "{:>12} {:>9.3}s {:>11.2} (sink {})", + del, + dt.as_secs_f64(), + dt.as_secs_f64() * 1e9 / iters, + sink + ); + } +} diff --git a/src/align/score.rs b/src/align/score.rs index 0e89881..408a2d2 100644 --- a/src/align/score.rs +++ b/src/align/score.rs @@ -619,20 +619,45 @@ impl MotifWindow { } /// Base encoding: A=0, C=1, G=2, T=3. + /// + /// A table lookup rather than a match on the four bases. The match compiled + /// to a chain of compares whose outcome is data-dependent and close to + /// unpredictable, which the scan pays at every junction position; the table + /// is one load at a computed index. `|` binds tighter than `>=` in Rust, so + /// the guard tests the union of the four bases and rejects anything that is + /// not `A`, `C`, `G` or `T` — `N`, the chromosome-boundary byte, and the + /// out-of-range sentinel all take that path, exactly as no match arm + /// covered them. #[inline] fn motif(&self) -> SpliceMotif { - match (self.d1, self.d2, self.a1, self.a2) { - (2, 3, 0, 2) => SpliceMotif::GtAg, - (2, 1, 0, 2) => SpliceMotif::GcAg, - (0, 3, 0, 1) => SpliceMotif::AtAc, - (1, 3, 0, 1) => SpliceMotif::CtAc, - (1, 3, 2, 1) => SpliceMotif::CtGc, - (2, 3, 0, 3) => SpliceMotif::GtAt, - _ => SpliceMotif::NonCanonical, + let (d1, d2, a1, a2) = (self.d1, self.d2, self.a1, self.a2); + if d1 | d2 | a1 | a2 >= 4 { + return SpliceMotif::NonCanonical; } + MOTIF_TABLE[motif_index(d1, d2, a1, a2)] } } +/// Every four-base combination, packed `d1 d2 a1 a2` at two bits each. +static MOTIF_TABLE: [SpliceMotif; 256] = build_motif_table(); + +/// Pack the four bases into the table index, two bits each. +const fn motif_index(d1: u8, d2: u8, a1: u8, a2: u8) -> usize { + ((d1 << 6) | (d2 << 4) | (a1 << 2) | a2) as usize +} + +const fn build_motif_table() -> [SpliceMotif; 256] { + // A=0, C=1, G=2, T=3. + let mut t = [SpliceMotif::NonCanonical; 256]; + t[motif_index(2, 3, 0, 2)] = SpliceMotif::GtAg; + t[motif_index(2, 1, 0, 2)] = SpliceMotif::GcAg; + t[motif_index(0, 3, 0, 1)] = SpliceMotif::AtAc; + t[motif_index(1, 3, 0, 1)] = SpliceMotif::CtAc; + t[motif_index(1, 3, 2, 1)] = SpliceMotif::CtGc; + t[motif_index(2, 3, 0, 3)] = SpliceMotif::GtAt; + t +} + /// Splice junction motif types #[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)] pub enum SpliceMotif { @@ -705,6 +730,47 @@ mod tests { use super::*; + /// The table has to agree with the match it replaced on every input the + /// scan can present, not merely on the six motifs. That includes `N` (4), + /// the chromosome-boundary byte (5) and the out-of-range sentinel, all of + /// which no match arm covered and which must therefore be `NonCanonical`. + #[test] + fn the_motif_table_agrees_with_the_original_match_on_every_input() { + fn original(d1: u8, d2: u8, a1: u8, a2: u8) -> SpliceMotif { + match (d1, d2, a1, a2) { + (2, 3, 0, 2) => SpliceMotif::GtAg, + (2, 1, 0, 2) => SpliceMotif::GcAg, + (0, 3, 0, 1) => SpliceMotif::AtAc, + (1, 3, 0, 1) => SpliceMotif::CtAc, + (1, 3, 2, 1) => SpliceMotif::CtGc, + (2, 3, 0, 3) => SpliceMotif::GtAt, + _ => SpliceMotif::NonCanonical, + } + } + + let values = [0u8, 1, 2, 3, 4, 5, OUT_OF_RANGE]; + for &d1 in &values { + for &d2 in &values { + for &a1 in &values { + for &a2 in &values { + let w = MotifWindow { + donor: 0, + d1, + d2, + a1, + a2, + }; + assert_eq!( + w.motif(), + original(d1, d2, a1, a2), + "disagreement at ({d1}, {d2}, {a1}, {a2})" + ); + } + } + } + } + } + fn make_test_genome(seq: &[u8]) -> Genome { // Create simple genome with one chromosome let n_genome = ((seq.len() as u64 + 1) / 64 + 1) * 64; // Pad to 64-byte boundary From 961366fc0cc04feee3c69dc634e5282f860bcd5a Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Thu, 30 Jul 2026 20:12:05 +0200 Subject: [PATCH 3/6] docs(changelog): record the junction-scan motif work Co-Authored-By: Claude Opus 5 (1M context) --- CHANGELOG.md | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 087aa49..fb43f0d 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -186,10 +186,16 @@ Sections commonly used: Features, Bug fixes, Other changes. `GenomeIndex::write` flow remains for tests and any caller that needs random access to the SA in RAM. -Initial release of Rust rewrite of STAR. -### Other changes - - Removed `Transcript::read_seq`, a public field that was filled with a copy of the read at every finalised alignment and never read. **API removal.** Output is unchanged. +- The splice-motif check in the junction scan carries a sliding window + and looks the motif up in a table, instead of re-reading four genome + bases and matching on them at every position. Consecutive junction + positions share two of the four bases, so the scan does two genome + reads per position rather than four. Output is unchanged (SAM + byte-identical on 200k reads); about 2.8% off the wall clock on a + 2M-read run, measured with `test/bench_ab.sh`. + +Initial release of Rust rewrite of STAR. From ba152cb1b16e89efa3b6ee17b992af484a991eef Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Fri, 31 Jul 2026 21:29:40 +0200 Subject: [PATCH 4/6] perf(align): reuse cluster_seeds' window-bin map across reads cluster_seeds built a fresh FxHashMap per read. The pre-sizing to anchor_indices.len() * 2 was only a floor: merging two windows re-keys every bin in the merged span, so a read with wide windows inserts far more entries than it has anchors and the map rehashes. Sampling a human 10x run put hashbrown::reserve_rehash at 2.1% of on-CPU time, all of it here. The map is now kept per thread and cleared per read, so its capacity settles at what the workload needs and the growth is paid once per thread. It is taken out of the thread-local rather than borrowed across the body, so a nested call would get a fresh map instead of panicking on an active borrow. Output-neutral: matrix.mtx, barcodes.tsv, features.tsv and SJ.out.tab are byte-identical before and after on 20 M read pairs of 10x pbmc_1k_v3 against GRCh38, compared decompressed. Interleaved A/B, 8 pairs: new faster in 8 of 8, median -7.0%. Co-Authored-By: Claude Opus 5 (1M context) --- CHANGELOG.md | 8 ++++++++ src/align/stitch.rs | 27 +++++++++++++++++++++++---- 2 files changed, 31 insertions(+), 4 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index fb43f0d..280323c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,14 @@ Sections commonly used: Features, Bug fixes, Other changes. ## [Unreleased] +### Other changes + +- `cluster_seeds` reuses its window-bin map across reads on a thread instead + of rebuilding it per read. Merging two windows re-keys every bin in the + merged span, so the per-read pre-sizing was only a floor and the map + rehashed; profiling a human 10x run put that rehash at 2.1% of on-CPU time. + Output-neutral, verified by an empty diff on 20 M read pairs. + ### Features - **STARsolo single-cell quantification (`--soloType`)** — the 10x diff --git a/src/align/stitch.rs b/src/align/stitch.rs index 349ff40..bfae6e2 100644 --- a/src/align/stitch.rs +++ b/src/align/stitch.rs @@ -461,7 +461,7 @@ pub fn cluster_seeds( ) -> Vec { // Integer-keyed maps in this hot path (the #1 align hotspot, ~19% self-time): // FxHash, not the default SipHash, and pre-sized to avoid rehashing. - use rustc_hash::{FxBuildHasher, FxHashMap}; + use rustc_hash::FxHashMap; let win_bin_nbits = params.win_bin_nbits; let win_anchor_dist_nbins = params.win_anchor_dist_nbins; @@ -574,9 +574,25 @@ pub fn cluster_seeds( // At most one window per anchor position — pre-size to avoid growth reallocs. let mut windows: Vec = Vec::with_capacity(anchor_indices.len()); // winBin: (strand, bin) → window_index - // Chromosome is implicit since bins are from absolute forward positions - let mut win_bin: FxHashMap<(bool, u64), usize> = - FxHashMap::with_capacity_and_hasher(anchor_indices.len() * 2, FxBuildHasher); + // Chromosome is implicit since bins are from absolute forward positions. + // + // Reused across reads on this thread. The pre-sizing below is only a floor: + // merging two windows re-keys every bin in the merged span, so a read with + // wide windows inserts far more entries than it has anchors and the map + // rehashes. Profiling a human 10x run put `hashbrown::reserve_rehash` at + // 2.1% of on-CPU time, all of it here. Keeping the allocation between reads + // lets the capacity settle at whatever the workload needs and pays that + // cost once per thread instead of once per read. + // + // Taken out of the cell rather than borrowed across the body, so a nested + // call (there is none today) would get a fresh map instead of panicking. + thread_local! { + static WIN_BIN: std::cell::RefCell> = + std::cell::RefCell::new(FxHashMap::default()); + } + let mut win_bin = WIN_BIN.with(|c| std::mem::take(&mut *c.borrow_mut())); + win_bin.clear(); + win_bin.reserve(anchor_indices.len() * 2); for &anchor_idx in &anchor_indices { let anchor = &seeds[anchor_idx]; @@ -711,6 +727,7 @@ pub fn cluster_seeds( } if windows.iter().all(|w| !w.alive) { + WIN_BIN.with(|c| *c.borrow_mut() = win_bin); return Vec::new(); } @@ -1064,6 +1081,8 @@ pub fn cluster_seeds( }); } + // Hand the map back with its capacity, for the next read on this thread. + WIN_BIN.with(|c| *c.borrow_mut() = win_bin); clusters } From 91a848a5816c4cba761cd69a9c037a3b61f5764c Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Fri, 31 Jul 2026 22:31:19 +0200 Subject: [PATCH 5/6] perf(align): move cluster_seeds' alignments into their cluster The per-window Vec was cloned into the SeedCluster while `windows` was dropped one statement later, so every read paid a full copy per window to throw the original away. std::mem::take moves it instead. Output-neutral, verified by an empty diff on 20 M read pairs of 10x pbmc_1k_v3 against GRCh38. Co-Authored-By: Claude Opus 5 (1M context) --- src/align/stitch.rs | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/src/align/stitch.rs b/src/align/stitch.rs index bfae6e2..8dcca5c 100644 --- a/src/align/stitch.rs +++ b/src/align/stitch.rs @@ -1065,13 +1065,17 @@ pub fn cluster_seeds( // Phase 5: Build SeedCluster output let mut clusters = Vec::with_capacity(windows.len()); - for window in &windows { + // `windows` is dropped at the end of this function, so each window's + // alignments move into its cluster rather than being cloned. The clone was + // a full Vec copy per window per read, thrown away one + // statement later. + for window in &mut windows { if !window.alive || window.alignments.is_empty() { continue; } clusters.push(SeedCluster { - alignments: window.alignments.clone(), + alignments: std::mem::take(&mut window.alignments), chr_idx: window.chr_idx, genome_start: window.actual_start, genome_end: window.actual_end, From da18148733de08e22ff6d41644abf5a6f9df2857 Mon Sep 17 00:00:00 2001 From: Psy-Fer Date: Fri, 7 Aug 2026 12:59:14 +1000 Subject: [PATCH 6/6] fix(test): bench_ab.sh could not run on Linux, or with relative paths --- test/bench_ab.sh | 24 +++++++++++++++++++++++- 1 file changed, 23 insertions(+), 1 deletion(-) diff --git a/test/bench_ab.sh b/test/bench_ab.sh index d10755b..d364c48 100755 --- a/test/bench_ab.sh +++ b/test/bench_ab.sh @@ -48,6 +48,11 @@ BIN_A=${1:?usage: bench_ab.sh [threads] BIN_B=${2:?usage: bench_ab.sh [threads] [rounds]} GENOME_DIR=${3:?usage: bench_ab.sh [threads] [rounds]} READS=${4:?usage: bench_ab.sh [threads] [rounds]} +# Each round runs from a private temp directory, so every path handed to the +# aligner has to be absolute or it will not resolve from there. +GENOME_DIR=$(cd "$GENOME_DIR" && pwd) +READS=$(cd "$(dirname "$READS")" && pwd)/$(basename "$READS") + THREADS=${5:-16} ROUNDS=${6:-6} @@ -69,7 +74,24 @@ cp "$BIN_A" "$WORK/a/rustar-aligner" cp "$BIN_B" "$WORK/b/rustar-aligner" cpu_idle() { # percent idle, sampled over one second - top -l 2 -n 0 2>/dev/null | awk '/CPU usage/ {gsub("%","",$(NF-1)); v=$(NF-1)} END {print v+0}' + # macOS and Linux disagree on how to ask. `top -l` is macOS-only; on Linux it + # exits non-zero, and under `set -o pipefail` that aborts the whole script + # with no message, so this must branch rather than rely on a fallback. + if [ "$(uname -s)" = "Darwin" ]; then + top -l 2 -n 0 2>/dev/null | awk '/CPU usage/ {gsub("%","",$(NF-1)); v=$(NF-1)} END {print v+0}' + return + fi + # Linux: two /proc/stat samples a second apart. Field 5 is idle; the total is + # every field, so the ratio is idle time over elapsed jiffies across all CPUs. + local s1 s2 + s1=$(awk '/^cpu /{t=0; for (j=2; j<=NF; j++) t+=$j; print $5, t; exit}' /proc/stat) + sleep 1 + s2=$(awk '/^cpu /{t=0; for (j=2; j<=NF; j++) t+=$j; print $5, t; exit}' /proc/stat) + awk -v a="$s1" -v b="$s2" 'BEGIN { + split(a, x, " "); split(b, y, " "); + di = y[1] - x[1]; dt = y[2] - x[2]; + print (dt > 0) ? int(100 * di / dt + 0.5) : 100 + }' } # Returns " ".