Belle II Software  release-05-02-19
ECLLocalMaximumFinderModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 - Belle II Collaboration *
4  * *
5  * Belle II Local Maximum Finder (LMF). Two stage procedure: First it *
6  * searches a local maximum (energy higher than all neighbours and above *
7  * threshold) and second it uses a method for different *
8  * particle types to decide if the local maximum is a candidate. Each *
9  * candidate will result in one shower/cluster in the end. *
10  * *
11  * Author: The Belle II Collaboration *
12  * Contributors: Torben Ferber (ferber@physics.ubc.ca) *
13  * *
14  * This software is provided "as is" without any warranty. *
15  **************************************************************************/
16 
17 // THIS MODULE
18 #include <ecl/modules/eclLocalMaximumFinder/ECLLocalMaximumFinderModule.h>
19 
20 // ROOT
21 #include "TFile.h"
22 #include "TTree.h"
23 
24 // FRAMEWORK
25 #include <framework/datastore/StoreArray.h>
26 #include <framework/logging/Logger.h>
27 
28 // MDST
29 #include <mdst/dataobjects/MCParticle.h>
30 
31 // ECL
32 #include <ecl/dataobjects/ECLHit.h>
33 #include <ecl/dataobjects/ECLDigit.h>
34 #include <ecl/dataobjects/ECLCalDigit.h>
35 #include <ecl/dataobjects/ECLLocalMaximum.h>
36 #include <ecl/dataobjects/ECLConnectedRegion.h>
37 #include <ecl/geometry/ECLNeighbours.h>
38 #include <ecl/geometry/ECLGeometryPar.h>
39 
40 // NAMESPACE(S)
41 using namespace Belle2;
42 using namespace ECL;
43 
44 //-----------------------------------------------------------------
45 // Register the Module
46 //-----------------------------------------------------------------
47 REG_MODULE(ECLLocalMaximumFinder)
48 REG_MODULE(ECLLocalMaximumFinderPureCsI)
49 
50 //-----------------------------------------------------------------
51 // Implementation
52 //-----------------------------------------------------------------
54  m_mcParticles(mcParticleArrayName()),
55  m_eclHits(eclHitArrayName()),
56  m_eclDigits(eclDigitArrayName()),
57  m_eclCalDigits(eclCalDigitArrayName()),
58  m_eclConnectedRegions(eclConnectedRegionArrayName()),
59  m_eclLocalMaximums(eclLocalMaximumArrayName())
60 {
61  // Set description.
62  setDescription("ECLLocalMaximumFinderModule");
63 
64  // Parallel processing certification.
65  setPropertyFlags(c_ParallelProcessingCertified);
66 
67  // Add module parameters.
68  addParam("energyCut", m_energyCut, "Seed energy cut [MeV], minimum is 5.0 MeV.", 10.0 * Belle2::Unit::MeV);
69  addParam("isTrainingMode", m_isTrainingMode,
70  "Run in training mode (i.e. fill file with MVA input variables and determine MC truth of LM.).", 0);
71  addParam("outfileName", m_outfileName, "Output file name for training file.", std::string("ECLLocalMaximumFinderOutput.root"));
72  addParam("method", m_method, "Method to determine the LM (cut, none, fastbdt).", std::string("none"));
73  addParam("truthFraction", m_truthFraction, "Minimum matched energy fraction truth/rec for the LM.", 0.51);
74  addParam("cutOffset", m_cutOffset, "Cut method specific: Offset. (BaBar: 2.5, high eff: 1.40)", 2.5);
75  addParam("cutSlope", m_cutSlope, "Cut method specific: Slope. (BaBar: 0.5, high eff: 3.0)", 0.5);
76  addParam("cutRatioCorrection", m_cutRatioCorrection, "Cut method specific: Ratio correction.", 0.0);
77 
78 }
79 
81 {
82  ;
83 }
84 
86 {
87  B2DEBUG(200, "ECLLocalMaximumFinderModule::initialize()");
88 
89  m_eclLocalMaximums.registerInDataStore(eclLocalMaximumArrayName());
90  m_eclCalDigits.registerInDataStore(eclCalDigitArrayName());
91  m_eclConnectedRegions.registerInDataStore(eclConnectedRegionArrayName());
92  m_eclConnectedRegions.registerRelationTo(m_eclLocalMaximums);
93  m_eclCalDigits.registerRelationTo(m_eclLocalMaximums);
94 
95  // Check user input.
96  if (m_energyCut < c_minEnergyCut) {
97  B2WARNING("ECLLocalMaximumFinderModule::initialize: Energy threshold too small, resetting to " << c_minEnergyCut << " GeV");
98  m_energyCut = c_minEnergyCut;
99  }
100 
101  if (m_truthFraction < 1e-6) {
102  B2WARNING("ECLLocalMaximumFinderModule::initialize: Matching fraction must be above 1e-6, input: " << m_truthFraction <<
103  " will be reset to 1e-6.");
104  m_truthFraction = 1e-6;
105  }
106 
107  if (m_truthFraction > 1.0) {
108  B2WARNING("ECLLocalMaximumFinderModule::initialize: Matching fraction must be below 1.0, input: " << m_truthFraction <<
109  " will be reset to 1.0.");
110  m_truthFraction = 1.0;
111  }
112 
113  // Geometry instance.
114  m_geom = ECLGeometryPar::Instance();
115 
116  // Initialize neighbour map.
117  m_neighbourMap = new ECLNeighbours("N", 1);
118 
119  // Reset all variables.
120  resetClassifierVariables();
121 
122  // Open output files and declare branches if in training mode. Each file will hold a flat ntuple of training data.
123  m_outfile = NULL;
124  m_tree = NULL;
125  if (m_isTrainingMode > 0) {
126  const int cBufferLength = 500;
127  char tmpBuffer[cBufferLength];
128  int tmpBufferLength = sprintf(tmpBuffer, "%s", m_outfileName.c_str());
129 
130  if (tmpBufferLength >= cBufferLength) {
131  B2FATAL("ECLLocalMaximumFinderModule::initialize(): Output training filename length too long!");
132  }
133 
134  m_outfile = new TFile(tmpBuffer, "RECREATE");
135  m_tree = new TTree("locmax", "locmax");
136  m_tree->Branch("energyRatioNeighbour0", &m_energyRatioNeighbour[0], "energyRatioNeighbour0/F");
137  m_tree->Branch("energyRatioNeighbour1", &m_energyRatioNeighbour[1], "energyRatioNeighbour1/F");
138  m_tree->Branch("energyRatioNeighbour2", &m_energyRatioNeighbour[2], "energyRatioNeighbour2/F");
139  m_tree->Branch("energyRatioNeighbour3", &m_energyRatioNeighbour[3], "energyRatioNeighbour3/F");
140  m_tree->Branch("energyRatioNeighbour4", &m_energyRatioNeighbour[4], "energyRatioNeighbour4/F");
141  m_tree->Branch("energyRatioNeighbour5", &m_energyRatioNeighbour[5], "energyRatioNeighbour5/F");
142  m_tree->Branch("energyRatioNeighbour6", &m_energyRatioNeighbour[6], "energyRatioNeighbour6/F");
143  m_tree->Branch("energyRatioNeighbour7", &m_energyRatioNeighbour[7], "energyRatioNeighbour7/F");
144  m_tree->Branch("energyRatioNeighbour8", &m_energyRatioNeighbour[8], "energyRatioNeighbour8/F");
145  m_tree->Branch("energyRatioNeighbour9", &m_energyRatioNeighbour[9], "energyRatioNeighbour9/F");
146  m_tree->Branch("energyRatioNeighbour10", &m_energyRatioNeighbour[10], "energyRatioNeighbour10/F");
147  m_tree->Branch("energyRatioNeighbour11", &m_energyRatioNeighbour[11], "energyRatioNeighbour11/F");
148  m_tree->Branch("energy", &m_energy, "energy/F");
149  m_tree->Branch("target", &m_target, "target/F");
150  m_tree->Branch("targetindex", &m_targetindex, "targetindex/F");
151  m_tree->Branch("targetpi0index", &m_targetpi0index, "targetpi0index/F");
152  m_tree->Branch("thetaId", &m_thetaId, "thetaId/F");
153  m_tree->Branch("phiId", &m_phiId, "phiId/F");
154  m_tree->Branch("cellId", &m_cellId, "cellId/F");
155  m_tree->Branch("maxNeighbourEnergy", &m_maxNeighbourEnergy, "maxNeighbourEnergy/F");
156  m_tree->Branch("nNeighbours10", &m_nNeighbours10, "nNeighbours10/F");
157  m_tree->Branch("time", &m_time, "time/F");
158  m_tree->Branch("timeResolution", &m_timeResolution, "timeResolution/F");
159  m_tree->Branch("timeFitFailed", &m_timeFitFailed, "timeFitFailed/F");
160  m_tree->Branch("CRId", &m_CRId, "CRId/F");
161  m_tree->Branch("LMId", &m_LMId, "LMId/F");
162  }
163 
164  // initialize the vector that gives the relation between cellid and store array position
165  m_StoreArrPosition.resize(8736 + 1);
166 
167 }
168 
170 {
171  ;
172 }
173 
175 {
176  B2DEBUG(200, "ECLLocalMaximumFinderModule::event()");
177 
178  // Fill a vector that can be used to map cellid -> store array position
179  std::fill_n(m_StoreArrPosition.begin(), m_StoreArrPosition.size(), -1);
180  for (int i = 0; i < m_eclCalDigits.getEntries(); i++) {
181  m_StoreArrPosition[m_eclCalDigits[i]->getCellId()] = i;
182  }
183 
184  // Vector with neighbour ids.
185  std::vector< double > vNeighourEnergies;
186  vNeighourEnergies.resize(c_nMaxNeighbours);
187 
188  // Loop over connected regions.
189  for (const ECLConnectedRegion& aCR : m_eclConnectedRegions) {
190  const int crId = aCR.getCRId();
191  int iLM = 1;
192 
193  // Loop over all entries in this CR.
194  for (const ECLCalDigit& aECLCalDigit : aCR.getRelationsTo<ECLCalDigit>()) {
195 
196  // Check seed energy cut.
197  if (aECLCalDigit.getEnergy() >= m_energyCut) {
198 
199  // Clean up for this candiate (MVA is trained per LM, regardless of CR)
200  std::fill_n(vNeighourEnergies.begin(), vNeighourEnergies.size(),
201  -999); // -999 means later: this digit is just not available in this neighbour definition.
202  resetTrainingVariables();
203  resetClassifierVariables();
204 
205  // Check neighbours: Must be a local energy maximum.
206  bool isLocMax = 1;
207  int neighbourCount = 0;
208  for (auto& neighbourId : m_neighbourMap->getNeighbours(aECLCalDigit.getCellId())) {
209  if (neighbourId == aECLCalDigit.getCellId()) continue; // Skip the center cell to avoid possible floating point issues.
210 
211  const int pos = m_StoreArrPosition[neighbourId]; // Get position in the store array for this digit.
212 
213  double energyNeighbour = 0.0;
214  if (pos >= 0) {
215  energyNeighbour = m_eclCalDigits[pos]->getEnergy(); // Get the energy directly from the store array.
216  vNeighourEnergies[neighbourCount] = energyNeighbour;
217  } else {
218  // Digit does not belong to this CR
219  vNeighourEnergies[neighbourCount] = 0.0;
220  }
221  ++neighbourCount;
222 
223  if (energyNeighbour > aECLCalDigit.getEnergy()) {
224  isLocMax = 0;
225  break;
226  }
227  }
228 
229  // It is a local maximum. Get all variables needed for classification.
230  if (isLocMax) {
231 
232  for (unsigned int npos = 0; npos < vNeighourEnergies.size(); ++npos) {
233  if (vNeighourEnergies[npos] >= 0) {
234  m_energyRatioNeighbour[npos] = static_cast < float >(vNeighourEnergies[npos] / aECLCalDigit.getEnergy());
235  if (vNeighourEnergies[npos] > m_maxNeighbourEnergy) m_maxNeighbourEnergy = vNeighourEnergies[npos];
236  if (vNeighourEnergies[npos] > 1.0 * Belle2::Unit::MeV) ++m_nNeighbours10;
237  } else m_energyRatioNeighbour[npos] = 0.0;
238  }
239  m_time = aECLCalDigit.getTime();
240  m_maxEnergyRatio = m_maxNeighbourEnergy / aECLCalDigit.getEnergy();
241 
242  // Fill training monitoring variables and MC truth information.
243  if (m_isTrainingMode > 0) {
244 
245  m_energy = static_cast < float >(aECLCalDigit.getEnergy());
246  m_cellId = static_cast < float >(aECLCalDigit.getCellId());
247  m_timeResolution = static_cast < float >(aECLCalDigit.getTimeResolution());
248  m_timeFitFailed = static_cast < float >(aECLCalDigit.isFailedFit());
249  m_CRId = static_cast < float >(crId);
250  m_LMId = static_cast < float >(iLM);
251 
252  m_geom->Mapping(m_cellId - 1);
253  m_thetaId = static_cast < float >(m_geom->GetThetaID());
254  m_phiId = static_cast < float >(m_geom->GetPhiID());
255  }
256 
257  if (m_isTrainingMode > 0) { // This requires MC matching before this stage!
258  int pi0index = -1;
259  int maxtype = 0;
260  int maxpos = 0;
261 
262  auto relatedParticlePairs = aECLCalDigit.getRelationsWith<MCParticle>();
263 
264  for (unsigned int irel = 0; irel < relatedParticlePairs.size(); irel++) {
265  const auto particle = relatedParticlePairs.object(irel);
266  const double weight = relatedParticlePairs.weight(irel);
267 
268  int motherpdg = -1;
269  int motherindex = -1;
270  pi0index = -1;
271  getEnteringMother(*particle, motherpdg, motherindex, pi0index);
272  addToSignalEnergy(motherpdg, motherindex, pi0index, weight);
273  }
274 
275  B2DEBUG(175, " -> digt energy: " << aECLCalDigit.getEnergy());
276  B2DEBUG(175, "photon: " << m_signalEnergy[0][0] << " " << m_signalEnergy[0][1]);
277  B2DEBUG(175, "pi0: " << m_signalEnergy[1][0] << " " << m_signalEnergy[1][1]);
278  B2DEBUG(175, "electron: " << m_signalEnergy[2][0] << " " << m_signalEnergy[2][1]);
279  B2DEBUG(175, "muon: " << m_signalEnergy[3][0] << " " << m_signalEnergy[3][1]);
280  B2DEBUG(175, "neutral hadron: " << m_signalEnergy[4][0] << " " << m_signalEnergy[4][1]);
281  B2DEBUG(175, "charged hadron: " << m_signalEnergy[5][0] << " " << m_signalEnergy[5][1]);
282  B2DEBUG(175, "other: " << m_signalEnergy[6][0]);
283 
284  maxtype = 0;
285  maxpos = 0;
286  getMax(maxtype, maxpos);
287 
288  if (maxtype >= 0) {
289  if (m_signalEnergy[maxtype][maxpos] >= m_truthFraction * aECLCalDigit.getEnergy()) {
290  m_target = maxtype;
291  m_targetindex = m_signalId[maxtype][maxpos];
292  m_targetpi0index = pi0index;
293  } else {
294  m_target = 7;
295  }
296  } else {
297  m_target = 7;
298  }
299 
300  m_tree->Fill();
301 
302  } // end training
303 
304  if (m_method == "cut") {
305 
306  B2DEBUG(200, "m_cutSlope: " << m_cutSlope << ", m_nNeighbours10: " << m_nNeighbours10 << ", m_cutOffset: " << m_cutOffset <<
307  ", m_maxNeighbourEnergy: " << m_maxNeighbourEnergy << ", m_cutRatioCorrection: " << m_cutRatioCorrection <<
308  ", aECLCalDigit.getEnergy(): " << aECLCalDigit.getEnergy());
309  B2DEBUG(200, "m_cutSlope * (m_nNeighbours10 - m_cutOffset): " << m_cutSlope * (m_nNeighbours10 - m_cutOffset));
310  B2DEBUG(200, "(m_maxNeighbourEnergy - m_cutRatioCorrection) / (aECLCalDigit.getEnergy() - m_cutRatioCorrection)" <<
311  (m_maxNeighbourEnergy - m_cutRatioCorrection) / (aECLCalDigit.getEnergy() - m_cutRatioCorrection) << "\n");
312 
313  if (m_cutSlope * (m_nNeighbours10 - m_cutOffset) >= (m_maxNeighbourEnergy - m_cutRatioCorrection) /
314  (aECLCalDigit.getEnergy() - m_cutRatioCorrection)) {
315  makeLocalMaximum(aCR, aECLCalDigit.getCellId(), iLM);
316  ++iLM;
317  }
318  } else if (m_method == "none") { // All energy local maximums will become local maximums.
319  makeLocalMaximum(aCR, aECLCalDigit.getCellId(), iLM);
320  ++iLM;
321  }
322 
323  }
324 
325  } // end check energy
326 
327  } // end CalDigit loop
328 
329  // Check if there is at least one local maximum in the CR. If not, make the highest energetic crystal one.
330  if (iLM == 1) {
331 
332  int highestEnergyCellId = -1;
333  double highestEnergy = 0.0;
334 
335  // Loop over all entries in this CR.
336  for (const ECLCalDigit& aECLCalDigit : aCR.getRelationsTo<ECLCalDigit>(eclCalDigitArrayName())) {
337  if (aECLCalDigit.getEnergy() > highestEnergy) {
338  highestEnergyCellId = aECLCalDigit.getCellId();
339  highestEnergy = aECLCalDigit.getEnergy();
340  }
341  } // end CalDigit loop
342 
343  makeLocalMaximum(aCR, highestEnergyCellId, 0);
344 
345  }
346 
347  } // end CR loop
348 
349 }
350 
351 
353 {
354  B2DEBUG(200, "ECLLocalMaximumFinderModule::endRun()");
355 }
356 
357 
359 {
360  B2DEBUG(200, "ECLLocalMaximumFinderModule::terminate()");
361 
362  // If run in trainingmode, write the training file.
363  if (m_outfile and m_tree) {
364  m_outfile->cd();
365  m_tree->Write();
366  m_outfile->Write();
367  m_outfile->Close();
368  delete m_outfile;
369  }
370 
371  if (m_neighbourMap) delete m_neighbourMap;
372 
373 }
374 
375 void ECLLocalMaximumFinderModule::makeLocalMaximum(const ECLConnectedRegion& aCR, const int cellId, const int lmId)
376 {
377  // Set the local maximum dataobject.
378  const auto aLocalMaximum = m_eclLocalMaximums.appendNew();
379 
380  B2DEBUG(175, "ECLLocalMaximumFinderModule::makeLocalMaximum(): local maximum cellid: " << cellId);
381 
382  // Set the id of this local maximum.
383  aLocalMaximum->setLMId(lmId);
384 
385  // Set the cell Id of the digit.
386  aLocalMaximum->setCellId(cellId);
387 
388  // Add relations to ECLConnectedRegion.
389  aCR.addRelationTo(aLocalMaximum);
390 
391 }
392 
394 {
395  m_target = -1;
396  m_targetindex = -1;
397  m_targetpi0index = -1;
398 
399  m_totalSignalEnergy = 0.;
400 
401  for (unsigned int i = 0; i < 10; ++i) {
402  for (unsigned int j = 0; j < 5; ++j) {
403  m_signalEnergy[i][j] = 0.;
404  m_signalId[i][j] = -1;
405  }
406  }
407 
408 }
409 
411 {
412 
413  for (unsigned i = 0; i < c_nMaxNeighbours; ++i) {
414  m_energyRatioNeighbour[i] = 0.;
415  }
416 
417  m_energy = 0.;
418  m_cellId = -1;
419  m_thetaId = -1;
420  m_phiId = -1;
421  m_nNeighbours10 = 0;
422  m_maxNeighbourEnergy = 0.0;
423  m_time = -9999;
424  m_timeResolution = -9999.;
425  m_timeFitFailed = -9999.;
426  m_CRId = -1;
427  m_LMId = -1;
428  m_maxEnergyRatio = 0.0;
429 }
430 
431 
432 void ECLLocalMaximumFinderModule::getEnteringMother(const MCParticle& particle, int& pdg, int& arrayindex, int& pi0arrayindex)
433 {
434 
435  int index = particle.getArrayIndex();
436  int pi0index = -1;
437 
438  while (!isEnteringECL(m_mcParticles[index]->getProductionVertex())) {
439  if (m_mcParticles[index]->getMother()) index = m_mcParticles[index]->getMother()->getArrayIndex();
440  else index = -1;
441  };
442 
443  // For photon mother: are they from a pi0? This can be used to improved overlap/merged pi0 reconstruction.
444  if (m_mcParticles[index]->getPDG() == 22) {
445  if (m_mcParticles[index]->getMother()->getPDG() == 111) {
446  pi0index = m_mcParticles[index]->getMother()->getArrayIndex();
447  }
448  }
449 
450  // Dont include mother if its energy is too low or if its a photon from a neutron interaction.
451  if ((m_mcParticles[index]->getEnergy() < 5.0 * Belle2::Unit::MeV)
452  or (m_mcParticles[index]->getPDG() == 22 and abs(m_mcParticles[index]->getMother()->getPDG()) == 2212)) {
453  pdg = -1;
454  arrayindex = -1;
455  pi0arrayindex = -1;
456  } else {
457  pdg = m_mcParticles[index]->getPDG();
458  arrayindex = index;
459  pi0arrayindex = pi0index;
460  }
461 }
462 
464 {
465 
466  const double theta = vertex.Theta();
467 
468  if (theta > 0.555015 and theta < 2.26369) { //barrel
469  double radius = vertex.Perp();
470  if (radius < 125 * Belle2::Unit::cm) return true;
471  } else if (theta <= 0.555015) { //fwd
472  if (vertex.Z() < 196.16 * Belle2::Unit::cm) return true;
473  } else if (theta >= 2.26369) { //bwd
474  if (vertex.Z() > -102.16 * Belle2::Unit::cm) return true;
475  }
476 
477  return false;
478 }
479 
480 int ECLLocalMaximumFinderModule::getIdPosition(const int type, const int id)
481 {
482 
483  for (int i = 0; i < 5; i++) {
484  if (m_signalId[type][i] == id) return i;
485  if (m_signalId[type][i] == -1) {
486  m_signalId[type][i] = id;
487  return i; //next free one
488  }
489  }
490 
491  return -1;
492 }
493 
494 // of all entries, get the maximum!
495 void ECLLocalMaximumFinderModule::getMax(int& type, int& id)
496 {
497  double maxe = 0.;
498  int maxtype = -1;
499  int maxid = -1;
500 
501  for (unsigned int i = 0; i < 10; ++i) {
502  for (unsigned int j = 0; j < 5; ++j) {
503  if (m_signalEnergy[i][j] > maxe) {
504  maxe = m_signalEnergy[i][j];
505  maxtype = i;
506  maxid = j;
507  }
508  }
509  }
510 
511  type = maxtype;
512  id = maxid;
513 
514 }
515 
516 void ECLLocalMaximumFinderModule::addToSignalEnergy(int& motherpdg, int& motherindex, int& pi0index, const double& weight)
517 {
518 
519  // for the LM training and CR/LM debugging
520  if (motherpdg == 22) {
521  if (pi0index >= 0) { // photon from pi0
522  int idpos = getIdPosition(1, motherindex);
523  m_signalEnergy[1][idpos] += weight;
524  } else {
525  int idpos = getIdPosition(0, motherindex);
526  m_signalEnergy[0][idpos] += weight; // photon from another source
527  }
528  } else if (abs(motherpdg) == 11) { // electron
529  int idpos = getIdPosition(2, motherindex);
530  m_signalEnergy[2][idpos] += weight;
531  } else if (abs(motherpdg) == 13) { // muon
532  int idpos = getIdPosition(3, motherindex);
533  m_signalEnergy[3][idpos] += weight;
534  } else if (abs(motherpdg) == 130 or abs(motherpdg) == 2112) { // neutral hadron
535  int idpos = getIdPosition(4, motherindex);
536  m_signalEnergy[4][idpos] += weight;
537  } else if (abs(motherpdg) == 211 or abs(motherpdg) == 321 or abs(motherpdg) == 2212) { // charged hadron
538  int idpos = getIdPosition(5, motherindex);
539  m_signalEnergy[5][idpos] += weight;
540  } else { // everything else
541  int idpos = getIdPosition(6, motherindex);
542  m_signalEnergy[6][idpos] += weight;
543  }
544 }
Belle2::ECLLocalMaximumFinderModule::getIdPosition
int getIdPosition(const int type, const int id)
Get Id position in the vector.
Definition: ECLLocalMaximumFinderModule.cc:480
Belle2::Unit::cm
static const double cm
Standard units with the value = 1.
Definition: Unit.h:57
Belle2::ECLCalDigit
Class to store calibrated ECLDigits: ECLCalDigits.
Definition: ECLCalDigit.h:38
Belle2::ECLLocalMaximumFinderModule::beginRun
virtual void beginRun() override
Begin.
Definition: ECLLocalMaximumFinderModule.cc:169
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::ECLLocalMaximumFinderModule::addToSignalEnergy
void addToSignalEnergy(int &motherpdg, int &motherindex, int &pi0index, const double &weight)
Add energy to vector.
Definition: ECLLocalMaximumFinderModule.cc:516
Belle2::RelationsInterface::addRelationTo
void addRelationTo(const RelationsInterface< BASE > *object, float weight=1.0, const std::string &namedRelation="") const
Add a relation from this object to another object (with caching).
Definition: RelationsObject.h:144
Belle2::RelationsInterface::getRelationsTo
RelationVector< TO > getRelationsTo(const std::string &name="", const std::string &namedRelation="") const
Get the relations that point from this object to another store array.
Definition: RelationsObject.h:199
Belle2::ECLLocalMaximumFinderModule::terminate
virtual void terminate() override
Terminate (close ROOT files here if you have opened any).
Definition: ECLLocalMaximumFinderModule.cc:358
Belle2::ECLLocalMaximumFinderModule::resetTrainingVariables
void resetTrainingVariables()
Reset Debug Variables.
Definition: ECLLocalMaximumFinderModule.cc:393
Belle2::ECLLocalMaximumFinderModule::event
virtual void event() override
Event.
Definition: ECLLocalMaximumFinderModule.cc:174
Belle2::B2Vector3< double >
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::ECL::ECLGeometryPar::Instance
static ECLGeometryPar * Instance()
Static method to get a reference to the ECLGeometryPar instance.
Definition: ECLGeometryPar.cc:151
Belle2::ECLLocalMaximumFinderModule::getEnteringMother
void getEnteringMother(const MCParticle &particle, int &pdg, int &arrayindex, int &pi0arrayindex)
Get enterging mother of this particle.
Definition: ECLLocalMaximumFinderModule.cc:432
Belle2::Unit::MeV
static const double MeV
[megaelectronvolt]
Definition: Unit.h:124
Belle2::ECL::ECLNeighbours
Class to get the neighbours for a given cell id.
Definition: ECLNeighbours.h:38
Belle2::ECLLocalMaximumFinderModule::isEnteringECL
bool isEnteringECL(const B2Vector3D &vec)
Check if particle is produced outside of the ECL.
Definition: ECLLocalMaximumFinderModule.cc:463
Belle2::ECLLocalMaximumFinderModule::~ECLLocalMaximumFinderModule
virtual ~ECLLocalMaximumFinderModule()
Destructor.
Definition: ECLLocalMaximumFinderModule.cc:80
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::ECLLocalMaximumFinderModule::resetClassifierVariables
void resetClassifierVariables()
Reset Classifier Variables.
Definition: ECLLocalMaximumFinderModule.cc:410
Belle2::ECLLocalMaximumFinderModule::initialize
virtual void initialize() override
Initialize.
Definition: ECLLocalMaximumFinderModule.cc:85
Belle2::ECLLocalMaximumFinderModule
Class to find connected regions.
Definition: ECLLocalMaximumFinderModule.h:50
Belle2::ECLLocalMaximumFinderModule::getMax
void getMax(int &maxtype, int &maxpos)
Get the highest energy deposition particle type.
Definition: ECLLocalMaximumFinderModule.cc:495
Belle2::ECLLocalMaximumFinderModule::makeLocalMaximum
void makeLocalMaximum(const ECLConnectedRegion &aCR, const int cellId, const int lmId)
Make local maximum for a given connected region.
Definition: ECLLocalMaximumFinderModule.cc:375
Belle2::MCParticle
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:43
Belle2::ECLConnectedRegion
Class to store connected regions (CRs)
Definition: ECLConnectedRegion.h:31
Belle2::ECLLocalMaximumFinderModule::endRun
virtual void endRun() override
End run.
Definition: ECLLocalMaximumFinderModule.cc:352