use crate::error::{Error, Result};
use crate::reader::BinaryReader;
use std::io::{Read, Seek, SeekFrom};
#[derive(Debug, Clone)]
pub struct Peak {
pub mz: f64,
pub abundance: f32,
}
#[derive(Debug, Clone, Copy)]
pub struct NoiseNode {
pub mz: f32,
pub noise: f32,
pub baseline: f32,
}
#[derive(Debug)]
pub struct ProfileChunk {
pub first_bin: u32,
pub signal: Vec<f32>,
pub fudge: Option<f32>,
}
#[derive(Debug)]
pub struct Profile {
pub first_value: f64,
pub step: f64,
pub peak_count: u32,
pub nbins: u32,
pub chunks: Vec<ProfileChunk>,
}
#[derive(Debug)]
pub struct PacketHeader {
pub profile_size: u32,
pub peak_list_size: u32,
pub layout: u32,
pub descriptor_list_size: u32,
pub unknown_stream_size: u32,
pub triplet_stream_size: u32,
pub low_mz: f32,
pub high_mz: f32,
}
#[derive(Debug)]
pub struct ScanDataPacket {
pub header: PacketHeader,
pub profile: Option<Profile>,
pub peaks: Vec<Peak>,
pub resolutions: Vec<f32>,
pub noise_nodes: Vec<NoiseNode>,
}
impl ScanDataPacket {
pub(crate) fn read<R: Read + Seek>(r: &mut BinaryReader<R>) -> Result<Self> {
let header = PacketHeader::read(r)?;
let profile = if header.profile_size > 0 {
Some(Profile::read(r, header.layout)?)
} else {
None
};
let peaks = Self::read_peaks(r, &header)?;
let (resolutions, noise_nodes) = Self::read_labels(r, &header, peaks.len())?;
Ok(Self {
header,
profile,
peaks,
resolutions,
noise_nodes,
})
}
pub(crate) fn read_skip_profile<R: Read + Seek>(r: &mut BinaryReader<R>) -> Result<Self> {
let header = PacketHeader::read(r)?;
if header.profile_size > 0 {
r.skip((header.profile_size as usize) * 4)?;
}
let peaks = Self::read_peaks(r, &header)?;
let (resolutions, noise_nodes) = Self::read_labels(r, &header, peaks.len())?;
Ok(Self {
header,
profile: None,
peaks,
resolutions,
noise_nodes,
})
}
pub(crate) fn read_peaks_only<R: Read + Seek>(r: &mut BinaryReader<R>) -> Result<Vec<Peak>> {
let header = PacketHeader::read(r)?;
if header.profile_size > 0 {
r.skip((header.profile_size as usize) * 4)?;
}
Self::read_peaks(r, &header)
}
fn read_peaks<R: Read + Seek>(
r: &mut BinaryReader<R>,
header: &PacketHeader,
) -> Result<Vec<Peak>> {
let wide_mz = header.layout & 0x10000 != 0;
if header.peak_list_size == 0 {
return Ok(Vec::new());
}
let count = r.read_u32()?;
let item_size: u64 = if wide_mz { 12 } else { 8 };
r.check_count(count as u64, item_size)?;
let mut peaks = Vec::with_capacity(count as usize);
for _ in 0..count {
let mz = if wide_mz {
r.read_f64()?
} else {
r.read_f32()? as f64
};
let abundance = r.read_f32()?;
peaks.push(Peak { mz, abundance });
}
Ok(peaks)
}
fn read_labels<R: Read + Seek>(
r: &mut BinaryReader<R>,
header: &PacketHeader,
n_peaks: usize,
) -> Result<(Vec<f32>, Vec<NoiseNode>)> {
if header.descriptor_list_size > 0 {
r.skip(header.descriptor_list_size as usize * 4)?;
}
let resolutions = if n_peaks > 0 && header.unknown_stream_size as usize == n_peaks + 1 {
let _count = r.read_u32()?;
r.check_count(n_peaks as u64, 4)?;
let mut res = Vec::with_capacity(n_peaks);
for _ in 0..n_peaks {
res.push(r.read_f32()?);
}
res
} else {
if header.unknown_stream_size > 0 {
r.skip(header.unknown_stream_size as usize * 4)?;
}
Vec::new()
};
let node_count = header.triplet_stream_size / 3;
r.check_count(node_count as u64, 12)?;
let mut noise_nodes = Vec::with_capacity(node_count as usize);
for _ in 0..node_count {
let mz = r.read_f32()?;
let noise = r.read_f32()?;
let baseline = r.read_f32()?;
noise_nodes.push(NoiseNode {
mz,
noise,
baseline,
});
}
let consumed = node_count * 3;
if header.triplet_stream_size > consumed {
r.skip((header.triplet_stream_size - consumed) as usize * 4)?;
}
Ok((resolutions, noise_nodes))
}
pub fn noise_at(&self, mz: f64) -> Option<(f32, f32)> {
let nodes = &self.noise_nodes;
if nodes.is_empty() {
return None;
}
let x = mz as f32;
if x <= nodes[0].mz {
return Some((nodes[0].noise, nodes[0].baseline));
}
for w in nodes.windows(2) {
let (lo, hi) = (&w[0], &w[1]);
if x <= hi.mz {
let span = hi.mz - lo.mz;
let f = if span > 0.0 { (x - lo.mz) / span } else { 0.0 };
let noise = lo.noise + f * (hi.noise - lo.noise);
let baseline = lo.baseline + f * (hi.baseline - lo.baseline);
return Some((noise, baseline));
}
}
let last = nodes[nodes.len() - 1];
Some((last.noise, last.baseline))
}
}
impl PacketHeader {
fn read<R: Read + Seek>(r: &mut BinaryReader<R>) -> Result<Self> {
let _unk1 = r.read_u32()?;
let profile_size = r.read_u32()?;
let peak_list_size = r.read_u32()?;
let layout = r.read_u32()?;
let descriptor_list_size = r.read_u32()?;
let unknown_stream_size = r.read_u32()?;
let triplet_stream_size = r.read_u32()?;
let _unk2 = r.read_u32()?;
let low_mz = r.read_f32()?;
let high_mz = r.read_f32()?;
Ok(Self {
profile_size,
peak_list_size,
layout,
descriptor_list_size,
unknown_stream_size,
triplet_stream_size,
low_mz,
high_mz,
})
}
}
impl Profile {
fn read<R: Read + Seek>(r: &mut BinaryReader<R>, layout: u32) -> Result<Self> {
let first_value = r.read_f64()?;
let step = r.read_f64()?;
let peak_count = r.read_u32()?;
let nbins = r.read_u32()?;
let has_fudge = layout & 0xFF != 0;
r.check_count(peak_count as u64, 8)?;
let mut chunks = Vec::with_capacity(peak_count as usize);
for _ in 0..peak_count {
let first_bin = r.read_u32()?;
let chunk_nbins = r.read_u32()?;
let fudge = if has_fudge { Some(r.read_f32()?) } else { None };
r.check_count(chunk_nbins as u64, 4)?;
let mut signal = Vec::with_capacity(chunk_nbins as usize);
for _ in 0..chunk_nbins {
signal.push(r.read_f32()?);
}
chunks.push(ProfileChunk {
first_bin,
signal,
fudge,
});
}
Ok(Self {
first_value,
step,
peak_count,
nbins,
chunks,
})
}
}
impl Profile {
pub fn to_mz_intensity(&self, coefficients: &[f64]) -> Vec<(f64, f64)> {
let cap: usize = self.chunks.iter().map(|c| c.signal.len()).sum();
let mut result = Vec::with_capacity(cap);
for chunk in &self.chunks {
for (i, &intensity) in chunk.signal.iter().enumerate() {
let bin_global = chunk.first_bin as f64 + i as f64;
let freq = self.first_value + bin_global * self.step;
let freq_adj = if let Some(fudge) = chunk.fudge {
freq + fudge as f64
} else {
freq
};
let mz = freq_to_mz(freq_adj, coefficients);
result.push((mz, intensity as f64));
}
}
result
}
}
pub fn freq_to_mz(freq: f64, coefficients: &[f64]) -> f64 {
if freq == 0.0 {
return 0.0;
}
match coefficients.len() {
0 => freq, 4 => {
let (a, b, c) = (coefficients[1], coefficients[2], coefficients[3]);
a + b / freq + c / (freq * freq)
}
5 => {
let (a, b, c) = (coefficients[2], coefficients[3], coefficients[4]);
let f2 = freq * freq;
a + b / f2 + c / (f2 * f2)
}
7 => {
let (a, b, c) = (coefficients[2], coefficients[3], coefficients[4]);
let f2 = freq * freq;
a + b / f2 + c / (f2 * f2)
}
_ => freq,
}
}
pub fn read_flat_peaks<R: Read + Seek>(
source: &mut R,
data_addr: u64,
cum_end: u64,
data_size: u32,
) -> Result<Vec<Peak>> {
if data_size <= 1 {
return Ok(Vec::new());
}
for subtract in [1u32, 2] {
if data_size <= subtract {
continue;
}
let peak_count = (data_size - subtract) as usize;
let peak_section_bytes = peak_count as u64 * 9;
if peak_section_bytes > cum_end {
continue;
}
let peaks_start = data_addr
.saturating_add(cum_end)
.saturating_sub(peak_section_bytes);
source.seek(SeekFrom::Start(peaks_start))?;
let mut r = BinaryReader::new(&mut *source);
r.check_count(peak_count as u64, 8)?;
let mut peaks = Vec::with_capacity(peak_count);
for _ in 0..peak_count {
let mz = r.read_f32()? as f64;
let abundance = r.read_f32()?;
peaks.push(Peak { mz, abundance });
}
let looks_valid = if let Some(first) = peaks.first() {
first.mz == 0.0 || (first.mz > 10.0 && first.mz < 10_000.0)
} else {
true
};
if looks_valid {
return Ok(peaks);
}
}
Err(Error::UnexpectedEof {
offset: data_addr.saturating_add(cum_end),
needed: 0,
})
}
pub fn read_scan_srm_v66<R: Read + Seek>(
source: &mut R,
data_addr: u64,
start_offset: u64,
_record_size: u32,
) -> Result<Vec<Peak>> {
let abs_start = data_addr.saturating_add(start_offset);
source.seek(SeekFrom::Start(abs_start))?;
let mut r = BinaryReader::new(source);
let n_peaks = r.read_u32()? as usize;
if n_peaks == 0 {
return Ok(Vec::new());
}
r.skip(28)?;
r.skip(n_peaks * 8)?;
r.check_count(n_peaks as u64, 12)?;
let mut peaks = Vec::with_capacity(n_peaks);
for _ in 0..n_peaks {
let _channel = r.read_u32()?;
let mz = r.read_f32()? as f64;
let abundance = r.read_f32()?;
peaks.push(Peak { mz, abundance });
}
Ok(peaks)
}
pub fn read_scan_srm_v66_windows<R: Read + Seek>(
source: &mut R,
data_addr: u64,
start_offset: u64,
) -> Result<Vec<(f32, f32)>> {
let abs_start = data_addr.saturating_add(start_offset);
source.seek(SeekFrom::Start(abs_start))?;
let mut r = BinaryReader::new(source);
let n_peaks = r.read_u32()? as usize;
if n_peaks == 0 {
return Ok(Vec::new());
}
r.skip(28)?;
r.check_count(n_peaks as u64, 8)?;
let mut windows = Vec::with_capacity(n_peaks);
for _ in 0..n_peaks {
let lo = r.read_f32()?;
let hi = r.read_f32()?;
windows.push((lo, hi));
}
Ok(windows)
}
pub fn search_v63_transition(data: &[u8], q3_center_target: f64) -> Option<(f64, f64, f64)> {
let end = data.len().saturating_sub(32);
for j in 8..end {
if j + 8 > data.len() {
break;
}
let v = f64::from_le_bytes(data[j..j + 8].try_into().ok()?);
if (v - q3_center_target).abs() > 0.002 {
continue;
}
let q1 = f64::from_le_bytes(data[j - 8..j].try_into().ok()?);
if !q1.is_finite() || !(50.0..=3000.0).contains(&q1) {
continue;
}
if j + 16 > data.len() {
continue;
}
let q3w = f64::from_le_bytes(data[j + 8..j + 16].try_into().ok()?);
if !q3w.is_finite() || !(0.01..=10.0).contains(&q3w) {
continue;
}
if j + 32 > data.len() {
continue;
}
let ce = f64::from_le_bytes(data[j + 24..j + 32].try_into().ok()?);
if !ce.is_finite() || !(0.1..=300.0).contains(&ce) {
continue;
}
return Some((q1, q3w, ce));
}
None
}
#[cfg(test)]
mod tests {
use super::*;
use std::io::Cursor;
fn header_with(
descriptor_list_size: u32,
unknown_stream_size: u32,
triplet_stream_size: u32,
) -> PacketHeader {
PacketHeader {
profile_size: 0,
peak_list_size: 0,
layout: 0,
descriptor_list_size,
unknown_stream_size,
triplet_stream_size,
low_mz: 0.0,
high_mz: 0.0,
}
}
#[test]
fn huge_descriptor_list_size_does_not_overflow() {
let header = header_with(u32::MAX / 4 + 1, 0, 0);
let mut r = BinaryReader::new(Cursor::new(vec![0u8; 8]));
let _ = ScanDataPacket::read_labels(&mut r, &header, 0);
}
#[test]
fn huge_unknown_stream_size_does_not_overflow() {
let header = header_with(0, u32::MAX / 4 + 1, 0);
let mut r = BinaryReader::new(Cursor::new(vec![0u8; 8]));
let _ = ScanDataPacket::read_labels(&mut r, &header, 0);
}
#[test]
fn read_labels_with_no_streams_returns_empty() {
let header = header_with(0, 0, 0);
let mut r = BinaryReader::new(Cursor::new(Vec::new()));
let (resolutions, noise_nodes) = ScanDataPacket::read_labels(&mut r, &header, 0).unwrap();
assert!(resolutions.is_empty());
assert!(noise_nodes.is_empty());
}
#[test]
fn read_labels_decodes_resolutions_when_size_matches() {
let header = header_with(0, 3, 0);
let mut bytes = 2u32.to_le_bytes().to_vec(); bytes.extend_from_slice(&1.5f32.to_le_bytes());
bytes.extend_from_slice(&2.5f32.to_le_bytes());
let mut r = BinaryReader::new(Cursor::new(bytes));
let (resolutions, _) = ScanDataPacket::read_labels(&mut r, &header, 2).unwrap();
assert_eq!(resolutions, vec![1.5, 2.5]);
}
#[test]
fn read_labels_decodes_noise_nodes() {
let header = header_with(0, 0, 3); let mut bytes = 100.0f32.to_le_bytes().to_vec();
bytes.extend_from_slice(&5.0f32.to_le_bytes());
bytes.extend_from_slice(&1.0f32.to_le_bytes());
let mut r = BinaryReader::new(Cursor::new(bytes));
let (_, noise_nodes) = ScanDataPacket::read_labels(&mut r, &header, 0).unwrap();
assert_eq!(noise_nodes.len(), 1);
assert_eq!(noise_nodes[0].mz, 100.0);
assert_eq!(noise_nodes[0].noise, 5.0);
assert_eq!(noise_nodes[0].baseline, 1.0);
}
#[test]
fn read_peaks_rejects_implausible_count_before_allocating() {
let mut header = header_with(0, 0, 0);
header.peak_list_size = 1;
let mut bytes = u32::MAX.to_le_bytes().to_vec(); bytes.extend_from_slice(&[0u8; 4]); let mut r = BinaryReader::new(Cursor::new(bytes));
let err = ScanDataPacket::read_peaks(&mut r, &header).unwrap_err();
assert!(matches!(err, Error::AllocationTooLarge { .. }));
}
#[test]
fn read_peaks_empty_when_peak_list_size_zero() {
let header = header_with(0, 0, 0);
let mut r = BinaryReader::new(Cursor::new(Vec::new()));
assert_eq!(
ScanDataPacket::read_peaks(&mut r, &header).unwrap().len(),
0
);
}
#[test]
fn read_flat_peaks_decodes_valid_record() {
let mut bytes = vec![0u8; 100]; let peak_start = bytes.len();
bytes.extend_from_slice(&500.0f32.to_le_bytes());
bytes.extend_from_slice(&10.0f32.to_le_bytes());
bytes.push(0); let cum_end = (bytes.len()) as u64;
let mut cursor = Cursor::new(bytes);
let peaks = read_flat_peaks(&mut cursor, 0, cum_end, 2).unwrap();
assert_eq!(peaks.len(), 1);
assert_eq!(peaks[0].mz, 500.0);
assert_eq!(peaks[0].abundance, 10.0);
let _ = peak_start; }
#[test]
fn read_flat_peaks_rejects_implausible_peak_count() {
let bytes = vec![0u8; 16];
let mut cursor = Cursor::new(bytes);
let err = read_flat_peaks(&mut cursor, 0, 16, u32::MAX).unwrap_err();
match err {
Error::AllocationTooLarge { .. } | Error::UnexpectedEof { .. } => {}
other => panic!("unexpected error variant: {other:?}"),
}
}
#[test]
fn read_scan_srm_v66_decodes_valid_record() {
let mut bytes = 1u32.to_le_bytes().to_vec();
bytes.extend_from_slice(&[0u8; 28]);
bytes.extend_from_slice(&100.0f32.to_le_bytes()); bytes.extend_from_slice(&110.0f32.to_le_bytes()); bytes.extend_from_slice(&7u32.to_le_bytes()); bytes.extend_from_slice(&105.0f32.to_le_bytes()); bytes.extend_from_slice(&42.0f32.to_le_bytes()); let mut cursor = Cursor::new(bytes);
let peaks = read_scan_srm_v66(&mut cursor, 0, 0, 0).unwrap();
assert_eq!(peaks.len(), 1);
assert_eq!(peaks[0].mz, 105.0);
assert_eq!(peaks[0].abundance, 42.0);
}
#[test]
fn read_scan_srm_v66_zero_peaks_is_empty() {
let bytes = 0u32.to_le_bytes().to_vec();
let mut cursor = Cursor::new(bytes);
assert!(read_scan_srm_v66(&mut cursor, 0, 0, 0).unwrap().is_empty());
}
#[test]
fn read_scan_srm_v66_windows_decodes_pairs() {
let mut bytes = 2u32.to_le_bytes().to_vec();
bytes.extend_from_slice(&[0u8; 28]);
bytes.extend_from_slice(&1.0f32.to_le_bytes());
bytes.extend_from_slice(&2.0f32.to_le_bytes());
bytes.extend_from_slice(&3.0f32.to_le_bytes());
bytes.extend_from_slice(&4.0f32.to_le_bytes());
let mut cursor = Cursor::new(bytes);
let windows = read_scan_srm_v66_windows(&mut cursor, 0, 0).unwrap();
assert_eq!(windows, vec![(1.0, 2.0), (3.0, 4.0)]);
}
#[test]
fn search_v63_transition_finds_matching_record() {
let mut data = vec![0u8; 8]; data.extend_from_slice(&1.0f64.to_le_bytes()); data.extend_from_slice(&0.0f64.to_le_bytes()); data.extend_from_slice(&500.0f64.to_le_bytes()); data.extend_from_slice(&300.0f64.to_le_bytes()); data.extend_from_slice(&1.0f64.to_le_bytes()); data.extend_from_slice(&0.02f64.to_le_bytes()); data.extend_from_slice(&25.0f64.to_le_bytes()); data.extend_from_slice(&[0u8; 16]);
let (q1, q3w, ce) = search_v63_transition(&data, 300.0).unwrap();
assert_eq!(q1, 500.0);
assert_eq!(q3w, 1.0);
assert_eq!(ce, 25.0);
}
#[test]
fn search_v63_transition_no_match_returns_none() {
let data = vec![0u8; 128];
assert!(search_v63_transition(&data, 300.0).is_none());
}
#[test]
fn to_mz_intensity_capacity_uses_actual_signal_length_not_untrusted_nbins() {
let profile = Profile {
first_value: 100.0,
step: 0.01,
peak_count: 1,
nbins: u32::MAX,
chunks: vec![ProfileChunk {
first_bin: 0,
signal: vec![1.0, 2.0, 3.0],
fudge: None,
}],
};
let result = profile.to_mz_intensity(&[]);
assert_eq!(result.len(), 3);
}
}