Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
18 changes: 18 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,24 @@ This project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.htm

## [Unreleased]

### Added

- mzML export now emits a `<chromatogramList>` via
`SpectrumSource::iter_chromatograms` on `TdfSource` / `OwnedTdfSource`,
which previously used the trait's default empty implementation. Two
whole-run traces are produced, one point per frame in acquisition order,
built directly from columns `Reader::frames` already SELECTs (no extra
SQL, no peak decode): a TIC ("total ion current chromatogram",
`MS:1000235`) from `Frames.SummedIntensities`, and a basepeak
chromatogram ("basepeak chromatogram", `MS:1000628`) from
`Frames.MaxIntensity`. `MaxIntensity` is a new additive column on the
`Frames` SELECT and a new `Frame::max_intensity` field (also exposed on
the Python `Frame`). A trace whose source column is absent on every frame
is omitted rather than emitted empty. SRM/PRM transition chromatograms
are intentionally left out: they need genuine per-target cross-frame
aggregation of the decoded peak stream, not just column wiring. Closes
#25. Contributed by @Nabejo.

### Fixed

- `openmassspec-core`'s declared minimum version in `Cargo.toml` was still
Expand Down
5 changes: 5 additions & 0 deletions crates/opentimstdf-py/src/lib.rs
Original file line number Diff line number Diff line change
Expand Up @@ -80,6 +80,8 @@ struct Frame {
accumulation_time: Option<f64>,
#[pyo3(get)]
summed_intensities: Option<u64>,
#[pyo3(get)]
max_intensity: Option<u64>,
}

