From ff6cbe851009a128b64f33caae45566bc0f79ec0 Mon Sep 17 00:00:00 2001 From: Benjamin Riley Date: Sun, 19 Jul 2026 22:01:09 +0000 Subject: [PATCH 1/2] feat: implement OpenTfRawSource::iter_chromatograms for TIC/BPC/SRM Wire the already-decoded per-scan chromatogram fields into openmassspec_core::ChromatogramRecords so TIC, BPC, and SRM traces reach mzML output through the shared writer. Previously OpenTfRawSource fell through to the trait's empty default. - TIC (MS:1000235) and BPC (MS:1000628) are built from the scan index's total_current / base_intensity fields (retention time minutes -> seconds). - For flat-peak SRM files (v63/v66), scans are grouped by scan_event into one selected reaction monitoring chromatogram (MS:1000627) per transition, with precursor_mz from the Q1 map and product_mz set only when a single Q3 window resolves unambiguously; transitions without a known Q1 are skipped. No new binary parsing is introduced. Added pure-function unit tests for the TIC/BPC/SRM mapping. Fixes Sigilweaver/OpenTFRaw#39 Claude-Session: https://claude.ai/code/session_01TLSifJpmaN7RvW3xi2uWgo --- CHANGELOG.md | 11 ++ CONTRIBUTORS.md | 7 + crates/opentfraw/src/mzml.rs | 239 +++++++++++++++++++++++++++++++++++ 3 files changed, 257 insertions(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index a89d208..8a5ad55 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -31,6 +31,17 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added +- `OpenTfRawSource::iter_chromatograms` (Sigilweaver/OpenTFRaw#39): wires the + already-decoded per-scan chromatogram fields into + `openmassspec_core::ChromatogramRecord`s so TIC, BPC, and SRM traces reach + mzML output via the shared `` writer. Emits a total ion + current chromatogram (`MS:1000235`) and a base-peak chromatogram + (`MS:1000628`) from the scan index's `total_current` / `base_intensity` + fields, and, for flat-peak SRM files (v63/v66), one selected reaction + monitoring chromatogram (`MS:1000627`) per transition, grouping scans by + `scan_event` with `precursor_mz` from the Q1 map and `product_mz` filled in + only when a single Q3 window resolves unambiguously. No new binary parsing + is involved. (contributed by @Nabejo) - Unit tests for the pure-decode functions in `bytes.rs`, `header.rs`, `audit_tag.rs`, `generic_data.rs`, and `error_log.rs`, using hand-built byte fixtures (valid, truncated, and out-of-range cases). Also added diff --git a/CONTRIBUTORS.md b/CONTRIBUTORS.md index 12393f9..62ef093 100644 --- a/CONTRIBUTORS.md +++ b/CONTRIBUTORS.md @@ -4,6 +4,13 @@ Thank you to everyone who has contributed to OpenTFRaw. ## Benjamin Riley ([@Nabejo](https://github.com/Nabejo)) +Contributed in [Unreleased]: + +- **`OpenTfRawSource::iter_chromatograms`** - wired the already-decoded + per-scan TIC / base-peak fields (and, for flat-peak SRM files, the Q1/Q3 + transition maps) into `openmassspec_core::ChromatogramRecord`s so TIC, BPC, + and SRM chromatograms reach mzML output. (Sigilweaver/OpenTFRaw#39) + Contributed in v1.3.0: - **`RawFile.sample_info`** - exposed the sample-sheet / sequence-row diff --git a/crates/opentfraw/src/mzml.rs b/crates/opentfraw/src/mzml.rs index f8990e0..2452e40 100644 --- a/crates/opentfraw/src/mzml.rs +++ b/crates/opentfraw/src/mzml.rs @@ -489,6 +489,142 @@ fn to_msc_record(rec: SpectrumRecord) -> msc::SpectrumRecord { } } +// -- Chromatogram records (TIC / BPC / SRM) -- +// +// These map already-decoded per-scan fields into vendor-neutral +// `openmassspec_core::ChromatogramRecord`s. No new binary parsing happens +// here: the TIC/BPC series come straight from the scan index +// (`total_current`, `base_intensity`), and SRM transitions are formed by +// grouping the same scans by their `scan_event` using the Q1/Q3 maps the +// reader already builds. Following the sibling Waters wiring +// (Sigilweaver/OpenWRaw#9), a trace is only emitted when it can be labeled +// with a real PSI-MS chromatogram-type term and its defining fields resolve; +// anything ambiguous is omitted rather than guessed. + +/// Build the total-ion-current chromatogram from `(rt_min, total_current)` +/// pairs, or `None` when the file has no scans. +/// +/// `ChromatogramRecord::time_sec` is seconds; the reader reports retention +/// time in minutes, hence the `* 60.0`. +fn tic_record(tic: &[(f64, f64)]) -> Option { + if tic.is_empty() { + return None; + } + Some(msc::ChromatogramRecord { + index: 0, + id: "TIC".to_string(), + chromatogram_type: Some(msc::CvTerm::new( + "MS:1000235", + "total ion current chromatogram", + )), + precursor_mz: None, + product_mz: None, + time_sec: tic.iter().map(|&(rt, _)| (rt * 60.0) as f32).collect(), + intensity: tic.iter().map(|&(_, i)| i as f32).collect(), + }) +} + +/// Build the base-peak chromatogram from `(rt_min, base_intensity, base_mz)` +/// triples, or `None` when the file has no scans. `base_mz` is not carried on +/// a chromatogram record, so only the intensity trace is emitted. +fn bpc_record(bpc: &[(f64, f64, f64)]) -> Option { + if bpc.is_empty() { + return None; + } + Some(msc::ChromatogramRecord { + index: 0, + id: "BPC".to_string(), + chromatogram_type: Some(msc::CvTerm::new("MS:1000628", "basepeak chromatogram")), + precursor_mz: None, + product_mz: None, + time_sec: bpc.iter().map(|&(rt, _, _)| (rt * 60.0) as f32).collect(), + intensity: bpc.iter().map(|&(_, bpi, _)| bpi as f32).collect(), + }) +} + +/// Group SRM/flat-peak scans into one chromatogram per transition. +/// +/// `scans` is `(scan_event, rt_min, total_current)` in scan (time) order. +/// Each `scan_event` is one Q1->Q3 transition; its points are the summed +/// product-ion current (`total_current`) already decoded per scan, so no +/// re-summing of peaks is needed. A transition is emitted only when its Q1 +/// precursor is known (`q1_by_event`); the product m/z is filled in only when +/// exactly one Q3 window resolves unambiguously, and left unset otherwise +/// rather than guessed. Events are visited in sorted order for deterministic +/// output. +fn srm_chromatograms( + scans: &[(u16, f64, f64)], + q1_by_event: &std::collections::HashMap, + q3_windows: &std::collections::HashMap>, +) -> Vec { + let mut points_by_event: std::collections::BTreeMap> = + std::collections::BTreeMap::new(); + for &(event, rt_min, current) in scans { + points_by_event + .entry(event) + .or_default() + .push(((rt_min * 60.0) as f32, current as f32)); + } + + let mut out = Vec::new(); + for (event, points) in points_by_event { + let Some(&q1) = q1_by_event.get(&event) else { + continue; + }; + let product_mz = match q3_windows.get(&event) { + Some(windows) if windows.len() == 1 => { + let (lo, hi) = windows[0]; + Some(((lo + hi) / 2.0) as f64) + } + _ => None, + }; + let id = match product_mz { + Some(q3) => format!("SRM Q1={q1:.4} Q3={q3:.4}"), + None => format!("SRM Q1={q1:.4}"), + }; + out.push(msc::ChromatogramRecord { + index: 0, + id, + chromatogram_type: Some(msc::CvTerm::new( + "MS:1000627", + "selected reaction monitoring chromatogram", + )), + precursor_mz: Some(q1), + product_mz, + time_sec: points.iter().map(|&(t, _)| t).collect(), + intensity: points.iter().map(|&(_, i)| i).collect(), + }); + } + out +} + +/// Assemble every chromatogram trace `raw` can supply: always a TIC and BPC +/// (when the file has scans), plus one SRM chromatogram per transition for +/// flat-peak SRM files. Records are indexed sequentially in emission order. +fn build_chromatograms(raw: &RawFileReader) -> Vec { + let mut out = Vec::new(); + out.extend(tic_record(&raw.tic_chromatogram())); + out.extend(bpc_record(&raw.bpc_chromatogram())); + + if raw.flat_peaks && !raw.srm_q1_by_event.is_empty() { + let scans: Vec<(u16, f64, f64)> = raw + .scan_index + .iter() + .map(|e| (e.scan_event, e.start_time, e.total_current)) + .collect(); + out.extend(srm_chromatograms( + &scans, + &raw.srm_q1_by_event, + &raw.srm_q3_windows, + )); + } + + for (i, rec) in out.iter_mut().enumerate() { + rec.index = i; + } + out +} + /// `SpectrumSource` adapter over a Thermo RAW reader. /// /// Use this when you want to feed Thermo data into any @@ -553,6 +689,12 @@ impl<'a, R: Read + Seek> msc::SpectrumSource for OpenTfRawSource<'a, R> { fn spectrum_count_hint(&self) -> Option { Some(self.raw.num_scans as usize) } + + fn iter_chromatograms<'s>( + &'s mut self, + ) -> Box + 's> { + Box::new(build_chromatograms(self.raw).into_iter()) + } } // -- Public mzML entry points (unchanged signatures) -- @@ -605,3 +747,100 @@ where msc::write_indexed_mzml(&mut src, out)?; Ok(()) } + +#[cfg(test)] +mod tests { + use super::*; + use std::collections::HashMap; + + #[test] + fn tic_record_maps_pairs_and_converts_minutes_to_seconds() { + // (rt_min, total_current) + let tic = [(0.0, 100.0), (0.5, 250.0), (1.0, 175.0)]; + let rec = tic_record(&tic).expect("non-empty TIC"); + assert_eq!(rec.id, "TIC"); + let cv = rec.chromatogram_type.as_ref().unwrap(); + assert_eq!(cv.accession, "MS:1000235"); + assert_eq!(cv.name, "total ion current chromatogram"); + assert!(rec.precursor_mz.is_none()); + assert!(rec.product_mz.is_none()); + assert_eq!(rec.time_sec, vec![0.0, 30.0, 60.0]); // rt_min * 60 + assert_eq!(rec.intensity, vec![100.0, 250.0, 175.0]); + } + + #[test] + fn bpc_record_uses_base_intensity_and_drops_base_mz() { + // (rt_min, base_intensity, base_mz) + let bpc = [(0.0, 90.0, 500.0), (0.5, 220.0, 501.0)]; + let rec = bpc_record(&bpc).expect("non-empty BPC"); + assert_eq!(rec.id, "BPC"); + let cv = rec.chromatogram_type.as_ref().unwrap(); + assert_eq!(cv.accession, "MS:1000628"); + assert_eq!(cv.name, "basepeak chromatogram"); + assert_eq!(rec.time_sec, vec![0.0, 30.0]); + assert_eq!(rec.intensity, vec![90.0, 220.0]); // base_mz not carried + } + + #[test] + fn tic_and_bpc_records_are_none_when_no_scans() { + assert!(tic_record(&[]).is_none()); + assert!(bpc_record(&[]).is_none()); + } + + #[test] + fn srm_groups_scans_by_event_and_sets_precursor_and_product() { + // Two transitions interleaved in scan order, one scan_event each. + let scans = [ + (1u16, 0.0f64, 10.0f64), + (2u16, 0.0, 20.0), + (1, 0.5, 12.0), + (2, 0.5, 22.0), + ]; + let mut q1 = HashMap::new(); + q1.insert(1u16, 500.0); + q1.insert(2u16, 600.0); + let mut q3 = HashMap::new(); + q3.insert(1u16, vec![(100.0f32, 101.0f32)]); // single window -> product m/z + q3.insert(2u16, vec![(200.0f32, 201.0f32)]); + + let recs = srm_chromatograms(&scans, &q1, &q3); + assert_eq!(recs.len(), 2); + + // BTreeMap order: event 1 first. + let r1 = &recs[0]; + assert_eq!(r1.precursor_mz, Some(500.0)); + assert_eq!(r1.product_mz, Some(100.5)); // (100 + 101) / 2 + assert_eq!(r1.id, "SRM Q1=500.0000 Q3=100.5000"); + let cv = r1.chromatogram_type.as_ref().unwrap(); + assert_eq!(cv.accession, "MS:1000627"); + assert_eq!(cv.name, "selected reaction monitoring chromatogram"); + assert_eq!(r1.time_sec, vec![0.0, 30.0]); + assert_eq!(r1.intensity, vec![10.0, 12.0]); + + let r2 = &recs[1]; + assert_eq!(r2.precursor_mz, Some(600.0)); + assert_eq!(r2.product_mz, Some(200.5)); + assert_eq!(r2.intensity, vec![20.0, 22.0]); + } + + #[test] + fn srm_skips_events_without_a_known_q1_and_omits_ambiguous_product() { + let scans = [ + (1u16, 0.0f64, 10.0f64), // no Q1 -> skipped + (2u16, 0.0, 20.0), // Q1 known, two Q3 windows -> product omitted + ]; + let mut q1 = HashMap::new(); + q1.insert(2u16, 600.0); + let mut q3 = HashMap::new(); + q3.insert(2u16, vec![(200.0f32, 201.0f32), (300.0f32, 301.0f32)]); + + let recs = srm_chromatograms(&scans, &q1, &q3); + assert_eq!(recs.len(), 1, "event without a Q1 must be skipped"); + assert_eq!(recs[0].precursor_mz, Some(600.0)); + assert!( + recs[0].product_mz.is_none(), + "ambiguous multi-window Q3 must not be guessed" + ); + assert_eq!(recs[0].id, "SRM Q1=600.0000"); + } +} From 39631b97ee7c34349f4263eb738c1fb9696607ee Mon Sep 17 00:00:00 2001 From: Nathan Riley Date: Tue, 21 Jul 2026 14:47:04 +0000 Subject: [PATCH 2/2] fix: disambiguate SRM chromatogram ids on Q1/Q3 collision Two scan_events can legitimately share the same Q1/Q3 (e.g. a scheduled SRM method re-monitoring the same transition across separate retention-time windows), which would otherwise produce two chromatograms with an identical id. mzML's chromatogram id is a plain xsd:string (not xsd:ID), so nothing upstream validates uniqueness, but openmassspec_core's indexed mzML writer keys its offset index by id, so a collision would leave one transition unreachable by id-based lookup. Append the scan_event to disambiguate any repeat. --- CHANGELOG.md | 6 ++++- crates/opentfraw/src/mzml.rs | 43 +++++++++++++++++++++++++++++++++++- 2 files changed, 47 insertions(+), 2 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 8a5ad55..6122b9a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -41,7 +41,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 monitoring chromatogram (`MS:1000627`) per transition, grouping scans by `scan_event` with `precursor_mz` from the Q1 map and `product_mz` filled in only when a single Q3 window resolves unambiguously. No new binary parsing - is involved. (contributed by @Nabejo) + is involved. SRM chromatogram ids are disambiguated with the `scan_event` + when two transitions share the same Q1/Q3 (e.g. a scheduled method + re-monitoring the same pair across separate retention-time windows), since + the indexed mzML writer keys its offset index by id. (contributed by + @Nabejo) - Unit tests for the pure-decode functions in `bytes.rs`, `header.rs`, `audit_tag.rs`, `generic_data.rs`, and `error_log.rs`, using hand-built byte fixtures (valid, truncated, and out-of-range cases). Also added diff --git a/crates/opentfraw/src/mzml.rs b/crates/opentfraw/src/mzml.rs index 2452e40..b4e6948 100644 --- a/crates/opentfraw/src/mzml.rs +++ b/crates/opentfraw/src/mzml.rs @@ -552,6 +552,16 @@ fn bpc_record(bpc: &[(f64, f64, f64)]) -> Option { /// exactly one Q3 window resolves unambiguously, and left unset otherwise /// rather than guessed. Events are visited in sorted order for deterministic /// output. +/// +/// Two transitions can legitimately share the same Q1/Q3 (e.g. a scheduled +/// method re-monitoring the same pair across separate retention-time +/// windows, each its own `scan_event`), which would otherwise collide on the +/// same `Q1={..} Q3={..}` id - `mzML`'s chromatogram `id` is a plain +/// `xsd:string`, not `xsd:ID`, so nothing upstream would catch a duplicate, +/// but the indexed-mzML writer keys its offset index by `id` +/// (`openmassspec_core::mzml::write_indexed_mzml`), so a collision there +/// would make one of the two transitions unreachable by id-based lookup. +/// Disambiguate any repeat by appending its `scan_event`. fn srm_chromatograms( scans: &[(u16, f64, f64)], q1_by_event: &std::collections::HashMap, @@ -567,6 +577,7 @@ fn srm_chromatograms( } let mut out = Vec::new(); + let mut seen_ids: std::collections::HashSet = std::collections::HashSet::new(); for (event, points) in points_by_event { let Some(&q1) = q1_by_event.get(&event) else { continue; @@ -578,10 +589,15 @@ fn srm_chromatograms( } _ => None, }; - let id = match product_mz { + let base_id = match product_mz { Some(q3) => format!("SRM Q1={q1:.4} Q3={q3:.4}"), None => format!("SRM Q1={q1:.4}"), }; + let id = if seen_ids.insert(base_id.clone()) { + base_id + } else { + format!("{base_id} event={event}") + }; out.push(msc::ChromatogramRecord { index: 0, id, @@ -843,4 +859,29 @@ mod tests { ); assert_eq!(recs[0].id, "SRM Q1=600.0000"); } + + #[test] + fn srm_disambiguates_ids_when_two_events_share_the_same_q1_and_q3() { + // Scheduled method re-monitoring the same transition in two separate + // time windows: same Q1/Q3, different scan_event. + let scans = [(1u16, 0.0f64, 10.0f64), (2u16, 5.0, 15.0)]; + let mut q1 = HashMap::new(); + q1.insert(1u16, 500.0); + q1.insert(2u16, 500.0); + let mut q3 = HashMap::new(); + q3.insert(1u16, vec![(100.0f32, 101.0f32)]); + q3.insert(2u16, vec![(100.0f32, 101.0f32)]); + + let recs = srm_chromatograms(&scans, &q1, &q3); + assert_eq!(recs.len(), 2); + assert_eq!(recs[0].id, "SRM Q1=500.0000 Q3=100.5000"); + assert_eq!( + recs[1].id, "SRM Q1=500.0000 Q3=100.5000 event=2", + "second event with an identical Q1/Q3 must get a disambiguated id" + ); + assert_ne!( + recs[0].id, recs[1].id, + "ids feed the indexed-mzML offset index and must stay unique" + ); + } }