Skip to main content

rustdf/sim/
dia.rs

1use mscore::simulation::noise_rng::noise_rng;
2use mscore::algorithm::isotope::{
3    calculate_precursor_transmission_factor,
4    calculate_transmission_dependent_fragment_ion_isotope_distribution,
5};
6use mscore::data::peptide::{PeptideIon, PeptideProductIonSeriesCollection};
7use mscore::data::spectrum::{IndexedMzSpectrum, MsType, MzSpectrum};
8use mscore::simulation::annotation::{
9    MzSpectrumAnnotated, PeakAnnotation, TimsFrameAnnotated, TimsSpectrumAnnotated,
10};
11use mscore::timstof::collision::{TimsTofCollisionEnergy, TimsTofCollisionEnergyDIA};
12use mscore::timstof::frame::TimsFrame;
13use mscore::timstof::quadrupole::{IonTransmission, TimsTransmissionDIA, WindowTransmission};
14use mscore::timstof::spectrum::TimsSpectrum;
15use std::collections::{BTreeMap, HashSet};
16use std::path::Path;
17use std::sync::Arc;
18
19use rand::Rng;
20use rayon::prelude::*;
21use rayon::ThreadPoolBuilder;
22
23use crate::sim::containers::{IsotopeTransmissionConfig, IsotopeTransmissionMode};
24use crate::sim::projector::{IntensityStage, MzCoordSpace, RenderedEvent, RenderedSpectrum};
25use crate::sim::scheme::{DataMode, InstrumentCapabilities, IsolationWindow};
26use crate::sim::handle::{FragmentIonsWithComplementary, TimsTofSyntheticsDataHandle};
27use crate::sim::precursor::TimsTofSyntheticsPrecursorFrameBuilder;
28
29/// Vendor-neutral per-fragment-series spectrum (P6a MS2 physics kernel): scale one
30/// fragment ion series into its `MzSpectrum` under the isotope-transmission mode,
31/// PRE m/z-noise. This is the fragment physics shared by every instrument — the
32/// transmission GATING (which scans/windows reach here) and the per-scan vs
33/// collapsed aggregation stay in the vendor adapter. `transmission_factor` is the
34/// precursor-scaling factor (PrecursorScaling), `frag_data` the per-(peptide,
35/// charge,CE) complementary data (PerFragment), `transmitted_indices` the
36/// quadrupole-transmitted precursor-isotope indices. Extracted verbatim from the
37/// inline DIA/DDA computation so rendered output is unchanged.
38#[allow(clippy::too_many_arguments)]
39pub(crate) fn fragment_series_spectrum(
40    mode: IsotopeTransmissionMode,
41    fragment_ion_series: &MzSpectrum,
42    series_idx: usize,
43    fraction_events: f32,
44    transmission_factor: f64,
45    frag_data: Option<&FragmentIonsWithComplementary>,
46    transmitted_indices: &HashSet<usize>,
47    max_isotopes: usize,
48) -> MzSpectrum {
49    match mode {
50        IsotopeTransmissionMode::None => fragment_ion_series.clone() * fraction_events as f64,
51        IsotopeTransmissionMode::PrecursorScaling => {
52            fragment_ion_series.clone() * (fraction_events as f64 * transmission_factor)
53        }
54        IsotopeTransmissionMode::PerFragment => {
55            if let Some(frag_data) = frag_data {
56                if series_idx < frag_data.per_fragment_data.len() {
57                    let series_data = &frag_data.per_fragment_data[series_idx];
58                    let mut aggregated_mz: Vec<f64> = Vec::new();
59                    let mut aggregated_intensity: Vec<f64> = Vec::new();
60                    for frag_ion_data in series_data {
61                        let adjusted_dist =
62                            calculate_transmission_dependent_fragment_ion_isotope_distribution(
63                                &frag_ion_data.fragment_distribution,
64                                &frag_ion_data.complementary_distribution,
65                                transmitted_indices,
66                                max_isotopes,
67                            );
68                        for (mz, abundance) in adjusted_dist {
69                            aggregated_mz.push(mz);
70                            aggregated_intensity
71                                .push(abundance * frag_ion_data.predicted_intensity * fraction_events as f64);
72                        }
73                    }
74                    if !aggregated_mz.is_empty() {
75                        MzSpectrum::new(aggregated_mz, aggregated_intensity)
76                    } else {
77                        fragment_ion_series.clone() * fraction_events as f64
78                    }
79                } else {
80                    fragment_ion_series.clone() * fraction_events as f64
81                }
82            } else {
83                fragment_ion_series.clone() * fraction_events as f64
84            }
85        }
86    }
87}
88
89pub struct TimsTofSyntheticsFrameBuilderDIA {
90    pub path: String,
91    pub precursor_frame_builder: TimsTofSyntheticsPrecursorFrameBuilder,
92    pub transmission_settings: TimsTransmissionDIA,
93    pub fragmentation_settings: TimsTofCollisionEnergyDIA,
94    pub fragment_ions:
95        Option<BTreeMap<(u32, i8, i32), (PeptideProductIonSeriesCollection, Vec<MzSpectrum>)>>,
96    pub fragment_ions_annotated: Option<
97        BTreeMap<(u32, i8, i32), (PeptideProductIonSeriesCollection, Vec<MzSpectrumAnnotated>)>,
98    >,
99    /// Configuration for quad-selection dependent isotope transmission (already
100    /// gated by the instrument capabilities at construction — P5e).
101    pub isotope_config: IsotopeTransmissionConfig,
102    /// Fragment ions with complementary data for transmission-dependent calculations
103    pub fragment_ions_with_transmission:
104        Option<BTreeMap<(u32, i8, i32), FragmentIonsWithComplementary>>,
105    /// Physical instrument capabilities (P5e). Default = Bruker timsTOF.
106    pub capabilities: InstrumentCapabilities,
107    /// Master seed for the simulation's noise. 0 keeps the legacy `thread_rng` behaviour
108    /// for API users that never set it; TimSim sets it from the run's `sample_seed`.
109    pub noise_seed: u64,
110}
111
112impl TimsTofSyntheticsFrameBuilderDIA {
113
114    /// Set the master seed for the m/z jitter, making a run's noise reproducible.
115    pub fn set_noise_seed(&mut self, seed: u64) {
116        self.noise_seed = seed;
117    }
118    pub fn new(path: &Path, with_annotations: bool, num_threads: usize) -> rusqlite::Result<Self> {
119        Self::new_with_config(path, with_annotations, num_threads, IsotopeTransmissionConfig::default())
120    }
121
122    pub fn new_with_config(
123        path: &Path,
124        with_annotations: bool,
125        num_threads: usize,
126        isotope_config: IsotopeTransmissionConfig,
127    ) -> rusqlite::Result<Self> {
128        Self::new_with_config_and_source(
129            path,
130            with_annotations,
131            num_threads,
132            isotope_config,
133            &crate::sim::projector::DistributionSource::Columns,
134        )
135    }
136
137    /// As [`new_with_config`], but the precursor builder's occurrence/abundance
138    /// distributions come from `source` (P4: `Columns` default, or the projector).
139    pub fn new_with_config_and_source(
140        path: &Path,
141        with_annotations: bool,
142        num_threads: usize,
143        isotope_config: IsotopeTransmissionConfig,
144        source: &crate::sim::projector::DistributionSource,
145    ) -> rusqlite::Result<Self> {
146        let synthetics = TimsTofSyntheticsPrecursorFrameBuilder::from_source(path, source)?;
147        let handle = TimsTofSyntheticsDataHandle::new(path)?;
148
149        // P5b: refuse to render fragments stored under an incompatible prediction
150        // set (CE encoding the render keying can't resolve). Bruker/legacy pass.
151        handle
152            .read_prediction_set()?
153            .assert_render_compatible()
154            .map_err(|_| rusqlite::Error::InvalidQuery)?;
155
156        let fragment_ions = handle.read_fragment_ions()?;
157
158        // P5e: gate the isotope-transmission config by instrument capabilities.
159        // Default = Bruker timsTOF (no-op); P6 threads real Astral capabilities.
160        let capabilities = InstrumentCapabilities::default();
161        let isotope_config = isotope_config.gated_by(capabilities);
162
163        // get collision energy settings per window group
164        let fragmentation_settings = handle.get_collision_energy_dia();
165        // get ion transmission settings per window group
166        let transmission_settings = handle.get_transmission_dia();
167
168        // Build transmission data if isotope transmission is enabled
169        let fragment_ions_with_transmission = if isotope_config.is_enabled() {
170            Some(TimsTofSyntheticsDataHandle::build_fragment_ions_with_transmission_data(
171                &synthetics.peptides,
172                &fragment_ions,
173                num_threads,
174            ))
175        } else {
176            None
177        };
178
179        match with_annotations {
180            true => {
181                let fragment_ions_annotated =
182                    Some(TimsTofSyntheticsDataHandle::build_fragment_ions_annotated(
183                        &synthetics.peptides,
184                        &fragment_ions,
185                        num_threads,
186                    ));
187                Ok(Self {
188                    noise_seed: 0,
189                    path: path.to_str().unwrap().to_string(),
190                    precursor_frame_builder: synthetics,
191                    transmission_settings,
192                    fragmentation_settings,
193                    fragment_ions: None,
194                    fragment_ions_annotated,
195                    isotope_config,
196                    fragment_ions_with_transmission,
197                    capabilities,
198                })
199            }
200
201            false => {
202                let fragment_ions = Some(TimsTofSyntheticsDataHandle::build_fragment_ions(
203                    &synthetics.peptides,
204                    &fragment_ions,
205                    num_threads,
206                ));
207                Ok(Self {
208                    noise_seed: 0,
209                    path: path.to_str().unwrap().to_string(),
210                    precursor_frame_builder: synthetics,
211                    transmission_settings,
212                    fragmentation_settings,
213                    fragment_ions,
214                    fragment_ions_annotated: None,
215                    isotope_config,
216                    fragment_ions_with_transmission,
217                    capabilities,
218                })
219            }
220        }
221    }
222
223    /// Build a frame for DIA synthetic experiment
224    ///
225    /// # Arguments
226    ///
227    /// * `frame_id` - The frame id
228    /// * `fragmentation` - A boolean indicating if fragmentation is enabled, if false, the frame has same mz distribution as the precursor frame but will be quadrupole filtered
229    ///
230    /// # Returns
231    ///
232    /// A TimsFrame
233    ///
234    pub fn build_frame(
235        &self,
236        frame_id: u32,
237        fragmentation: bool,
238        mz_noise_precursor: bool,
239        uniform: bool,
240        precursor_noise_ppm: f64,
241        mz_noise_fragment: bool,
242        fragment_noise_ppm: f64,
243        right_drag: bool,
244    ) -> TimsFrame {
245        // determine if the frame is a precursor frame
246        match self
247            .precursor_frame_builder
248            .precursor_frame_id_set
249            .contains(&frame_id)
250        {
251            true => self.build_ms1_frame(
252                frame_id,
253                mz_noise_precursor,
254                uniform,
255                precursor_noise_ppm,
256                right_drag,
257            ),
258            false => self.build_ms2_frame(
259                frame_id,
260                fragmentation,
261                mz_noise_fragment,
262                uniform,
263                fragment_noise_ppm,
264                right_drag,
265            ),
266        }
267    }
268
269    pub fn build_frame_annotated(
270        &self,
271        frame_id: u32,
272        fragmentation: bool,
273        mz_noise_precursor: bool,
274        uniform: bool,
275        precursor_noise_ppm: f64,
276        mz_noise_fragment: bool,
277        fragment_noise_ppm: f64,
278        right_drag: bool,
279    ) -> TimsFrameAnnotated {
280        match self
281            .precursor_frame_builder
282            .precursor_frame_id_set
283            .contains(&frame_id)
284        {
285            true => self.build_ms1_frame_annotated(
286                frame_id,
287                mz_noise_precursor,
288                uniform,
289                precursor_noise_ppm,
290                right_drag,
291            ),
292            false => self.build_ms2_frame_annotated(
293                frame_id,
294                fragmentation,
295                mz_noise_fragment,
296                uniform,
297                fragment_noise_ppm,
298                right_drag,
299            ),
300        }
301    }
302
303    pub fn get_fragment_ion_ids(&self, precursor_frame_ids: Vec<u32>) -> Vec<u32> {
304        let mut peptide_ids: HashSet<u32> = HashSet::new();
305        // get all peptide ids for the precursor frame ids
306        for frame_id in precursor_frame_ids {
307            for (peptide_id, peptide) in self.precursor_frame_builder.peptides.iter() {
308                if peptide.frame_start <= frame_id && peptide.frame_end >= frame_id {
309                    peptide_ids.insert(*peptide_id);
310                }
311            }
312        }
313        // get all ion ids for the peptide ids
314        let mut result: Vec<u32> = Vec::new();
315        for item in peptide_ids {
316            let ions = self.precursor_frame_builder.ions.get(&item).unwrap();
317            for ion in ions.iter() {
318                result.push(ion.ion_id);
319            }
320        }
321        result
322    }
323
324    pub fn build_frames(
325        &self,
326        frame_ids: Vec<u32>,
327        fragmentation: bool,
328        mz_noise_precursor: bool,
329        uniform: bool,
330        precursor_noise_ppm: f64,
331        mz_noise_fragment: bool,
332        fragment_noise_ppm: f64,
333        right_drag: bool,
334        num_threads: usize,
335    ) -> Vec<TimsFrame> {
336        // Use global thread pool with custom parallelism instead of creating new pool each call
337        let pool = rayon::ThreadPoolBuilder::new()
338            .num_threads(num_threads)
339            .build()
340            .unwrap();
341
342        pool.install(|| {
343            // Use indexed parallel iteration to maintain order, avoiding post-sort
344            let mut tims_frames: Vec<TimsFrame> = Vec::with_capacity(frame_ids.len());
345            // Safety: we're about to fill all elements
346            unsafe { tims_frames.set_len(frame_ids.len()); }
347
348            frame_ids.par_iter().enumerate().for_each(|(idx, frame_id)| {
349                let frame = self.build_frame(
350                    *frame_id,
351                    fragmentation,
352                    mz_noise_precursor,
353                    uniform,
354                    precursor_noise_ppm,
355                    mz_noise_fragment,
356                    fragment_noise_ppm,
357                    right_drag,
358                );
359                // Safety: each index is unique due to enumerate
360                unsafe {
361                    let ptr = tims_frames.as_ptr() as *mut TimsFrame;
362                    std::ptr::write(ptr.add(idx), frame);
363                }
364            });
365
366            tims_frames
367        })
368    }
369
370    pub fn build_frames_annotated(
371        &self,
372        frame_ids: Vec<u32>,
373        fragmentation: bool,
374        mz_noise_precursor: bool,
375        uniform: bool,
376        precursor_noise_ppm: f64,
377        mz_noise_fragment: bool,
378        fragment_noise_ppm: f64,
379        right_drag: bool,
380        num_threads: usize,
381    ) -> Vec<TimsFrameAnnotated> {
382        // Use thread pool with custom parallelism
383        let pool = rayon::ThreadPoolBuilder::new()
384            .num_threads(num_threads)
385            .build()
386            .unwrap();
387
388        pool.install(|| {
389            // Use indexed parallel iteration to maintain order, avoiding post-sort
390            let mut tims_frames: Vec<TimsFrameAnnotated> = Vec::with_capacity(frame_ids.len());
391            unsafe { tims_frames.set_len(frame_ids.len()); }
392
393            frame_ids.par_iter().enumerate().for_each(|(idx, frame_id)| {
394                let frame = self.build_frame_annotated(
395                    *frame_id,
396                    fragmentation,
397                    mz_noise_precursor,
398                    uniform,
399                    precursor_noise_ppm,
400                    mz_noise_fragment,
401                    fragment_noise_ppm,
402                    right_drag,
403                );
404                unsafe {
405                    let ptr = tims_frames.as_ptr() as *mut TimsFrameAnnotated;
406                    std::ptr::write(ptr.add(idx), frame);
407                }
408            });
409
410            tims_frames
411        })
412    }
413
414    fn build_ms1_frame(
415        &self,
416        frame_id: u32,
417        mz_noise_precursor: bool,
418        uniform: bool,
419        precursor_ppm: f64,
420        right_drag: bool,
421    ) -> TimsFrame {
422        let mut tims_frame = self.precursor_frame_builder.build_precursor_frame(
423            frame_id,
424            mz_noise_precursor,
425            uniform,
426            precursor_ppm,
427            right_drag,
428        );
429        let intensities_rounded = tims_frame
430            .ims_frame
431            .intensity
432            .iter()
433            .map(|x| x.round())
434            .collect::<Vec<_>>();
435        tims_frame.ims_frame.intensity = Arc::new(intensities_rounded);
436        tims_frame
437    }
438
439    fn build_ms1_frame_annotated(
440        &self,
441        frame_id: u32,
442        mz_noise_precursor: bool,
443        uniform: bool,
444        precursor_ppm: f64,
445        right_drag: bool,
446    ) -> TimsFrameAnnotated {
447        let mut tims_frame = self
448            .precursor_frame_builder
449            .build_precursor_frame_annotated(
450                frame_id,
451                mz_noise_precursor,
452                uniform,
453                precursor_ppm,
454                right_drag,
455            );
456        let intensities_rounded = tims_frame
457            .intensity
458            .iter()
459            .map(|x| x.round())
460            .collect::<Vec<_>>();
461        tims_frame.intensity = intensities_rounded;
462        tims_frame
463    }
464
465    fn build_ms2_frame(
466        &self,
467        frame_id: u32,
468        fragmentation: bool,
469        mz_noise_fragment: bool,
470        uniform: bool,
471        fragment_ppm: f64,
472        right_drag: bool,
473    ) -> TimsFrame {
474        match fragmentation {
475            false => {
476                let mut frame = self.transmission_settings.transmit_tims_frame(
477                    &self.build_ms1_frame(
478                        frame_id,
479                        mz_noise_fragment,
480                        uniform,
481                        fragment_ppm,
482                        right_drag,
483                    ),
484                    None,
485                );
486                let intensities_rounded = frame
487                    .ims_frame
488                    .intensity
489                    .iter()
490                    .map(|x| x.round())
491                    .collect::<Vec<_>>();
492                frame.ims_frame.intensity = Arc::new(intensities_rounded);
493                frame.ms_type = MsType::FragmentDia;
494                frame
495            }
496            true => {
497                let mut frame = self.build_fragment_frame(
498                    frame_id,
499                    &self.fragment_ions.as_ref().unwrap(),
500                    mz_noise_fragment,
501                    uniform,
502                    fragment_ppm,
503                    None,
504                    None,
505                    None,
506                    Some(right_drag),
507                );
508                let intensities_rounded = frame
509                    .ims_frame
510                    .intensity
511                    .iter()
512                    .map(|x| x.round())
513                    .collect::<Vec<_>>();
514                frame.ims_frame.intensity = Arc::new(intensities_rounded);
515                frame
516            }
517        }
518    }
519
520    fn build_ms2_frame_annotated(
521        &self,
522        frame_id: u32,
523        fragmentation: bool,
524        mz_noise_fragment: bool,
525        uniform: bool,
526        fragment_ppm: f64,
527        right_drag: bool,
528    ) -> TimsFrameAnnotated {
529        match fragmentation {
530            false => {
531                let mut frame = self.transmission_settings.transmit_tims_frame_annotated(
532                    &self.build_ms1_frame_annotated(
533                        frame_id,
534                        mz_noise_fragment,
535                        uniform,
536                        fragment_ppm,
537                        right_drag,
538                    ),
539                    None,
540                );
541                let intensities_rounded = frame
542                    .intensity
543                    .iter()
544                    .map(|x| x.round())
545                    .collect::<Vec<_>>();
546                frame.intensity = intensities_rounded;
547                frame.ms_type = MsType::FragmentDia;
548                frame
549            }
550            true => {
551                let mut frame = self.build_fragment_frame_annotated(
552                    frame_id,
553                    &self.fragment_ions_annotated.as_ref().unwrap(),
554                    mz_noise_fragment,
555                    uniform,
556                    fragment_ppm,
557                    None,
558                    None,
559                    None,
560                    Some(right_drag),
561                );
562                let intensities_rounded = frame
563                    .intensity
564                    .iter()
565                    .map(|x| x.round())
566                    .collect::<Vec<_>>();
567                frame.intensity = intensities_rounded;
568                frame
569            }
570        }
571    }
572
573    /// Build a fragment frame
574    ///
575    /// # Arguments
576    ///
577    /// * `frame_id` - The frame id
578    /// * `mz_min` - The minimum m/z value in fragment spectrum
579    /// * `mz_max` - The maximum m/z value in fragment spectrum
580    /// * `intensity_min` - The minimum intensity value in fragment spectrum
581    ///
582    /// # Returns
583    ///
584    /// A TimsFrame
585    ///
586    fn build_fragment_frame(
587        &self,
588        frame_id: u32,
589        fragment_ions: &BTreeMap<
590            (u32, i8, i32),
591            (PeptideProductIonSeriesCollection, Vec<MzSpectrum>),
592        >,
593        mz_noise_fragment: bool,
594        uniform: bool,
595        fragment_ppm: f64,
596        mz_min: Option<f64>,
597        mz_max: Option<f64>,
598        intensity_min: Option<f64>,
599        right_drag: Option<bool>,
600    ) -> TimsFrame {
601        // One RNG per frame, keyed by the master seed and the frame id, so the m/z jitter
602        // does not depend on which thread builds the frame. See mscore::simulation::noise_rng.
603        let mut rng = noise_rng(self.noise_seed, &[frame_id as u64]);
604        // Cache frame-level lookups once
605        let ms_type = if self.precursor_frame_builder.precursor_frame_id_set.contains(&frame_id) {
606            MsType::Unknown
607        } else {
608            MsType::FragmentDia
609        };
610
611        let rt = *self.precursor_frame_builder.frame_to_rt.get(&frame_id).unwrap() as f64;
612        let right_drag_val = right_drag.unwrap_or(false);
613        let mz_min_val = mz_min.unwrap_or(100.0);
614        let mz_max_val = mz_max.unwrap_or(1700.0);
615        let intensity_min_val = intensity_min.unwrap_or(1.0);
616
617        // Use single lookup instead of contains_key + get
618        let Some((peptide_ids, frame_abundances)) = self
619            .precursor_frame_builder
620            .frame_to_abundances
621            .get(&frame_id)
622        else {
623            return TimsFrame::new(frame_id as i32, ms_type, rt, vec![], vec![], vec![], vec![], vec![]);
624        };
625
626        // Preallocate with estimated capacity
627        let estimated_capacity = peptide_ids.len() * 4;
628        let mut tims_spectra: Vec<TimsSpectrum> = Vec::with_capacity(estimated_capacity);
629
630        // Go over all peptides in the frame with their respective abundances
631        for (peptide_id, frame_abundance) in peptide_ids.iter().zip(frame_abundances.iter()) {
632            // Single lookup instead of contains_key + get
633            let Some((ion_abundances, scan_occurrences, scan_abundances, charges, spectra)) = self
634                .precursor_frame_builder
635                .peptide_to_ions
636                .get(peptide_id)
637            else {
638                continue;
639            };
640
641            // Cache peptide-level lookup
642            let total_events = *self.precursor_frame_builder.peptide_to_events.get(peptide_id).unwrap();
643
644            for (index, ion_abundance) in ion_abundances.iter().enumerate() {
645                let all_scan_occurrence = &scan_occurrences[index];
646                let all_scan_abundance = &scan_abundances[index];
647                let spectrum = &spectra[index];
648                let charge_state = charges[index];
649
650                for (scan, scan_abundance) in all_scan_occurrence.iter().zip(all_scan_abundance.iter()) {
651                    // Get transmitted isotope indices based on config mode
652                    let (is_transmitted, transmitted_indices) = match self.isotope_config.mode {
653                        IsotopeTransmissionMode::None => {
654                            // Standard check without indices
655                            let any = self.transmission_settings.any_transmitted(
656                                frame_id as i32,
657                                *scan as i32,
658                                &spectrum.mz,
659                                None,
660                            );
661                            (any, HashSet::new())
662                        },
663                        IsotopeTransmissionMode::PrecursorScaling | IsotopeTransmissionMode::PerFragment => {
664                            let indices = self.transmission_settings.get_transmission_set(
665                                frame_id as i32,
666                                *scan as i32,
667                                &spectrum.mz,
668                                Some(self.isotope_config.min_probability),
669                            );
670                            (!indices.is_empty(), indices)
671                        },
672                    };
673
674                    if !is_transmitted {
675                        continue;
676                    }
677
678                    // Calculate abundance factor (total_events cached above)
679                    let fraction_events = frame_abundance * scan_abundance * ion_abundance * total_events;
680
681                    // Get collision energy for the ion
682                    let collision_energy = self.fragmentation_settings.get_collision_energy(frame_id as i32, *scan as i32);
683                    // Resolve the fragment CE key, tolerant to ~0.1 eV quantization
684                    // noise (no-op for DIA's pre-rounded window CE; shared with DDA).
685                    let Some(collision_energy_quantized) = crate::sim::handle::resolve_fragment_ce_key(
686                        fragment_ions, *peptide_id, charge_state, collision_energy,
687                    ) else {
688                        // Fail loud if fragments exist for this ion but none near the
689                        // applied CE (prediction set doesn't cover this CE, P5b); a
690                        // precursor with no predicted fragments at all is a legit skip.
691                        if crate::sim::handle::fragment_prefix_exists(fragment_ions, *peptide_id, charge_state) {
692                            panic!(
693                                "DIA fragment lookup miss: peptide {} charge {} applied CE {:.4} eV \
694                                 has predicted fragments, but none within 0.1 eV — the prediction \
695                                 set does not cover this instrument's collision energy",
696                                *peptide_id, charge_state, collision_energy,
697                            );
698                        }
699                        continue;
700                    };
701                    let (_, fragment_series_vec) = fragment_ions
702                        .get(&(*peptide_id, charge_state, collision_energy_quantized))
703                        .expect("resolve_fragment_ce_key returned a present key");
704
705                    // Cache scan mobility lookup
706                    let scan_mobility = *self.precursor_frame_builder.scan_to_mobility.get(scan).unwrap() as f64;
707
708                    // Calculate transmission factor for PrecursorScaling mode
709                    let transmission_factor = if self.isotope_config.mode == IsotopeTransmissionMode::PrecursorScaling {
710                        if let Some(comp_data) = self.fragment_ions_with_transmission.as_ref() {
711                            if let Some(frag_data) = comp_data.get(&(*peptide_id, charge_state, collision_energy_quantized)) {
712                                calculate_precursor_transmission_factor(
713                                    &frag_data.precursor_isotope_distribution,
714                                    &transmitted_indices,
715                                )
716                            } else {
717                                1.0
718                            }
719                        } else {
720                            1.0
721                        }
722                    } else {
723                        1.0
724                    };
725
726                    // Complementary per-(peptide,charge,CE) data for PerFragment;
727                    // looked up once (same for all series) and passed to the kernel.
728                    let frag_data = self
729                        .fragment_ions_with_transmission
730                        .as_ref()
731                        .and_then(|c| c.get(&(*peptide_id, charge_state, collision_energy_quantized)));
732
733                    for (series_idx, fragment_ion_series) in fragment_series_vec.iter().enumerate() {
734                        let final_spectrum = fragment_series_spectrum(
735                            self.isotope_config.mode,
736                            fragment_ion_series,
737                            series_idx,
738                            fraction_events,
739                            transmission_factor,
740                            frag_data,
741                            &transmitted_indices,
742                            self.isotope_config.max_isotopes,
743                        );
744
745                        let mz_spectrum = if mz_noise_fragment {
746                            if uniform {
747                                final_spectrum.add_mz_noise_uniform_with_rng(fragment_ppm, right_drag_val, &mut rng)
748                            } else {
749                                final_spectrum.add_mz_noise_normal_with_rng(fragment_ppm, &mut rng)
750                            }
751                        } else {
752                            final_spectrum
753                        };
754
755                        let spectrum_len = mz_spectrum.mz.len();
756                        tims_spectra.push(TimsSpectrum::new(
757                            frame_id as i32,
758                            *scan as i32,
759                            rt,
760                            scan_mobility,
761                            ms_type.clone(),
762                            IndexedMzSpectrum::from_mz_spectrum(
763                                vec![0; spectrum_len],
764                                mz_spectrum,
765                            ).filter_ranged(100.0, 1700.0, 1.0, 1e9),
766                        ));
767                    }
768
769                    // Add unfragmented precursor ions (survival) if configured
770                    if self.isotope_config.has_precursor_survival() {
771                        let mut rng = rand::thread_rng();
772                        let survival_fraction = rng.gen_range(
773                            self.isotope_config.precursor_survival_min
774                            ..=self.isotope_config.precursor_survival_max
775                        );
776
777                        if survival_fraction > 0.0 {
778                            // Transmit the precursor spectrum through the quadrupole
779                            let precursor_transmitted = self.transmission_settings.transmit_spectrum(
780                                frame_id as i32,
781                                *scan as i32,
782                                spectrum.clone(),
783                                Some(self.isotope_config.min_probability),
784                            );
785
786                            if !precursor_transmitted.mz.is_empty() {
787                                // Scale by survival fraction and event count
788                                let precursor_scaled = precursor_transmitted * (fraction_events as f64 * survival_fraction);
789
790                                let precursor_mz_spectrum = if mz_noise_fragment {
791                                    if uniform {
792                                        precursor_scaled.add_mz_noise_uniform_with_rng(fragment_ppm, right_drag_val, &mut rng)
793                                    } else {
794                                        precursor_scaled.add_mz_noise_normal_with_rng(fragment_ppm, &mut rng)
795                                    }
796                                } else {
797                                    precursor_scaled
798                                };
799
800                                let precursor_len = precursor_mz_spectrum.mz.len();
801                                tims_spectra.push(TimsSpectrum::new(
802                                    frame_id as i32,
803                                    *scan as i32,
804                                    rt,
805                                    scan_mobility,
806                                    ms_type.clone(),
807                                    IndexedMzSpectrum::from_mz_spectrum(
808                                        vec![0; precursor_len],
809                                        precursor_mz_spectrum,
810                                    ).filter_ranged(100.0, 1700.0, 1.0, 1e9),
811                                ));
812                            }
813                        }
814                    }
815                }
816            }
817        }
818
819        if tims_spectra.is_empty() {
820            return TimsFrame::new(frame_id as i32, ms_type, rt, vec![], vec![], vec![], vec![], vec![]);
821        }
822
823        let tims_frame = TimsFrame::from_tims_spectra(tims_spectra);
824        tims_frame.filter_ranged(
825            mz_min_val,
826            mz_max_val,
827            0,
828            1000,
829            0.0,
830            10.0,
831            intensity_min_val,
832            1e9,
833            0,
834            i32::MAX,
835        )
836    }
837
838    /// Render a fragment (MS2) frame as a vendor-neutral [`RenderedEvent::Scan`]
839    /// for a non-IMS instrument (P6e): collapse ONE isolation window's fragment
840    /// signal into the single MS2 spectrum an Astral/Orbitrap records.
841    ///
842    /// Astral has no ion-mobility axis, so (a) each precursor ion is gated into the
843    /// window ONCE via [`WindowTransmission::any_transmitted`] on its isotope
844    /// envelope (an m/z window, NOT a scan range), and (b) every transmitted ion
845    /// contributes its fragment series at the FULL mobility marginal
846    /// `frame_abundance × ion_abundance × total_events` (scan factor folded to 1.0)
847    /// — the same marginal contract as the MS1 render
848    /// ([`TimsTofSyntheticsPrecursorFrameBuilder::precursor_scan_marginal_spectrum`]),
849    /// NOT the Bruker per-scan `scan_abundance` sum. The isotope-transmission mode
850    /// is `None` (Astral capabilities gate it off), so the shared
851    /// [`fragment_series_spectrum`] kernel does the per-series scaling. `nce` is the
852    /// window's normalized collision energy: both the CE key the stored fragments
853    /// were predicted at AND the applied CE (the keying is unit-agnostic). A
854    /// precursor whose fragments are not stored near `nce` is SKIPPED (predicted at
855    /// a different window's CE) — compatibility is enforced once at registration
856    /// (P6d), so the render core never panics on a per-precursor miss. Returns an
857    /// empty-spectrum MS2 `Scan` when nothing is transmitted / no fragments exist
858    /// — the writer must still consume and clear that template slot (zero residual).
859    ///
860    /// Fragments ONLY: precursor-survival signal is intentionally not modelled here
861    /// (it is stochastic, which would break this deterministic render), and an
862    /// Astral run that configures `precursor_survival_*` is rejected at config load
863    /// rather than silently dropping it.
864    pub fn render_fragment_scan(
865        &self,
866        frame_id: u32,
867        window: &WindowTransmission,
868        nce: f64,
869        data_mode: DataMode,
870    ) -> RenderedEvent {
871        let rt = *self
872            .precursor_frame_builder
873            .frame_to_rt
874            .get(&frame_id)
875            .unwrap_or(&0.0) as f64;
876        let isolation = Some(IsolationWindow {
877            center_mz: window.center_mz,
878            width_mz: window.width_mz,
879        });
880        let make_scan = |peaks: MzSpectrum| RenderedEvent::Scan {
881            ms_level: 2,
882            retention_time_s: rt,
883            isolation,
884            spectrum: RenderedSpectrum {
885                mz: (*peaks.mz).clone(),
886                intensity: (*peaks.intensity).clone(),
887                coords: MzCoordSpace::Physical,
888                mode: data_mode,
889                detector_applied: false,
890                stage: IntensityStage::Transmitted,
891            },
892        };
893
894        let (Some(fragment_ions), Some((peptide_ids, frame_abundances))) = (
895            self.fragment_ions.as_ref(),
896            self.precursor_frame_builder.frame_to_abundances.get(&frame_id),
897        ) else {
898            return make_scan(MzSpectrum::from_collection(vec![]));
899        };
900
901        let mut specs: Vec<MzSpectrum> = Vec::new();
902        for (peptide_id, frame_abundance) in peptide_ids.iter().zip(frame_abundances.iter()) {
903            let Some((ion_abundances, _scan_occ, _scan_abu, charges, spectra)) =
904                self.precursor_frame_builder.peptide_to_ions.get(peptide_id)
905            else {
906                continue;
907            };
908            let total_events = *self
909                .precursor_frame_builder
910                .peptide_to_events
911                .get(peptide_id)
912                .unwrap();
913            for (index, ion_abundance) in ion_abundances.iter().enumerate() {
914                let spectrum = &spectra[index];
915                let charge_state = charges[index];
916                // Quadrupole window gating (m/z), once per ion — no IMS scans.
917                if !window.any_transmitted(&spectrum.mz, None) {
918                    continue;
919                }
920                // Resolve the fragment CE key at the window NCE (unit-agnostic).
921                // A miss SKIPS this precursor for this window — unlike the Bruker
922                // same-instrument frame builder (which panics on a coverage gap),
923                // the Astral render gates precursors by m/z ALONE (no mobility), so
924                // a precursor whose fragments were predicted at a different window's
925                // CE legitimately does not contribute here. Prediction-set/instrument
926                // compatibility is enforced once, up front, at registration (P6d).
927                let Some(ce_key) = crate::sim::handle::resolve_fragment_ce_key(
928                    fragment_ions,
929                    *peptide_id,
930                    charge_state,
931                    nce,
932                ) else {
933                    continue;
934                };
935                let (_, fragment_series_vec) = fragment_ions
936                    .get(&(*peptide_id, charge_state, ce_key))
937                    .expect("resolve_fragment_ce_key returned a present key");
938
939                // Full mobility marginal: scan factor folded to 1.0 (no grid).
940                let fraction_events = frame_abundance * ion_abundance * total_events;
941                for (series_idx, fragment_ion_series) in fragment_series_vec.iter().enumerate() {
942                    specs.push(fragment_series_spectrum(
943                        IsotopeTransmissionMode::None,
944                        fragment_ion_series,
945                        series_idx,
946                        fraction_events,
947                        1.0,
948                        None,
949                        &HashSet::new(),
950                        self.isotope_config.max_isotopes,
951                    ));
952                }
953            }
954        }
955        make_scan(MzSpectrum::from_collection(specs))
956    }
957
958    pub fn build_fragment_frame_annotated(
959        &self,
960        frame_id: u32,
961        fragment_ions: &BTreeMap<
962            (u32, i8, i32),
963            (PeptideProductIonSeriesCollection, Vec<MzSpectrumAnnotated>),
964        >,
965        mz_noise_fragment: bool,
966        uniform: bool,
967        fragment_ppm: f64,
968        mz_min: Option<f64>,
969        mz_max: Option<f64>,
970        intensity_min: Option<f64>,
971        right_drag: Option<bool>,
972    ) -> TimsFrameAnnotated {
973        // One RNG per frame, keyed by the master seed and the frame id, so the m/z jitter
974        // does not depend on which thread builds the frame. See mscore::simulation::noise_rng.
975        let mut rng = noise_rng(self.noise_seed, &[frame_id as u64]);
976        // Cache frame-level lookups
977        let ms_type = if self.precursor_frame_builder.precursor_frame_id_set.contains(&frame_id) {
978            MsType::Unknown
979        } else {
980            MsType::FragmentDia
981        };
982
983        let rt = *self.precursor_frame_builder.frame_to_rt.get(&frame_id).unwrap() as f64;
984        let right_drag_val = right_drag.unwrap_or(false);
985        let mz_min_val = mz_min.unwrap_or(100.0);
986        let mz_max_val = mz_max.unwrap_or(1700.0);
987        let intensity_min_val = intensity_min.unwrap_or(1.0);
988
989        // Single lookup instead of contains_key + get
990        let Some((peptide_ids, frame_abundances)) = self
991            .precursor_frame_builder
992            .frame_to_abundances
993            .get(&frame_id)
994        else {
995            return TimsFrameAnnotated::new(frame_id as i32, rt, ms_type, vec![], vec![], vec![], vec![], vec![], vec![]);
996        };
997
998        // Preallocate with estimated capacity
999        let estimated_capacity = peptide_ids.len() * 4;
1000        let mut tims_spectra: Vec<TimsSpectrumAnnotated> = Vec::with_capacity(estimated_capacity);
1001
1002        for (peptide_id, frame_abundance) in peptide_ids.iter().zip(frame_abundances.iter()) {
1003            // Single lookup
1004            let Some((ion_abundances, scan_occurrences, scan_abundances, charges, _)) = self
1005                .precursor_frame_builder
1006                .peptide_to_ions
1007                .get(peptide_id)
1008            else {
1009                continue;
1010            };
1011
1012            // Cache peptide-level lookups
1013            let total_events = *self.precursor_frame_builder.peptide_to_events.get(peptide_id).unwrap();
1014            let peptide = self.precursor_frame_builder.peptides.get(peptide_id).unwrap();
1015
1016            for (index, ion_abundance) in ion_abundances.iter().enumerate() {
1017                let all_scan_occurrence = &scan_occurrences[index];
1018                let all_scan_abundance = &scan_abundances[index];
1019                let charge_state = charges[index];
1020
1021                let ion = PeptideIon::new(
1022                    peptide.sequence.sequence.clone(),
1023                    charge_state as i32,
1024                    *ion_abundance as f64,
1025                    Some(*peptide_id as i32),
1026                );
1027                // TODO: make this configurable
1028                let spectrum = ion.calculate_isotopic_spectrum_annotated(1e-3, 1e-8, 200, 1e-4);
1029
1030                for (scan, scan_abundance) in all_scan_occurrence.iter().zip(all_scan_abundance.iter()) {
1031                    if !self.transmission_settings.any_transmitted(
1032                        frame_id as i32,
1033                        *scan as i32,
1034                        &spectrum.mz,
1035                        None,
1036                    ) {
1037                        continue;
1038                    }
1039
1040                    let fraction_events = frame_abundance * scan_abundance * ion_abundance * total_events;
1041
1042                    let collision_energy = self.fragmentation_settings.get_collision_energy(frame_id as i32, *scan as i32);
1043                    // Resolve the fragment CE key, tolerant to ~0.1 eV quantization noise.
1044                    let Some(collision_energy_quantized) = crate::sim::handle::resolve_fragment_ce_key(
1045                        fragment_ions, *peptide_id, charge_state, collision_energy,
1046                    ) else {
1047                        if crate::sim::handle::fragment_prefix_exists(fragment_ions, *peptide_id, charge_state) {
1048                            panic!(
1049                                "DIA (annotated) fragment lookup miss: peptide {} charge {} applied CE \
1050                                 {:.4} eV has predicted fragments, but none within 0.1 eV — the \
1051                                 prediction set does not cover this instrument's collision energy",
1052                                *peptide_id, charge_state, collision_energy,
1053                            );
1054                        }
1055                        continue;
1056                    };
1057                    let (_, fragment_series_vec) = fragment_ions
1058                        .get(&(*peptide_id, charge_state, collision_energy_quantized))
1059                        .expect("resolve_fragment_ce_key returned a present key");
1060
1061                    // Cache scan mobility
1062                    let scan_mobility = *self.precursor_frame_builder.scan_to_mobility.get(scan).unwrap() as f64;
1063
1064                    for fragment_ion_series in fragment_series_vec.iter() {
1065                        let scaled_spec = fragment_ion_series.clone() * fraction_events as f64;
1066
1067                        let mz_spectrum = if mz_noise_fragment {
1068                            if uniform {
1069                                scaled_spec.add_mz_noise_uniform_with_rng(fragment_ppm, right_drag_val, &mut rng)
1070                            } else {
1071                                scaled_spec.add_mz_noise_normal_with_rng(fragment_ppm, &mut rng)
1072                            }
1073                        } else {
1074                            scaled_spec
1075                        };
1076
1077                        let spectrum_len = mz_spectrum.mz.len();
1078                        tims_spectra.push(TimsSpectrumAnnotated::new(
1079                            frame_id as i32,
1080                            *scan,
1081                            rt,
1082                            scan_mobility,
1083                            ms_type.clone(),
1084                            vec![0; spectrum_len],
1085                            mz_spectrum,
1086                        ));
1087                    }
1088
1089                    // Add unfragmented precursor ions (survival) if configured
1090                    if self.isotope_config.has_precursor_survival() {
1091                        let mut rng = rand::thread_rng();
1092                        let survival_fraction = rng.gen_range(
1093                            self.isotope_config.precursor_survival_min
1094                            ..=self.isotope_config.precursor_survival_max
1095                        );
1096
1097                        if survival_fraction > 0.0 {
1098                            // Create a non-annotated spectrum for transmission
1099                            let precursor_mz_spectrum = MzSpectrum::new(spectrum.mz.clone(), spectrum.intensity.clone());
1100
1101                            // Transmit through the quadrupole
1102                            let precursor_transmitted = self.transmission_settings.transmit_spectrum(
1103                                frame_id as i32,
1104                                *scan as i32,
1105                                precursor_mz_spectrum,
1106                                Some(self.isotope_config.min_probability),
1107                            );
1108
1109                            if !precursor_transmitted.mz.is_empty() {
1110                                // Scale by survival fraction and event count
1111                                let precursor_scaled = precursor_transmitted * (fraction_events as f64 * survival_fraction);
1112
1113                                let precursor_final = if mz_noise_fragment {
1114                                    if uniform {
1115                                        precursor_scaled.add_mz_noise_uniform_with_rng(fragment_ppm, right_drag_val, &mut rng)
1116                                    } else {
1117                                        precursor_scaled.add_mz_noise_normal_with_rng(fragment_ppm, &mut rng)
1118                                    }
1119                                } else {
1120                                    precursor_scaled
1121                                };
1122
1123                                // Convert to annotated spectrum (with precursor annotations)
1124                                let annotations: Vec<PeakAnnotation> = precursor_final.mz.iter()
1125                                    .map(|_| PeakAnnotation { contributions: vec![] })
1126                                    .collect();
1127                                let precursor_annotated = MzSpectrumAnnotated::new(
1128                                    precursor_final.mz.to_vec(),
1129                                    precursor_final.intensity.to_vec(),
1130                                    annotations,
1131                                );
1132
1133                                let precursor_len = precursor_annotated.mz.len();
1134                                tims_spectra.push(TimsSpectrumAnnotated::new(
1135                                    frame_id as i32,
1136                                    *scan,
1137                                    rt,
1138                                    scan_mobility,
1139                                    ms_type.clone(),
1140                                    vec![0; precursor_len],
1141                                    precursor_annotated,
1142                                ));
1143                            }
1144                        }
1145                    }
1146                }
1147            }
1148        }
1149
1150        if tims_spectra.is_empty() {
1151            return TimsFrameAnnotated::new(frame_id as i32, rt, ms_type, vec![], vec![], vec![], vec![], vec![], vec![]);
1152        }
1153
1154        TimsFrameAnnotated::from_tims_spectra_annotated(tims_spectra).filter_ranged(
1155            mz_min_val, mz_max_val, 0.0, 10.0, 0, 1000, intensity_min_val, 1e9,
1156        )
1157    }
1158
1159    pub fn get_ion_transmission_matrix(
1160        &self,
1161        peptide_id: u32,
1162        charge: i8,
1163        include_precursor_frames: bool,
1164    ) -> Vec<Vec<f32>> {
1165        let maybe_peptide_sim = self.precursor_frame_builder.peptides.get(&peptide_id);
1166
1167        let mut frame_ids = match maybe_peptide_sim {
1168            Some(maybe_peptide_sim) => maybe_peptide_sim.frame_distribution.occurrence.clone(),
1169            _ => vec![],
1170        };
1171
1172        if !include_precursor_frames {
1173            frame_ids = frame_ids
1174                .iter()
1175                .filter(|frame_id| {
1176                    !self
1177                        .precursor_frame_builder
1178                        .precursor_frame_id_set
1179                        .contains(frame_id)
1180                })
1181                .cloned()
1182                .collect();
1183        }
1184
1185        let ion = self
1186            .precursor_frame_builder
1187            .ions
1188            .get(&peptide_id)
1189            .unwrap()
1190            .iter()
1191            .find(|ion| ion.charge == charge)
1192            .unwrap();
1193        let spectrum = ion.simulated_spectrum.clone();
1194        let scan_distribution = &ion.scan_distribution;
1195
1196        let mut transmission_matrix =
1197            vec![vec![0.0; frame_ids.len()]; scan_distribution.occurrence.len()];
1198
1199        for (frame_index, frame) in frame_ids.iter().enumerate() {
1200            for (scan_index, scan) in scan_distribution.occurrence.iter().enumerate() {
1201                if self.transmission_settings.all_transmitted(
1202                    *frame as i32,
1203                    *scan as i32,
1204                    &spectrum.mz,
1205                    None,
1206                ) {
1207                    transmission_matrix[scan_index][frame_index] = 1.0;
1208                } else if self.transmission_settings.any_transmitted(
1209                    *frame as i32,
1210                    *scan as i32,
1211                    &spectrum.mz,
1212                    None,
1213                ) {
1214                    let transmitted_spectrum = self.transmission_settings.transmit_spectrum(
1215                        *frame as i32,
1216                        *scan as i32,
1217                        spectrum.clone(),
1218                        None,
1219                    );
1220                    let percentage_transmitted = transmitted_spectrum.intensity.iter().sum::<f64>()
1221                        / spectrum.intensity.iter().sum::<f64>();
1222                    transmission_matrix[scan_index][frame_index] = percentage_transmitted as f32;
1223                }
1224            }
1225        }
1226
1227        transmission_matrix
1228    }
1229
1230    pub fn count_number_transmissions(&self, peptide_id: u32, charge: i8) -> (usize, usize) {
1231        let frame_ids: Vec<_> = self
1232            .precursor_frame_builder
1233            .peptides
1234            .get(&peptide_id)
1235            .unwrap()
1236            .frame_distribution
1237            .occurrence
1238            .clone()
1239            .iter()
1240            .filter(|frame_id| {
1241                !self
1242                    .precursor_frame_builder
1243                    .precursor_frame_id_set
1244                    .contains(frame_id)
1245            })
1246            .cloned()
1247            .collect();
1248        let ion = self
1249            .precursor_frame_builder
1250            .ions
1251            .get(&peptide_id)
1252            .unwrap()
1253            .iter()
1254            .find(|ion| ion.charge == charge)
1255            .unwrap();
1256        let spectrum = ion.simulated_spectrum.clone();
1257        let scan_distribution = &ion.scan_distribution;
1258
1259        let mut frame_count = 0;
1260        let mut scan_count = 0;
1261
1262        for frame in frame_ids.iter() {
1263            let mut frame_transmitted = false;
1264            for scan in scan_distribution.occurrence.iter() {
1265                if self.transmission_settings.any_transmitted(
1266                    *frame as i32,
1267                    *scan as i32,
1268                    &spectrum.mz,
1269                    None,
1270                ) {
1271                    frame_transmitted = true;
1272                    scan_count += 1;
1273                }
1274            }
1275            if frame_transmitted {
1276                frame_count += 1;
1277            }
1278        }
1279
1280        (frame_count, scan_count)
1281    }
1282
1283    pub fn count_number_transmissions_parallel(
1284        &self,
1285        peptide_ids: Vec<u32>,
1286        charge: Vec<i8>,
1287        num_threads: usize,
1288    ) -> Vec<(usize, usize)> {
1289        let thread_pool = ThreadPoolBuilder::new()
1290            .num_threads(num_threads)
1291            .build()
1292            .unwrap();
1293        let result: Vec<(usize, usize)> = thread_pool.install(|| {
1294            peptide_ids
1295                .par_iter()
1296                .zip(charge.par_iter())
1297                .map(|(peptide_id, charge)| self.count_number_transmissions(*peptide_id, *charge))
1298                .collect()
1299        });
1300
1301        result
1302    }
1303}
1304
1305impl TimsTofCollisionEnergy for TimsTofSyntheticsFrameBuilderDIA {
1306    fn get_collision_energy(&self, frame_id: i32, scan_id: i32) -> f64 {
1307        self.fragmentation_settings
1308            .get_collision_energy(frame_id, scan_id)
1309    }
1310}
1311
1312#[cfg(test)]
1313mod p6a_fragment_kernel_tests {
1314    use super::*;
1315    use std::collections::HashSet;
1316
1317    // The scaling branches of the shared fragment-series kernel. (PerFragment is a
1318    // verbatim move of the inline computation; no DB simulated with quad-transmission
1319    // data is available to integration-test it, so it is covered by the byte-parity
1320    // gate on the None path + code review.)
1321    #[test]
1322    fn fragment_series_spectrum_scaling_branches() {
1323        let series = MzSpectrum::new(vec![200.0, 300.0], vec![10.0, 20.0]);
1324        let empty: HashSet<usize> = HashSet::new();
1325
1326        // None: scale by fraction_events only.
1327        let none = fragment_series_spectrum(
1328            IsotopeTransmissionMode::None, &series, 0, 2.0, 99.0 /*ignored*/, None, &empty, 10,
1329        );
1330        assert_eq!(*none.mz, vec![200.0, 300.0]);
1331        assert_eq!(*none.intensity, vec![20.0, 40.0]);
1332
1333        // PrecursorScaling: scale by fraction_events * transmission_factor.
1334        let ps = fragment_series_spectrum(
1335            IsotopeTransmissionMode::PrecursorScaling, &series, 0, 2.0, 0.5, None, &empty, 10,
1336        );
1337        assert_eq!(*ps.intensity, vec![10.0, 20.0]); // 10*(2*0.5), 20*(2*0.5)
1338
1339        // PerFragment with no complementary data falls back to None-style scaling.
1340        let pf = fragment_series_spectrum(
1341            IsotopeTransmissionMode::PerFragment, &series, 0, 2.0, 0.5, None, &empty, 10,
1342        );
1343        assert_eq!(*pf.intensity, vec![20.0, 40.0]);
1344    }
1345
1346    use crate::sim::precursor::TimsTofSyntheticsPrecursorFrameBuilder;
1347
1348    /// Hand-built single-frame DIA builder: frame 2 holds peptide 10 (charge 2,
1349    /// precursor m/z 600) at frame_abundance 0.8, total_events 1000, captured on a
1350    /// scan grid Σ scan_abundance = 0.5. Its fragments (one b/y series, a single
1351    /// peak at m/z 200) are stored at NCE 27 (map key round(27*10)=270).
1352    fn one_fragment_frame_dia() -> TimsTofSyntheticsFrameBuilderDIA {
1353        let mut frame_to_abundances = BTreeMap::new();
1354        frame_to_abundances.insert(2u32, (vec![10u32], vec![0.8f32]));
1355        let mut peptide_to_ions = BTreeMap::new();
1356        peptide_to_ions.insert(
1357            10u32,
1358            (
1359                vec![1.0f32],
1360                vec![vec![100u32, 101u32]],
1361                vec![vec![0.3f32, 0.2f32]], // Σ = 0.5 captured (ignored by the marginal)
1362                vec![2i8],
1363                vec![MzSpectrum::new(vec![600.0], vec![1.0])], // precursor envelope
1364            ),
1365        );
1366        let mut peptide_to_events = BTreeMap::new();
1367        peptide_to_events.insert(10u32, 1000.0f32);
1368        let mut frame_to_rt = BTreeMap::new();
1369        frame_to_rt.insert(2u32, 30.0f32);
1370
1371        let precursor = TimsTofSyntheticsPrecursorFrameBuilder {
1372            noise_seed: 0,
1373            ions: BTreeMap::new(),
1374            peptides: BTreeMap::new(),
1375            scans: Vec::new(),
1376            frames: Vec::new(),
1377            precursor_frame_id_set: HashSet::new(),
1378            frame_to_abundances,
1379            peptide_to_ions,
1380            frame_to_rt,
1381            scan_to_mobility: BTreeMap::new(),
1382            peptide_to_events,
1383            ion_id_to_peptide_charge: BTreeMap::new(),
1384        };
1385
1386        let mut fragment_ions = BTreeMap::new();
1387        fragment_ions.insert(
1388            (10u32, 2i8, 270i32), // NCE 27 -> round(27*10)
1389            (
1390                PeptideProductIonSeriesCollection::new(vec![]),
1391                vec![MzSpectrum::new(vec![200.0], vec![1.0])],
1392            ),
1393        );
1394
1395        TimsTofSyntheticsFrameBuilderDIA {
1396            noise_seed: 0,
1397            path: String::new(),
1398            precursor_frame_builder: precursor,
1399            transmission_settings: TimsTransmissionDIA::new(
1400                vec![], vec![], vec![], vec![], vec![], vec![], vec![], None,
1401            ),
1402            fragmentation_settings: TimsTofCollisionEnergyDIA::new(
1403                vec![], vec![], vec![], vec![], vec![], vec![],
1404            ),
1405            fragment_ions: Some(fragment_ions),
1406            fragment_ions_annotated: None,
1407            isotope_config: IsotopeTransmissionConfig::default(),
1408            fragment_ions_with_transmission: None,
1409            capabilities: InstrumentCapabilities::astral(),
1410        }
1411    }
1412
1413    #[test]
1414    fn astral_ms2_render_full_marginal_and_window_gating() {
1415        let b = one_fragment_frame_dia();
1416
1417        // Window covering the precursor (m/z 600): the ion is transmitted once and
1418        // contributes its fragment series at the FULL mobility marginal
1419        // 0.8 * 1.0 * 1000 * 1.0 = 800 (NOT the captured Σ scan_abundance=0.5 -> 400).
1420        let win = WindowTransmission::new(600.0, 50.0, 15.0);
1421        let RenderedEvent::Scan { ms_level, isolation, spectrum, .. } =
1422            b.render_fragment_scan(2, &win, 27.0, DataMode::Centroid)
1423        else {
1424            panic!("MS2 render must be a Scan");
1425        };
1426        assert_eq!(ms_level, 2);
1427        let iso = isolation.expect("MS2 scan carries an isolation window");
1428        assert!((iso.center_mz - 600.0).abs() < 1e-9 && (iso.width_mz - 50.0).abs() < 1e-9);
1429        assert_eq!(spectrum.mz, vec![200.0]);
1430        assert!((spectrum.intensity[0] - 800.0).abs() < 1e-2, "full marginal, got {}", spectrum.intensity[0]);
1431
1432        // A window that does NOT cover the precursor transmits nothing -> empty scan
1433        // (the slot is still authored/cleared by the writer).
1434        let off = WindowTransmission::new(900.0, 50.0, 15.0);
1435        let RenderedEvent::Scan { spectrum: empty, .. } =
1436            b.render_fragment_scan(2, &off, 27.0, DataMode::Centroid)
1437        else {
1438            panic!("must be a Scan");
1439        };
1440        assert!(empty.mz.is_empty(), "no precursor transmitted -> empty MS2 scan");
1441
1442        // Determinism.
1443        let RenderedEvent::Scan { spectrum: s2, .. } =
1444            b.render_fragment_scan(2, &win, 27.0, DataMode::Centroid)
1445        else { panic!() };
1446        assert_eq!(s2.mz, spectrum.mz);
1447        assert_eq!(s2.intensity, spectrum.intensity);
1448    }
1449}