#[pymethods]
Expand Down Expand Up @@ -117,6 +119,7 @@ impl From<RsFrame> for Frame {
mz_calibration_id: f.mz_calibration_id,
accumulation_time: f.accumulation_time,
summed_intensities: f.summed_intensities,
max_intensity: f.max_intensity,
}
}
}
Expand Down Expand Up @@ -535,6 +538,7 @@ impl Reader {
mz_calibration_id: frame.mz_calibration_id,
accumulation_time: frame.accumulation_time,
summed_intensities: frame.summed_intensities,
max_intensity: frame.max_intensity,
};
Ok(self
.locked_inner()?
Expand Down Expand Up @@ -565,6 +569,7 @@ impl Reader {
mz_calibration_id: frame.mz_calibration_id,
accumulation_time: frame.accumulation_time,
summed_intensities: frame.summed_intensities,
max_intensity: frame.max_intensity,
};
let peaks = guard.decode_peaks(&rs_frame).map_err(to_py_err)?;
let n = peaks.len();
Expand Down
3 changes: 2 additions & 1 deletion crates/opentimstdf/src/codec.rs
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ use crate::types::{Frame, Peak};
///
/// Column order must match the SELECT used in `Reader::frame()` and `Reader::frames()`:
/// Id, Time, NumScans, NumPeaks, TimsId, ScanMode, MsMsType,
/// MzCalibration, AccumulationTime, SummedIntensities
/// MzCalibration, AccumulationTime, SummedIntensities, MaxIntensity
pub(crate) fn frame_from_row(row: &rusqlite::Row<'_>) -> rusqlite::Result<Frame> {
Ok(Frame {
id: row.get(0)?,
Expand All @@ -17,6 +17,7 @@ pub(crate) fn frame_from_row(row: &rusqlite::Row<'_>) -> rusqlite::Result<Frame>
mz_calibration_id: row.get::<_, Option<i64>>(7)?.unwrap_or(1) as u32,
accumulation_time: row.get(8)?,
summed_intensities: row.get::<_, Option<i64>>(9)?.map(|v| v as u64),
max_intensity: row.get::<_, Option<i64>>(10)?.map(|v| v as u64),
})
}

Expand Down
184 changes: 184 additions & 0 deletions crates/opentimstdf/src/mzml.rs
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,11 @@
//! * [`write_mzml`] / [`write_indexed_mzml`] - convenience entry points
//! that wrap the canonical writer in `openmassspec-core`.
//!
//! Chromatograms: [`iter_chromatograms`](msc::SpectrumSource::iter_chromatograms)
//! emits whole-run TIC and basepeak (BPC) traces built directly from the
//! per-frame `Frames.SummedIntensities` / `Frames.MaxIntensity` columns (one
//! point per frame). See [`chromatogram_records_for`].
//!
//! Frame -> spectrum projection:
//!
//! * **MS1 frames** (`msms_type == 0`): all mobility scans pooled into one
Expand Down Expand Up @@ -569,13 +574,92 @@ fn run_metadata_for(meta: &Metadata, bundle_name: &str) -> msc::RunMetadata {
}
}

/// Build the whole-run chromatogram traces from the already-decoded per-frame
/// `Frames` columns carried on `frames`.
///
/// Two traces are emitted, each a time series with one point per frame (in
/// acquisition order) built purely from columns `Reader::frames` already
/// SELECTs - no extra SQL, no peak decode:
///
/// * **TIC** - "total ion current chromatogram" (`MS:1000235`) from
/// [`Frame::summed_intensities`] (`Frames.SummedIntensities`, the
/// accumulation-normalized per-frame intensity sum documented in
/// `docs/docs/format/01-tdf-sqlite-schema.md` and cross-checked against the
/// decoded peak sum in `tests/roundtrip.rs`).
/// * **BPC** - "basepeak chromatogram" (`MS:1000628`) from
/// [`Frame::max_intensity`] (`Frames.MaxIntensity`, the most intense peak in
/// the frame).
///
/// `Frame::time` is retention time in seconds, which is what
/// [`msc::ChromatogramRecord::time_sec`] wants (the canonical writer converts
/// to minutes for the mzML `time array` itself). A frame whose source column
/// is absent (`None`) contributes no point; a trace with no points at all is
/// omitted rather than emitted empty, so a bundle missing a column yields no
/// misleading zero-length chromatogram - matching the writer's contract of
/// only emitting `<chromatogramList>` when the source yields something.
///
/// SRM/PRM transition chromatograms are intentionally not produced here: that
/// needs genuine per-target cross-frame aggregation of the decoded peak stream
/// (see Sigilweaver/OpenTimsTDF#25), not just wiring of existing columns.
fn chromatogram_records_for(frames: &[Frame]) -> Vec<msc::ChromatogramRecord> {
fn trace(frames: &[Frame], value_of: impl Fn(&Frame) -> Option<u64>) -> (Vec<f32>, Vec<f32>) {
let mut time_sec = Vec::new();
let mut intensity = Vec::new();
for f in frames {
if let Some(v) = value_of(f) {
time_sec.push(f.time as f32);
intensity.push(v as f32);
}
}
(time_sec, intensity)
}

let mut out: Vec<msc::ChromatogramRecord> = Vec::new();

let (tic_time, tic_int) = trace(frames, |f| f.summed_intensities);
if !tic_time.is_empty() {
out.push(msc::ChromatogramRecord {
index: out.len(),
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_time,
intensity: tic_int,
});
}

let (bpc_time, bpc_int) = trace(frames, |f| f.max_intensity);
if !bpc_time.is_empty() {
out.push(msc::ChromatogramRecord {
index: out.len(),
id: "BPC".to_string(),
chromatogram_type: Some(msc::CvTerm::new("MS:1000628", "basepeak chromatogram")),
precursor_mz: None,
product_mz: None,
time_sec: bpc_time,
intensity: bpc_int,
});
}

out
}

impl<'a> msc::SpectrumSource for TdfSource<'a> {
fn run_metadata(&self) -> msc::RunMetadata {
run_metadata_for(&self.metadata, &self.bundle_name)
}
fn iter_spectra<'s>(&'s mut self) -> Box<dyn Iterator<Item = msc::SpectrumRecord> + 's> {
Box::new(frame_iter(self.reader, &self.frames, &self.calibration))
}
fn iter_chromatograms<'s>(
&'s mut self,
) -> Box<dyn Iterator<Item = msc::ChromatogramRecord> + 's> {
Box::new(chromatogram_records_for(&self.frames).into_iter())
}
}

impl msc::SpectrumSource for OwnedTdfSource {
Expand All @@ -585,6 +669,11 @@ impl msc::SpectrumSource for OwnedTdfSource {
fn iter_spectra<'s>(&'s mut self) -> Box<dyn Iterator<Item = msc::SpectrumRecord> + 's> {
Box::new(frame_iter(&self.reader, &self.frames, &self.calibration))
}
fn iter_chromatograms<'s>(
&'s mut self,
) -> Box<dyn Iterator<Item = msc::ChromatogramRecord> + 's> {
Box::new(chromatogram_records_for(&self.frames).into_iter())
}
}

