Belle II Software  release-05-02-19
KLMClustersReconstructorModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Kirill Chilikin *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 /* Own header. */
12 #include <klm/modules/KLMClustersReconstructor/KLMClustersReconstructorModule.h>
13 
14 /* KLM headers. */
15 #include <klm/dataobjects/bklm/BKLMElementNumbers.h>
16 #include <klm/dataobjects/eklm/EKLMElementNumbers.h>
17 #include <klm/modules/KLMClustersReconstructor/KLMHit2d.h>
18 
19 /* C++ headers. */
20 #include <algorithm>
21 
22 using namespace Belle2;
23 
24 REG_MODULE(KLMClustersReconstructor)
25 
27  m_PositionMode(c_FirstLayer),
28  m_ClusterMode(c_AnyHit)
29 {
30  setDescription("Unified BKLM/EKLM module for the reconstruction of KLMClusters.");
31  setPropertyFlags(c_ParallelProcessingCertified);
32  addParam("ClusteringAngle", m_ClusteringAngle, "Clustering angle (rad).",
33  0.26);
34  addParam("PositionMode", m_PositionModeString,
35  "Vertex position calculation mode ('FullAverage' or 'FirstLayer').",
36  std::string("FirstLayer"));
37  addParam("ClusterMode", m_ClusterModeString,
38  "Clusterization mode ('AnyHit' or 'FirstHit').",
39  std::string("AnyHit"));
40 }
41 
43 {
44 }
45 
47 {
48  m_KLMClusters.registerInDataStore();
49  m_BKLMHit2ds.isRequired();
50  m_EKLMHit2ds.isRequired();
51  m_KLMClusters.registerRelationTo(m_BKLMHit2ds);
52  m_KLMClusters.registerRelationTo(m_EKLMHit2ds);
53  if (m_PositionModeString == "FullAverage")
55  else if (m_PositionModeString == "FirstLayer")
57  else
58  B2FATAL("Incorrect PositionMode argument.");
59  if (m_ClusterModeString == "AnyHit")
61  else if (m_ClusterModeString == "FirstHit")
63  else
64  B2FATAL("Incorrect ClusterMode argument.");
65 }
66 
68 {
69 }
70 
71 static bool compareDistance(const KLMHit2d& hit1, const KLMHit2d& hit2)
72 {
73  return hit1.getPosition().Mag() < hit2.getPosition().Mag();
74 }
75 
77 {
78  //static double mass = TDatabasePDG::Instance()->GetParticle(130)->Mass();
79  int i, nLayers, innermostLayer, nHits;
82  int* layerHitsBKLM, *layerHitsEKLM;
83  float minTime = -1;
84  double p;//, v;
85  std::vector<KLMHit2d> klmHit2ds, klmClusterHits;
86  std::vector<KLMHit2d>::iterator it, it0, it2;
87  KLMCluster* klmCluster;
88  TVector3 hitPos;
89  layerHitsBKLM = new int[nLayersBKLM];
90  layerHitsEKLM = new int[nLayersEKLM];
91  /* Fill vector of 2d hits. */
92  int nHitsBKLM = m_BKLMHit2ds.getEntries();
93  for (i = 0; i < nHitsBKLM; i++) {
94  if (m_BKLMHit2ds[i]->isOutOfTime())
95  continue;
96  klmHit2ds.push_back(KLMHit2d(m_BKLMHit2ds[i]));
97  }
98  int nHitsEKLM = m_EKLMHit2ds.getEntries();
99  for (i = 0; i < nHitsEKLM; i++) {
100  klmHit2ds.push_back(KLMHit2d(m_EKLMHit2ds[i]));
101  }
102  /* Sort by the distance from center. */
103  sort(klmHit2ds.begin(), klmHit2ds.end(), compareDistance);
104  /* Clustering. */
105  while (klmHit2ds.size() > 0) {
106  klmClusterHits.clear();
107  it = klmHit2ds.begin();
108  klmClusterHits.push_back(*it);
109  it = klmHit2ds.erase(it);
110  while (it != klmHit2ds.end()) {
111  it2 = klmClusterHits.begin();
112  switch (m_ClusterMode) {
113  case c_AnyHit:
114  while (it2 != klmClusterHits.end()) {
115  if (it->getPosition().Angle(it2->getPosition()) <
117  klmClusterHits.push_back(*it);
118  it = klmHit2ds.erase(it);
119  goto clusterFound;
120  } else
121  ++it2;
122  }
123  break;
124  case c_FirstHit:
125  if (it->getPosition().Angle(it2->getPosition()) <
127  klmClusterHits.push_back(*it);
128  it = klmHit2ds.erase(it);
129  goto clusterFound;
130  }
131  break;
132  }
133  ++it;
134 clusterFound:;
135  }
136  hitPos.SetX(0);
137  hitPos.SetY(0);
138  hitPos.SetZ(0);
139  for (i = 0; i < nLayersBKLM; i++)
140  layerHitsBKLM[i] = 0;
141  for (i = 0; i < nLayersEKLM; i++)
142  layerHitsEKLM[i] = 0;
143  /* Find minimal time, fill layer array, find hit position. */
144  it0 = klmClusterHits.begin();
145  nHits = 0;
146  for (it = klmClusterHits.begin(); it != klmClusterHits.end(); ++it) {
147  if ((it->getLayer() == it0->getLayer() &&
148  it->inBKLM() == it0->inBKLM()) || m_PositionMode == c_FullAverage) {
149  hitPos = hitPos + it->getPosition();
150  nHits++;
151  }
152  if (minTime < 0 || it->getTime() < minTime)
153  minTime = it->getTime();
154  if (it->inBKLM())
155  layerHitsBKLM[it->getLayer() - 1]++;
156  else
157  layerHitsEKLM[it->getLayer() - 1]++;
158  }
159  hitPos = hitPos * (1.0 / nHits);
160  /* Find innermost layer. */
161  nLayers = 0;
162  innermostLayer = -1;
163  for (i = 0; i < nLayersBKLM; i++) {
164  if (layerHitsBKLM[i] > 0) {
165  nLayers++;
166  if (innermostLayer < 0)
167  innermostLayer = i + 1;
168  }
169  }
170  for (i = 0; i < nLayersEKLM; i++) {
171  if (layerHitsEKLM[i] > 0) {
172  nLayers++;
173  if (innermostLayer < 0)
174  innermostLayer = i + 1;
175  }
176  }
177  /* Calculate energy. */
178  //if (it0->inBKLM()) {
179  /*
180  * TODO: The constant is from BKLM K0L reconstructor,
181  * it must be recalculated.
182  */
183  p = klmClusterHits.size() * 0.215;
184  /* FIXME: Reimplement time calculation after completion of time calibration.
185  } else {
186  v = hitPos.Mag() / minTime / Const::speedOfLight;
187  if (v < 0.999999)
188  p = mass * v / sqrt(1.0 - v * v);
189  else
190  p = 0;
191  }*/
192  klmCluster = m_KLMClusters.appendNew(
193  hitPos.x(), hitPos.y(), hitPos.z(), minTime, nLayers,
194  innermostLayer, p);
195  for (it = klmClusterHits.begin(); it != klmClusterHits.end(); ++it) {
196  if (it->inBKLM())
197  klmCluster->addRelationTo(it->getBKLMHit2d());
198  else
199  klmCluster->addRelationTo(it->getEKLMHit2d());
200  }
201  }
202  delete[] layerHitsBKLM;
203  delete[] layerHitsEKLM;
204 }
205 
207 {
208 }
209 
211 {
212 }
213 
Belle2::KLMClustersReconstructorModule::m_PositionModeString
std::string m_PositionModeString
Vertex position calculation mode.
Definition: KLMClustersReconstructorModule.h:105
Belle2::KLMClustersReconstructorModule::c_FullAverage
@ c_FullAverage
Full average.
Definition: KLMClustersReconstructorModule.h:50
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::KLMClustersReconstructorModule::m_KLMClusters
StoreArray< KLMCluster > m_KLMClusters
KLM clusters.
Definition: KLMClustersReconstructorModule.h:117
Belle2::KLMClustersReconstructorModule::m_ClusterMode
enum ClusterMode m_ClusterMode
Clusterization mode.
Definition: KLMClustersReconstructorModule.h:114
Belle2::KLMClustersReconstructorModule::c_FirstHit
@ c_FirstHit
Angle from first hit.
Definition: KLMClustersReconstructorModule.h:59
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::KLMClustersReconstructorModule::~KLMClustersReconstructorModule
~KLMClustersReconstructorModule()
Destructor.
Definition: KLMClustersReconstructorModule.cc:42
Belle2::KLMHit2d
Class for simultaneous handling of the BKLM and EKLM 2d hits.
Definition: KLMHit2d.h:34
Belle2::KLMHit2d::getPosition
TVector3 getPosition() const
Get hit position.
Definition: KLMHit2d.cc:34
Belle2::KLMClustersReconstructorModule::beginRun
void beginRun() override
Called when entering a new run.
Definition: KLMClustersReconstructorModule.cc:67
Belle2::KLMClustersReconstructorModule::m_EKLMHit2ds
StoreArray< EKLMHit2d > m_EKLMHit2ds
EKLM 2d hits.
Definition: KLMClustersReconstructorModule.h:123
Belle2::EKLMElementNumbers::getMaximalLayerNumber
static constexpr int getMaximalLayerNumber()
Get maximal layer number.
Definition: EKLMElementNumbers.h:291
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::BKLMElementNumbers::getMaximalLayerNumber
static constexpr int getMaximalLayerNumber()
Get maximal layer number (1-based).
Definition: BKLMElementNumbers.h:251
Belle2::KLMClustersReconstructorModule::m_ClusterModeString
std::string m_ClusterModeString
Clusterization mode.
Definition: KLMClustersReconstructorModule.h:111
Belle2::KLMClustersReconstructorModule::terminate
void terminate() override
This method is called at the end of the event processing.
Definition: KLMClustersReconstructorModule.cc:210
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::KLMClustersReconstructorModule
Module KLMClustersReconstructorModule.
Definition: KLMClustersReconstructorModule.h:44
Belle2::KLMClustersReconstructorModule::m_ClusteringAngle
double m_ClusteringAngle
Clustering angle.
Definition: KLMClustersReconstructorModule.h:102
Belle2::KLMCluster
KLM cluster data.
Definition: KLMCluster.h:38
Belle2::KLMClustersReconstructorModule::event
void event() override
This method is called for each event.
Definition: KLMClustersReconstructorModule.cc:76
Belle2::KLMClustersReconstructorModule::c_AnyHit
@ c_AnyHit
Angle from any hit.
Definition: KLMClustersReconstructorModule.h:58
Belle2::KLMClustersReconstructorModule::initialize
void initialize() override
Initializer.
Definition: KLMClustersReconstructorModule.cc:46
Belle2::KLMClustersReconstructorModule::m_PositionMode
enum PositionMode m_PositionMode
Vertex position calculation mode.
Definition: KLMClustersReconstructorModule.h:108
Belle2::KLMClustersReconstructorModule::endRun
void endRun() override
This method is called if the current run ends.
Definition: KLMClustersReconstructorModule.cc:206
Belle2::KLMClustersReconstructorModule::c_FirstLayer
@ c_FirstLayer
First layer only.
Definition: KLMClustersReconstructorModule.h:51
Belle2::KLMClustersReconstructorModule::m_BKLMHit2ds
StoreArray< BKLMHit2d > m_BKLMHit2ds
BKLM 2d hits.
Definition: KLMClustersReconstructorModule.h:120