Skip to main content

gondola_core/events/
mc_event.rs

1// The following file is part of gaps-online-software and published 
2// under the GPLv3 license
3
4use crate::prelude::*;
5
6#[derive(Debug, Clone)]
7#[cfg_attr(feature="pybindings", pyclass)] 
8pub struct McEvent {
9  pub run_id   : u32,
10  pub event_id : u32,
11  pub hits     : Vec<McHit>,
12  pub mctree   : McTree, 
13  pub primary  : Tracklet
14}
15
16impl McEvent {
17
18  pub fn new() -> Self {
19    Self {
20      run_id   : 0,
21      event_id : 0,
22      hits     : Vec::<McHit>::new(),
23      mctree   : McTree::new(),
24      primary  : Tracklet::new()
25    }
26  }
27
28  fn create_mc_tree(&mut self) {
29    self.mctree.assemble(&mut self.hits);
30  }
31}
32
33impl Serialization for McEvent {
34  fn from_bytestream(stream : &Vec<u8>,
35                     pos    : &mut usize) 
36    -> Result<Self, SerializationError> {
37    let mut event = Self::new();
38    let head      = parse_u16(stream, pos);
39    if head != Self::HEAD {
40      error!("Decoding of HEAD failed! Got {} instead!", head);
41      return Err(SerializationError::HeadInvalid);
42    }
43    event.run_id   = parse_u32(stream, pos);
44    event.event_id = parse_u32(stream, pos);
45    let nhits      = parse_u16(stream, pos);
46    for _ in 0..nhits {
47      event.hits.push(McHit::from_bytestream(stream, pos)?);
48    }
49    let tail = parse_u16(stream, pos);
50    if tail != Self::TAIL {
51      //error!("Decoding of TAIL failed for version {}! Got {} instead!", version, tail);
52      error!("Decoding of TAIL failed! Got {} instead!", tail);
53      return Err(SerializationError::TailInvalid);
54    }
55    // assemble and sort the mc tree 
56    //event.mctree.assemble(event.hits);
57    event.create_mc_tree();
58    // fill primary properties from the tree
59    if event.mctree.trackmap.contains_key(&0) {
60      if event.mctree.trackmap[&0].len() > 0 {
61        //self.primary.is_infinite   
62
63        event.primary.vertex_mom_x  = event.mctree.trackmap[&0][0].vertex_mom_x;
64        event.primary.vertex_mom_y  = event.mctree.trackmap[&0][0].vertex_mom_y;
65        event.primary.vertex_mom_z  = event.mctree.trackmap[&0][0].vertex_mom_z;
66        event.primary.vertex_x      = event.mctree.trackmap[&0][0].vertex_pos_x;
67        event.primary.vertex_y      = event.mctree.trackmap[&0][0].vertex_pos_y;
68        event.primary.vertex_z      = event.mctree.trackmap[&0][0].vertex_pos_z;
69        event.primary.vertex_energy = event.mctree.trackmap[&0][0].vertex_kin_E;
70      }
71    }
72    Ok(event) 
73  }
74}
75
76impl fmt::Display for McEvent {
77  fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
78    let mut repr = format!("<McEvent ");
79    repr += &(format!("\n run   id : {}",self.run_id));
80    repr += &(format!("\n event id : {}",self.event_id));
81    repr += "\n === PRIMARY ===";
82    repr += &(format!("\n  {}", self.primary));
83    repr += "\n === HITS ===";
84    for h in &self.hits {
85      repr += &(format!("\n {}", h));
86    }
87    repr += "\n === MCTREE ===";
88    repr += &(format!("\n {}", self.mctree));
89    repr += ">";
90    write!(f, "{}", repr)
91  }
92}
93
94impl Frameable for McEvent {
95  const CRFRAMEOBJECT_TYPE : CRFrameObjectType = CRFrameObjectType::McTree;
96}
97
98#[cfg(feature="pybindings")] 
99#[pymethods] 
100impl McEvent {
101
102  /// The properties of the injected primary 
103  #[getter] 
104  #[pyo3(name="primary")]
105  fn get_primary_py(&self) -> Tracklet {
106    self.primary.clone() 
107  }
108
109  /// Get McTruth Hits in the format 
110  /// (x,y,z,time,edep,volumeid,hardware_id)
111  #[getter]
112  fn get_mctruth_hits(&self) -> Vec<(f32, f32, f32, f32, f32, u32, u32)> {
113    let mut hits = Vec::<(f32,f32,f32,f32,f32,u32,u32)>::new();
114    for track_id in self.mctree.trackmap.keys() {
115      for h in &self.mctree.trackmap[track_id] {
116        let mut hp = (h.pos_x, h.pos_y, h.pos_z, h.glob_time,h.kin_E,h.hw_id, h.volume_id);
117        hits.push(hp);
118      }
119    }
120    return hits;
121  }
122  
123  /// Combine all energy depositions in a certain volume 
124  #[getter]
125  fn get_volid_edep_map(&self) -> HashMap<u32, f32> {
126    let mut edep_map = HashMap::<u32,f32>::new();
127    for track_id in self.mctree.trackmap.keys() {
128      for h in &self.mctree.trackmap[track_id] {
129        if edep_map.contains_key(&h.volume_id) {
130          *edep_map.get_mut(&h.volume_id).unwrap() += h.step_edep; 
131        } else {
132          edep_map.insert(h.volume_id, h.step_edep);
133        }
134      }
135    }
136    return edep_map;
137  }
138  
139  /// Combine all energy depositions in a certain volume 
140  #[getter]
141  fn get_hwid_edep_map(&self) -> HashMap<u32, f32> {
142    let mut edep_map = HashMap::<u32,f32>::new();
143    for track_id in self.mctree.trackmap.keys() {
144      for h in &self.mctree.trackmap[track_id] {
145        if edep_map.contains_key(&h.hw_id) {
146          *edep_map.get_mut(&h.hw_id).unwrap() += h.step_edep; 
147        } else {
148          edep_map.insert(h.hw_id, h.step_edep);
149        }
150      }
151    }
152    return edep_map;
153  }
154}
155
156#[cfg(feature="pybindings")]
157pythonize!(McEvent);
158