Belle II Software  release-05-01-25
ECLSplitterN2Module.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 - Belle II Collaboration *
4  * *
5  * See .h file for a description. *
6  * *
7  * Author: The Belle II Collaboration *
8  * Contributors: Torben Ferber (ferber@physics.ubc.ca) (TF) *
9  * *
10  * This software is provided "as is" without any warranty. *
11  **************************************************************************/
12 
13 // THIS MODULE
14 #include <ecl/modules/eclSplitterN2/ECLSplitterN2Module.h>
15 
16 // FRAMEWORK
17 #include <framework/datastore/RelationArray.h>
18 #include <framework/logging/Logger.h>
19 
20 // ECL
21 #include <ecl/utility/Position.h>
22 #include <ecl/dataobjects/ECLCalDigit.h>
23 #include <ecl/dataobjects/ECLConnectedRegion.h>
24 #include <ecl/dataobjects/ECLLocalMaximum.h>
25 #include <ecl/dataobjects/ECLShower.h>
26 
27 // MDST
28 #include <mdst/dataobjects/ECLCluster.h>
29 
30 // OTHER
31 #include <string>
32 
33 // NAMESPACES
34 using namespace Belle2;
35 using namespace ECL;
36 
37 //-----------------------------------------------------------------
38 // Register the Module(s)
39 //-----------------------------------------------------------------
40 REG_MODULE(ECLSplitterN2)
41 REG_MODULE(ECLSplitterN2PureCsI)
42 
43 //-----------------------------------------------------------------
44 // Implementation
45 //-----------------------------------------------------------------
46 
47 ECLSplitterN2Module::ECLSplitterN2Module() : Module(), m_eclCalDigits(eclCalDigitArrayName()),
48  m_eclConnectedRegions(eclConnectedRegionArrayName()),
49  m_eclLocalMaximums(eclLocalMaximumArrayName()),
50  m_eclShowers(eclShowerArrayName())
51 {
52  // Set description.
53  setDescription("ECLSplitterN2Module: Baseline reconstruction splitter code for the neutral hadron hypothesis.");
54 
55  // Set module parameters.
56  addParam("positionMethod", m_positionMethod, "Position determination method.", std::string("lilo"));
57  addParam("liloParameterA", m_liloParameterA, "Position determination linear-log. parameter A.", 4.0);
58  addParam("liloParameterB", m_liloParameterB, "Position determination linear-log. parameter B.", 0.0);
59  addParam("liloParameterC", m_liloParameterC, "Position determination linear-log. parameter C.", 0.0);
60 
61  // Set parallel processing flag.
62  setPropertyFlags(c_ParallelProcessingCertified);
63 }
64 
66 {
67  // do not delete objects here, do it in terminate()!
68 }
69 
71 {
72  // Check user input.
73  m_liloParameters.resize(3);
74  m_liloParameters.at(0) = m_liloParameterA;
75  m_liloParameters.at(1) = m_liloParameterB;
76  m_liloParameters.at(2) = m_liloParameterC;
77 
78  // ECL dataobjects.
79  m_eclCalDigits.registerInDataStore(eclCalDigitArrayName());
80  m_eclConnectedRegions.registerInDataStore(eclConnectedRegionArrayName());
81  m_eclShowers.registerInDataStore(eclShowerArrayName());
82 
83  // Register relations (we probably dont need all, but keep them for now for debugging).
84  m_eclShowers.registerRelationTo(m_eclConnectedRegions);
85  m_eclShowers.registerRelationTo(m_eclCalDigits);
86  m_eclShowers.registerRelationTo(m_eclLocalMaximums);
87 
88 }
89 
91 {
92  ;
93 }
94 
96 {
97  B2DEBUG(175, "ECLCRSplitterN2Module::event()");
98 
99  // Loop over all connected regions (CR_.
100  for (auto& aCR : m_eclConnectedRegions) {
101  unsigned int iShower = 1;
102 
103  const auto aECLShower = m_eclShowers.appendNew();
104 
105  // Add relation to the CR.
106  aECLShower->addRelationTo(&aCR);
107 
108  // Loop over all local maximums (LM).
109  for (auto& aLM : aCR.getRelationsWith<ECLLocalMaximum>(eclLocalMaximumArrayName())) {
110  // Add relation to the CR.
111  aECLShower->addRelationTo(&aLM);
112  }
113 
114  // Prepare shower variables.
115  double highestEnergy = 0.0;
116  double highestEnergyTime = 0.;
117  double highestEnergyTimeResolution = 0.;
118  double weightSum = 0.0;
119  double energySum = 0.0;
120  unsigned int highestEnergyID = 0;
121  std::vector< ECLCalDigit > digits;
122  std::vector< double > weights;
123 
124  // Loop over all digits that are related to the CR, they can be weighted (in the future?).
125  auto relatedDigitsPairs = aCR.getRelationsTo<ECLCalDigit>(eclCalDigitArrayName());
126  for (unsigned int iRel = 0; iRel < relatedDigitsPairs.size(); iRel++) {
127  const auto aECLCalDigit = relatedDigitsPairs.object(iRel);
128  const auto weight = relatedDigitsPairs.weight(iRel);
129 
130  // Add Relation to ECLCalDigits.
131  aECLShower->addRelationTo(aECLCalDigit, weight);
132 
133  // Find highest energetic crystal, its time, and its time resolution. This is not neceessarily the LM!
134  const double energyDigit = aECLCalDigit->getEnergy();
135  if (energyDigit > highestEnergy) {
136  highestEnergy = energyDigit * weight;
137  highestEnergyID = aECLCalDigit->getCellId();
138  highestEnergyTime = aECLCalDigit->getTime();
139  highestEnergyTimeResolution = aECLCalDigit->getTimeResolution();
140  }
141 
142  digits.push_back(*aECLCalDigit);
143  weights.push_back(weight);
144 
145  weightSum += weight;
146  energySum += energyDigit * weight;
147 
148  }
149 
150  const TVector3& showerposition = Belle2::ECL::computePositionLiLo(digits, weights, m_liloParameters);
151  aECLShower->setTheta(showerposition.Theta());
152  aECLShower->setPhi(showerposition.Phi());
153  aECLShower->setR(showerposition.Mag());
154 
155  aECLShower->setEnergy(energySum);
156  aECLShower->setEnergyRaw(energySum);
157  aECLShower->setEnergyHighestCrystal(highestEnergy);
158  aECLShower->setCentralCellId(highestEnergyID);
159  aECLShower->setTime(highestEnergyTime);
160  aECLShower->setDeltaTime99(highestEnergyTimeResolution);
161  aECLShower->setNumberOfCrystals(weightSum);
162 
163  aECLShower->setShowerId(iShower);
164  aECLShower->setHypothesisId(Belle2::ECLShower::c_neutralHadron);
165  aECLShower->setConnectedRegionId(aCR.getCRId());
166 
167  B2DEBUG(175, "N2 shower " << iShower);
168  B2DEBUG(175, " theta = " << aECLShower->getTheta());
169  B2DEBUG(175, " phi = " << aECLShower->getPhi());
170  B2DEBUG(175, " R = " << aECLShower->getR());
171  B2DEBUG(175, " energy = " << aECLShower->getEnergy());
172  B2DEBUG(175, " time = " << aECLShower->getTime());
173  B2DEBUG(175, " time resolution = " << aECLShower->getDeltaTime99());
174 
175  } // end auto& aCR
176 
177 }
178 
179 
181 {
182  ;
183 }
184 
185 
187 {
188  ;
189 }
Belle2::ECLCalDigit
Class to store calibrated ECLDigits: ECLCalDigits.
Definition: ECLCalDigit.h:38
Belle2::ECLSplitterN2Module::initialize
virtual void initialize() override
Initialize.
Definition: ECLSplitterN2Module.cc:70
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::ECLSplitterN2Module::endRun
virtual void endRun() override
End run.
Definition: ECLSplitterN2Module.cc:180
Belle2::ECLSplitterN2Module
Class to perform the shower correction.
Definition: ECLSplitterN2Module.h:45
Belle2::ECLShower::c_neutralHadron
@ c_neutralHadron
CR is reconstructed as a neutral hadron (N2)
Definition: ECLShower.h:56
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::ECLSplitterN2Module::beginRun
virtual void beginRun() override
Begin run.
Definition: ECLSplitterN2Module.cc:90
Belle2::ECLSplitterN2Module::event
virtual void event() override
Event.
Definition: ECLSplitterN2Module.cc:95
Belle2::ECLLocalMaximum
Class to store local maxima (LM)
Definition: ECLLocalMaximum.h:41
Belle2::ECLSplitterN2Module::terminate
virtual void terminate() override
Terminate.
Definition: ECLSplitterN2Module.cc:186
Belle2::ECLSplitterN2Module::~ECLSplitterN2Module
~ECLSplitterN2Module()
Destructor.
Definition: ECLSplitterN2Module.cc:65