/// Write a `.d/` bundle to mzML in one call.
Expand Down Expand Up @@ -630,6 +719,7 @@ mod tests {
mz_calibration_id: 1,
accumulation_time: None,
summed_intensities: None,
max_intensity: None,
}
}

Expand Down Expand Up @@ -878,4 +968,98 @@ mod tests {
assert_eq!(cv.name, term_name, "wrong CV name for {name:?}");
}
}

fn frame_with_intensities(
id: u32,
time: f64,
summed_intensities: Option<u64>,
max_intensity: Option<u64>,
) -> Frame {
Frame {
id,
time,
num_scans: 10,
num_peaks: 3,
tims_id: 0,
scan_mode: 0,
msms_type: 0,
mz_calibration_id: 1,
accumulation_time: None,
summed_intensities,
max_intensity,
}
}

#[test]
fn chromatogram_records_for_builds_tic_and_bpc_from_frame_columns() {
let frames = [
frame_with_intensities(1, 1.5, Some(100), Some(40)),
frame_with_intensities(2, 3.0, Some(250), Some(90)),
];
let recs = chromatogram_records_for(&frames);
assert_eq!(recs.len(), 2, "expected a TIC and a BPC record");

let tic = &recs[0];
assert_eq!(tic.index, 0);
assert_eq!(tic.id, "TIC");
let tic_cv = tic.chromatogram_type.as_ref().expect("TIC CV term");
assert_eq!(tic_cv.accession, "MS:1000235");
assert_eq!(tic_cv.name, "total ion current chromatogram");
assert!(tic.precursor_mz.is_none());
assert!(tic.product_mz.is_none());
// Retention time carried through in seconds (not minutes).
assert_eq!(tic.time_sec, vec![1.5, 3.0]);
assert_eq!(tic.intensity, vec![100.0, 250.0]);

let bpc = &recs[1];
assert_eq!(bpc.index, 1);
assert_eq!(bpc.id, "BPC");
let bpc_cv = bpc.chromatogram_type.as_ref().expect("BPC CV term");
assert_eq!(bpc_cv.accession, "MS:1000628");
assert_eq!(bpc_cv.name, "basepeak chromatogram");
assert_eq!(bpc.time_sec, vec![1.5, 3.0]);
assert_eq!(bpc.intensity, vec![40.0, 90.0]);
}

#[test]
fn chromatogram_records_for_skips_frames_missing_a_column_keeping_arrays_parallel() {
// Middle frame has no SummedIntensities and no MaxIntensity: it must
// contribute no point to either trace, and the two arrays stay
// parallel (one time per intensity).
let frames = [
frame_with_intensities(1, 1.0, Some(10), Some(5)),
frame_with_intensities(2, 2.0, None, None),
frame_with_intensities(3, 3.0, Some(30), Some(15)),
];
let recs = chromatogram_records_for(&frames);
assert_eq!(recs.len(), 2);
assert_eq!(recs[0].time_sec, vec![1.0, 3.0]);
assert_eq!(recs[0].intensity, vec![10.0, 30.0]);
assert_eq!(recs[1].time_sec, vec![1.0, 3.0]);
assert_eq!(recs[1].intensity, vec![5.0, 15.0]);
}

#[test]
fn chromatogram_records_for_omits_a_trace_with_no_points() {
// SummedIntensities present on every frame but MaxIntensity absent
// (e.g. an older schema without the column): only the TIC is emitted,
// and it keeps index 0 - no empty BPC record is produced.
let frames = [
frame_with_intensities(1, 1.0, Some(10), None),
frame_with_intensities(2, 2.0, Some(20), None),
];
let recs = chromatogram_records_for(&frames);
assert_eq!(recs.len(), 1);
assert_eq!(recs[0].id, "TIC");
assert_eq!(recs[0].index, 0);
}

