Belle II Software  release-05-01-25
FANGSDigitizerModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2015 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Igal Jaegle *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <beast/fangs/modules/FANGSDigitizerModule.h>
12 #include <beast/fangs/dataobjects/FANGSSimHit.h>
13 
14 #include <mdst/dataobjects/MCParticle.h>
15 #include <framework/logging/Logger.h>
16 #include <framework/gearbox/GearDir.h>
17 #include <framework/core/RandomNumbers.h>
18 
19 //c++
20 #include <cmath>
21 #include <vector>
22 
23 using namespace Belle2;
24 using namespace fangs;
25 
26 
27 //-----------------------------------------------------------------
28 // Register the Module
29 //-----------------------------------------------------------------
30 REG_MODULE(FANGSDigitizer)
31 
32 //-----------------------------------------------------------------
33 // Implementation
34 //-----------------------------------------------------------------
35 
37 {
38  // Set module properties
39  setDescription("FANGS digitizer module");
40 
41  //Default values are set here. New values can be in FANGS.xml.
42  addParam("LowerTimingCut", m_lowerTimingCut, "Lower timing cut", 0.);
43  addParam("UpperTimingCut", m_upperTimingCut, "Upper timing cut", 1000000.);
44 }
45 
46 FANGSDigitizerModule::~FANGSDigitizerModule()
47 {
48 }
49 
51 {
52  B2INFO("FANGSDigitizer: Initializing");
53  m_fangsHit.registerInDataStore();
54 
55  //get xml data
56  getXMLData();
57 
58  //converter: electron number to TOT part I
59  fctToT_Calib1 = new TF1("fctToT_Calib1", "[0]*(x/[3]+[1])/(x/[3]+[2])", 0., 100000.);
60  fctToT_Calib1->SetParameters(m_TOTA1, m_TOTB1, m_TOTC1, m_TOTQ1);
61  //converter: electron number to TOT part II
62  fctToT_Calib2 = new TF1("fctToT_Calib2", "[0]*(x/[3]+[1])/(x/[3]+[2])", 0., 100000.);
63  fctToT_Calib2->SetParameters(m_TOTA2, m_TOTB2, m_TOTC2, m_TOTQ2);
64 }
65 
67 {
68 }
69 
71 {
72  StoreArray<MCParticle> mcParticles;
73  StoreArray<FANGSSimHit> FANGSSimHits;
74  m_nFANGS = 15;
75  std::vector<double> T0(m_nFANGS,
77 
78  //Find the start time
79  for (const auto& FANGSSimHit : FANGSSimHits) {
80  const int lad = FANGSSimHit.getLadder();
81  const int sen = FANGSSimHit.getSensor();
82  const int detNb = (lad - 1) * 5 + sen - 1;
83  const TVector3 trackPosition = FANGSSimHit.getLocalPosEntry();
84  const double z = trackPosition.Y() + m_sensor_width / 2.; //cm
85  if (z < T0[detNb]) {
86  T0[detNb] = z;
87  }
88  }
89 
90  for (auto& val : T0) {
91  if (m_lowerTimingCut < val && val < m_upperTimingCut) {
92  val = val / m_v_sensor; //cm / (cm / ns) = ns
93  } else {
94  val = -1.;
95  }
96  }
97 
98  //loop on all entries to store in 3D the ionization for each FANGS
99  for (const auto& FANGSSimHit : FANGSSimHits) {
100 
101  const int lad = FANGSSimHit.getLadder();
102  const int sen = FANGSSimHit.getSensor();
103  const int detNb = (lad - 1) * 5 + sen - 1;
104  const TVector3 simHitPosition = FANGSSimHit.getLocalPosEntry();
105  const double edep = FANGSSimHit.getEnergyDep() * 1e9; //GeV to eV
106 
107  const TVector3 chipPosition(//cm
108  simHitPosition.X(),
109  simHitPosition.Z(),
110  simHitPosition.Y() + m_sensor_width / 2.);
111 
112  //If new detector filled the chip
113  if (olddetNb != detNb && m_dchip_map.size() > 0) {
114  Pixelization();
115  olddetNb = detNb;
116  m_dchip_map.clear();
117  m_dchip.clear();
118  m_dchip_detNb_map.clear();
119  m_dchip_pdg_map.clear();
120  m_dchip_trkID_map.clear();
121  }
122 
123  // There is some issue with these cuts! please check, remove for now!
124  //check if ionization within sensitive volume
125  /*if ((-m_ChipColumnX < chipPosition.X() && chipPosition.X() < m_ChipColumnX) &&
126  (-m_ChipRowY < chipPosition.Y() && chipPosition.Y() < m_ChipRowY) &&
127  (0. < chipPosition.Z() && chipPosition.Z() < m_sensor_width) &&
128  (m_lowerTimingCut < T0[detNb] && T0[detNb] < m_upperTimingCut)) {*/
129  if (1) {
130 
131  if (edep < m_Workfct) break;
133  // check if enough energy to ionize
134  else if (edep > m_Workfct) {
135 
136  const double meanEl = edep / m_Workfct;
137  const double sigma = sqrt(m_Fanofac * meanEl);
138  const int NbEle = (int)gRandom->Gaus(meanEl, sigma);
139 
140  int col = (int)((chipPosition.X() + m_ChipColumnX) / (2. * m_ChipColumnX / (double)m_ChipColumnNb));
141  int row = (int)((chipPosition.Y() + m_ChipRowY) / (2. * m_ChipRowY / (double)m_ChipRowNb));
142  int pix = col + m_ChipColumnNb * row;
143  int quT = gRandom->Uniform(-1, 1);
144  int bci = (int)((chipPosition.Z() / m_v_sensor - T0[detNb]) / (double)m_PixelTimeBin) + quT;
145  if (bci < 0)bci = 0;
146 
147  //check if amplified drifted electron are within pixel boundaries
148  if ((0 <= col && col < m_ChipColumnNb) &&
149  (0 <= row && row < m_ChipRowNb) &&
150  (0 <= pix && pix < m_ChipColumnNb * m_ChipRowNb) &&
151  (0 <= bci && bci < MAXtSIZE)) {
152  //store info into 3D array for each FANGSs
153  m_dchip_map[std::tuple<int, int>(col, row)] = 1;
154  m_dchip_detNb_map[std::tuple<int, int>(col, row)] = detNb;
155  m_dchip_pdg_map[std::tuple<int, int>(col, row)] = FANGSSimHit.getPDG();
156  m_dchip_trkID_map[std::tuple<int, int>(col, row)] = FANGSSimHit.gettrkID();
157  m_dchip[std::tuple<int, int, int>(col, row, bci)] += (int)(NbEle);
158  }
159  }
160  }
161  }
162 
163  if (m_dchip_map.size() > 0) Pixelization();
164  m_dchip_map.clear();
165  m_dchip.clear();
166  m_dchip_detNb_map.clear();
167  m_dchip_pdg_map.clear();
168  m_dchip_trkID_map.clear();
169 }
170 
172 {
173  std::vector<int> t0;
174  std::vector<int> col;
175  std::vector<int> row;
176  std::vector<int> ToT;
177  std::vector<int> bci;
178 
179  StoreArray<FANGSHit> FANGSHits;
180 
181  for (auto& keyValuePair : m_dchip_map) {
182  const auto& key = keyValuePair.first;
183  //column
184  int i = std::get<0>(key);
185  //raw
186  int j = std::get<1>(key);
187 
188  if (m_dchip_map[std::tuple<int, int>(i, j)] == 1) {
189 
190  int k0 = 1e9;
191  const int quE = gRandom->Uniform(0, 2);
192  const double thresEl = m_PixelThreshold + gRandom->Uniform(-1.*m_PixelThresholdRMS, 1.*m_PixelThresholdRMS);
193  int kcounter = 0;
194  //determined t0 ie first time above pixel threshold
195  for (auto& keyValuePair2 : m_dchip) {
196  const auto& key2 = keyValuePair2.first;
197  int k = std::get<2>(key2);
198  if (m_dchip[std::tuple<int, int, int>(i, j, k)] > thresEl) {
199  if (k0 > k)k0 = k;
200  kcounter ++;
201  }
202  }
203  //determined nb of bc per pixel
204  //if good t0
205  if (k0 != 1e9) {
206  int ik = 0;
207  int NbOfEl = 0;
208  for (auto& keyValuePair2 : m_dchip) {
209  const auto& key2 = keyValuePair2.first;
210  int k = std::get<2>(key2);
211  //sum up charge with 16 cycles
212  if (ik < 16) {
213  NbOfEl += m_dchip[std::tuple<int, int, int>(i, j, k)];
214  } else {
215  //calculate ToT
216  int tot = -1;
217  if (NbOfEl > thresEl && NbOfEl <= 45.*m_TOTQ1) {
218  tot = (int)fctToT_Calib1->Eval((double)NbOfEl) + quE;
219  } else if (NbOfEl > 45.*m_TOTQ1 && NbOfEl <= 900.*m_TOTQ1) {
220  tot = (int)fctToT_Calib2->Eval((double)NbOfEl);
221  } else if (NbOfEl > 800.*m_TOTQ1) {
222  tot = 14;
223  }
224  if (tot > 13) {
225  tot = 14;
226  }
227  if (tot >= 0) {
228  ToT.push_back(tot);
229  t0.push_back(k0);
230  col.push_back(i);
231  row.push_back(j);
232  bci.push_back(k0);
233  }
234  ik = 0;
235  NbOfEl = 0;
236  }
237  ik++;
238  }
239 
240  if (kcounter < 16) {
241  //calculate ToT
242  int tot = -1;
243  if (NbOfEl > thresEl && NbOfEl <= 45.*m_TOTQ1) {
244  tot = (int)fctToT_Calib1->Eval((double)NbOfEl) + quE;
245  } else if (NbOfEl > 45.*m_TOTQ1 && NbOfEl <= 900.*m_TOTQ1) {
246  tot = (int)fctToT_Calib2->Eval((double)NbOfEl);
247  } else if (NbOfEl > 800.*m_TOTQ1) {
248  tot = 14;
249  }
250  if (tot > 13) {
251  tot = 14;
252  }
253  if (tot >= 0) {
254  ToT.push_back(tot);
255  t0.push_back(k0);
256  col.push_back(i);
257  row.push_back(j);
258  bci.push_back(k0);
259  }
260  }
261 
262  }
263  } //end loop on row
264  } // end loop on col
265 
266  //bool PixHit = false;
267  //if entry
268  if (bci.size() > 0) {
269  //PixHit = true;
270  //find start time
271  sort(t0.begin(), t0.end());
272 
273  //loop on nb of hit
274  for (int j = 0; j < (int)bci.size(); j++) {
275  if ((bci[j] - t0[0]) > (m_PixelTimeBinNb - 1)) {
276  continue;
277  }
278  //create FANGSHit
279  FANGSHits.appendNew(FANGSHit(col[j], row[j], bci[j] - t0[0], ToT[j],
280  m_dchip_detNb_map[std::tuple<int, int>(col[j], row[j])],
281  m_dchip_pdg_map[std::tuple<int, int>(col[j], row[j])],
282  m_dchip_trkID_map[std::tuple<int, int>(col[j], row[j])]));
283  } //end if entry
284  }
285  //return PixHit;
286 }
287 //read tube centers, impulse response, and garfield drift data filename from FANGS.xml
289 {
290  //const GearDir& content;
291  GearDir content = GearDir("/Detector/DetectorComponent[@name=\"FANGS\"]/Content/");
292 
293  //get the location of the tubes
294 
295  m_PixelThreshold = content.getInt("PixelThreshold");
296  m_PixelThresholdRMS = content.getInt("PixelThresholdRMS");
297  m_PixelTimeBinNb = content.getInt("PixelTimeBinNb");
298  m_PixelTimeBin = content.getDouble("PixelTimeBin");
299  m_ChipColumnNb = content.getInt("ChipColumnNb");
300  m_ChipRowNb = content.getInt("ChipRowNb");
301  m_ChipColumnX = content.getDouble("ChipColumnX");
302  m_ChipRowY = content.getDouble("ChipRowY");
303  m_TOTA1 = content.getDouble("TOTA1");
304  m_TOTB1 = content.getDouble("TOTB1");
305  m_TOTC1 = content.getDouble("TOTC1");
306  m_TOTQ1 = content.getDouble("TOTQ1");
307  m_TOTA2 = content.getDouble("TOTA2");
308  m_TOTB2 = content.getDouble("TOTB2");
309  m_TOTC2 = content.getDouble("TOTC2");
310  m_TOTQ2 = content.getDouble("TOTQ2");
311  m_v_sensor = content.getDouble("v_sensor");
312  m_sensor_width = content.getDouble("sensor_width");
313  m_Workfct = content.getDouble("Workfct");
314  m_Fanofac = content.getDouble("Fanofac");
315 
316  B2INFO("FANGSDigitizer: Aquired FANGS locations and gas parameters");
317  B2INFO(" from FANGS.xml. There are " << m_nFANGS << " FANGSs implemented");
318 
319 }
320 
322 {
323 }
324 
326 {
327 }
328 
329 
Belle2::StoreArray::appendNew
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:256
Belle2::fangs::FANGSDigitizerModule::beginRun
virtual void beginRun() override
Called when entering a new run.
Definition: FANGSDigitizerModule.cc:66
Belle2::FANGSSimHit
Class FANGSSimHit - Geant4 simulated hit for the FANGS detector.
Definition: FANGSSimHit.h:38
Belle2::fangs::FANGSDigitizerModule::m_TOTA2
double m_TOTA2
TOT factor A2.
Definition: FANGSDigitizerModule.h:113
Belle2::fangs::FANGSDigitizerModule::m_PixelTimeBin
double m_PixelTimeBin
Pixel time bin.
Definition: FANGSDigitizerModule.h:103
Belle2::fangs::FANGSDigitizerModule::m_ChipRowY
double m_ChipRowY
Chip row y dimension.
Definition: FANGSDigitizerModule.h:99
Belle2::fangs::FANGSDigitizerModule::m_ChipRowNb
int m_ChipRowNb
Chip row number.
Definition: FANGSDigitizerModule.h:95
Belle2::fangs::FANGSDigitizerModule::m_PixelThresholdRMS
int m_PixelThresholdRMS
Pixel threshold RMS.
Definition: FANGSDigitizerModule.h:91
Belle2::fangs::FANGSDigitizerModule::initialize
virtual void initialize() override
Initialize the Module.
Definition: FANGSDigitizerModule.cc:50
Belle2::fangs::FANGSDigitizerModule::m_Fanofac
double m_Fanofac
Fano factor.
Definition: FANGSDigitizerModule.h:123
Belle2::fangs::FANGSDigitizerModule::fctToT_Calib1
TF1 * fctToT_Calib1
Define ToT calib 1.
Definition: FANGSDigitizerModule.h:85
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::fangs::FANGSDigitizerModule::m_ChipColumnNb
int m_ChipColumnNb
Chip column number.
Definition: FANGSDigitizerModule.h:93
Belle2::FANGSSimHit::gettrkID
int gettrkID() const
Return track ID.
Definition: FANGSSimHit.h:66
Belle2::fangs::FANGSDigitizerModule::m_dchip
std::map< std::tuple< int, int, int >, int > m_dchip
chip store arrays
Definition: FANGSDigitizerModule.h:129
Belle2::fangs::FANGSDigitizerModule::m_dchip_map
std::map< std::tuple< int, int >, int > m_dchip_map
chip map arrays
Definition: FANGSDigitizerModule.h:131
Belle2::fangs::FANGSDigitizerModule::m_TOTC1
double m_TOTC1
TOT factor C1.
Definition: FANGSDigitizerModule.h:109
Belle2::fangs::FANGSDigitizerModule::m_Workfct
double m_Workfct
Work function.
Definition: FANGSDigitizerModule.h:121
Belle2::FANGSHit
ClassFANGSHit - digitization simulated hit for the FANGS detector.
Definition: FANGSHit.h:36
Belle2::fangs::FANGSDigitizerModule::m_lowerTimingCut
double m_lowerTimingCut
Lower timing cut.
Definition: FANGSDigitizerModule.h:142
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::fangs::FANGSDigitizerModule::event
virtual void event() override
This method is the core of the module.
Definition: FANGSDigitizerModule.cc:70
Belle2::fangs::FANGSDigitizerModule::m_dchip_trkID_map
std::map< std::tuple< int, int >, int > m_dchip_trkID_map
chip track ID map arrays
Definition: FANGSDigitizerModule.h:135
Belle2::fangs::FANGSDigitizerModule::m_dchip_pdg_map
std::map< std::tuple< int, int >, int > m_dchip_pdg_map
chip pdg map arrays
Definition: FANGSDigitizerModule.h:133
Belle2::fangs::FANGSDigitizerModule::Pixelization
void Pixelization()
Produces the pixelization.
Definition: FANGSDigitizerModule.cc:171
Belle2::fangs::FANGSDigitizerModule::olddetNb
int olddetNb
Old detector counter.
Definition: FANGSDigitizerModule.h:147
Belle2::fangs::FANGSDigitizerModule::m_ChipColumnX
double m_ChipColumnX
Chip column x dimension.
Definition: FANGSDigitizerModule.h:97
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::fangs::FANGSDigitizerModule::m_sensor_width
double m_sensor_width
sensor width
Definition: FANGSDigitizerModule.h:125
Belle2::GearDir
GearDir is the basic class used for accessing the parameter store.
Definition: GearDir.h:41
Belle2::fangs::FANGSDigitizerModule::m_TOTC2
double m_TOTC2
TOT factor C2.
Definition: FANGSDigitizerModule.h:117
Belle2::FANGSSimHit::getSensor
int getSensor() const
Return the Sensor number (starting at 1, increasing with decreasing z)
Definition: FANGSSimHit.h:70
Belle2::fangs::FANGSDigitizerModule::m_TOTB1
double m_TOTB1
TOT factor B1.
Definition: FANGSDigitizerModule.h:107
Belle2::fangs::FANGSDigitizerModule::fctToT_Calib2
TF1 * fctToT_Calib2
Define ToT calib 2.
Definition: FANGSDigitizerModule.h:87
Belle2::fangs::FANGSDigitizerModule::m_PixelTimeBinNb
int m_PixelTimeBinNb
Pixel time number of bin.
Definition: FANGSDigitizerModule.h:101
Belle2::fangs::FANGSDigitizerModule::m_dchip_detNb_map
std::map< std::tuple< int, int >, int > m_dchip_detNb_map
chip Nb map arrays
Definition: FANGSDigitizerModule.h:137
Belle2::FANGSSimHit::getLocalPosEntry
TVector3 getLocalPosEntry() const
Return the local entry track position.
Definition: FANGSSimHit.h:80
Belle2::FANGSSimHit::getEnergyDep
float getEnergyDep() const
Return the energy deposition in electrons.
Definition: FANGSSimHit.h:76
Belle2::fangs::FANGSDigitizerModule::terminate
virtual void terminate() override
This method is called at the end of the event processing.
Definition: FANGSDigitizerModule.cc:325
Belle2::fangs::FANGSDigitizerModule::m_TOTB2
double m_TOTB2
TOT factor B2.
Definition: FANGSDigitizerModule.h:115
Belle2::fangs::FANGSDigitizerModule::m_TOTQ1
double m_TOTQ1
TOT factor Q1.
Definition: FANGSDigitizerModule.h:111
Belle2::FANGSSimHit::getLadder
int getLadder() const
Return the Ladder number (starting at 1, increasing with phi)
Definition: FANGSSimHit.h:68
Belle2::FANGSSimHit::getPDG
int getPDG() const
Return the PDG number of the track.
Definition: FANGSSimHit.h:72
Belle2::fangs::FANGSDigitizerModule::m_upperTimingCut
double m_upperTimingCut
Upper timing cut.
Definition: FANGSDigitizerModule.h:144
Belle2::StoreArray< MCParticle >
Belle2::fangs::FANGSDigitizerModule::m_TOTA1
double m_TOTA1
TOT factor A1.
Definition: FANGSDigitizerModule.h:105
Belle2::fangs::FANGSDigitizerModule::m_nFANGS
int m_nFANGS
number of detectors.
Definition: FANGSDigitizerModule.h:139
Belle2::fangs::FANGSDigitizerModule::m_v_sensor
double m_v_sensor
Drift velocity in sensor.
Definition: FANGSDigitizerModule.h:127
Belle2::fangs::FANGSDigitizerModule::endRun
virtual void endRun() override
This method is called if the current run ends.
Definition: FANGSDigitizerModule.cc:321
Belle2::fangs::FANGSDigitizerModule::m_TOTQ2
double m_TOTQ2
TOT factor Q2.
Definition: FANGSDigitizerModule.h:119
Belle2::fangs::FANGSDigitizerModule::getXMLData
void getXMLData()
array for FANGSHit
Definition: FANGSDigitizerModule.cc:288
Belle2::fangs::FANGSDigitizerModule::m_PixelThreshold
int m_PixelThreshold
Pixel threshold.
Definition: FANGSDigitizerModule.h:89
Belle2::fangs::FANGSDigitizerModule
FANGS ATLAS FE-I4B ASIC chip digitizer.
Definition: FANGSDigitizerModule.h:40