#[test]
fn chromatogram_records_for_empty_when_no_columns_populated() {
let frames = [
frame_with_intensities(1, 1.0, None, None),
frame_with_intensities(2, 2.0, None, None),
];
assert!(chromatogram_records_for(&frames).is_empty());
}
}
4 changes: 2 additions & 2 deletions crates/opentimstdf/src/reader.rs
Original file line number Diff line number Diff line change
Expand Up @@ -236,7 +236,7 @@ impl Reader {
let conn = self.conn.lock().map_err(|_| Error::LockPoisoned)?;
let frame = conn.query_row(
"SELECT Id, Time, NumScans, NumPeaks, TimsId, ScanMode, MsMsType,
MzCalibration, AccumulationTime, SummedIntensities
MzCalibration, AccumulationTime, SummedIntensities, MaxIntensity
FROM Frames WHERE Id = ?1",
[frame_id],
frame_from_row,
Expand All @@ -249,7 +249,7 @@ impl Reader {
let conn = self.conn.lock().map_err(|_| Error::LockPoisoned)?;
let mut stmt = conn.prepare(
"SELECT Id, Time, NumScans, NumPeaks, TimsId, ScanMode, MsMsType,
MzCalibration, AccumulationTime, SummedIntensities
MzCalibration, AccumulationTime, SummedIntensities, MaxIntensity
FROM Frames ORDER BY Id ASC",
)?;
let rows = stmt.query_map([], frame_from_row)?;
Expand Down
3 changes: 3 additions & 0 deletions crates/opentimstdf/src/types.rs
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,9 @@ pub struct Frame {
/// Raw sum of decoded intensities divided by (AccumulationTime/100).
/// See SPEC §2.2 for the normalization formula.
pub summed_intensities: Option<u64>,
/// Most intense peak in the frame (`Frames.MaxIntensity`). `None` when the
/// column is absent. Used to build a basepeak chromatogram.
pub max_intensity: Option<u64>,
}

/// One decoded peak.
Expand Down
3 changes: 2 additions & 1 deletion docs/docs/format/01-tdf-sqlite-schema.md
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,8 @@ Columns that matter for the binary stream:
| `Time` | REAL | Retention time (seconds from start of acquisition). |
| `AccumulationTime` | REAL | Ion accumulation time in milliseconds (sv >= 3.5; absent in older schemas). |
| `RampTime` | REAL | TIMS ramp (elution) time in milliseconds. Equal to `AccumulationTime` in all observed bundles. |
| `SummedIntensities` | INT | Sum of all peak intensities in the frame, **normalised to 100 ms accumulation**: `SummedIntensities = sum(raw_intensity) * 100.0 / AccumulationTime_ms`. For bundles where `AccumulationTime ~ 100` the factor is ~ 1.0. Verified on PXD022216 (AT=108.46 ms): ratio decoded_sum / SummedIntensities = 1.08460 +/- 0.00001 across all 57,886 frames. |
| `SummedIntensities` | INT | Sum of all peak intensities in the frame, **normalised to 100 ms accumulation**: `SummedIntensities = sum(raw_intensity) * 100.0 / AccumulationTime_ms`. For bundles where `AccumulationTime ~ 100` the factor is ~ 1.0. Verified on PXD022216 (AT=108.46 ms): ratio decoded_sum / SummedIntensities = 1.08460 +/- 0.00001 across all 57,886 frames. Consumed as the per-frame point of the whole-run TIC chromatogram in mzML export. |
| `MaxIntensity` | INT | Most intense single peak in the frame. Consumed as the per-frame point of the whole-run basepeak chromatogram (BPC) in mzML export. |
| `MzCalibration` | INT | Foreign key into `MzCalibration` table. |
| `TimsCalibration` | INT | Foreign key into `TimsCalibration` table. |

Expand Down
1 change: 1 addition & 0 deletions docs/docs/reference/python-api.md
Original file line number Diff line number Diff line change
Expand Up @@ -76,6 +76,7 @@ stream.
| `mz_calibration_id` | `int` | Calibration id used for this frame |
| `accumulation_time` | `float \| None` | TIMS accumulation time (ms) |
| `summed_intensities` | `int \| None` | Sum of all peak intensities in the frame |
| `max_intensity` | `int \| None` | Most intense peak in the frame (`Frames.MaxIntensity`) |
| `polarity` | `str` | `"positive"` or `"negative"`, derived from `mz_calibration_id` |

### Metadata
Expand Down
Loading