Belle II Software  release-05-01-25
SimpleBoxDivisionHoughTree3D.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2015 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Nils Braun, Oliver Frost, Dmitrii Neverov *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #pragma once
11 #include <tracking/trackFindingCDC/hough/trees/BoxDivisionHoughTree.h>
12 #include <tracking/trackFindingCDC/eventdata/hits/CDCRecoHit3D.h>
13 
14 #include <TGraph.h>
15 #include <TF1.h>
16 #include <TCanvas.h>
17 #include <TAxis.h>
18 
19 namespace Belle2 {
24  namespace TrackFindingCDC {
25 
28  template<class AHitPtr, class AInBoxAlgorithm, size_t divisionX, size_t divisionY, size_t divisionZ>
29  class SimpleBoxDivisionHoughTree3D : public
30  BoxDivisionHoughTree<AHitPtr, typename AInBoxAlgorithm::HoughBox, divisionX, divisionY, divisionZ> {
31 
32  private:
34  using Super = BoxDivisionHoughTree<AHitPtr, typename AInBoxAlgorithm::HoughBox, divisionX, divisionY, divisionZ>;
35 
37  using HoughBox = typename AInBoxAlgorithm::HoughBox;
38 
40  template <size_t I>
41  using Width = typename HoughBox::template Width<I>;
42 
43  public:
46  float maximumY,
47  float maximumZ,
48  Width<0> overlapX = 0,
49  Width<1> overlapY = 0,
50  Width<2> overlapZ = 0)
51  : Super(0)
52  , m_maximumX(maximumX)
53  , m_maximumY(maximumY)
54  , m_maximumZ(maximumZ)
55  , m_overlapX(overlapX)
56  , m_overlapY(overlapY)
57  , m_overlapZ(overlapZ)
58  {
59  }
60 
62  void initialize()
63  {
64  Super::template constructArray<0>(-getMaximumX(), getMaximumX(), getOverlapX());
65  Super::template constructArray<1>(-getMaximumY(), getMaximumY(), getOverlapY());
66  Super::template constructArray<2>(-getMaximumZ(), getMaximumZ(), getOverlapZ());
67 
69  }
70 
72  std::vector<std::pair<HoughBox, std::vector<AHitPtr>>>
73  findSingleBest(const Weight& minWeight)
74  {
75  AInBoxAlgorithm inBoxAlgorithm;
76  auto skipLowWeightNode = [minWeight](const typename Super::Node * node) {
77  return not(node->getWeight() >= minWeight);
78  };
79  auto found = this->getTree()->findHeaviestLeafSingle(inBoxAlgorithm, this->getMaxLevel(), skipLowWeightNode);
80 
81  std::vector<std::pair<HoughBox, std::vector<AHitPtr>>> result;
82  if (found) {
83  // Move the found content over. unique_ptr still destroys the left overs.
84  result.push_back(std::move(*found));
85  }
86  return result;
87  }
88 
90  void writeDebugInfoToFile(const std::string& filename __attribute__((unused)))
91  {
92  //do nothing;
93  }
94 
100  void drawDebugPlot(const std::vector<CDCRecoHit3D>& allHits,
101  const std::vector<CDCRecoHit3D>& foundHits,
102  const typename AInBoxAlgorithm::HoughBox& node)
103  {
104  TGraph* allHitsGraph = new TGraph();
105  allHitsGraph->SetLineWidth(2);
106  allHitsGraph->SetLineColor(9);
107 
108  for (const CDCRecoHit3D& recoHit3D : allHits) {
109  const Vector3D& recoPos3D = recoHit3D.getRecoPos3D();
110  const double R = std::sqrt(recoPos3D.x() * recoPos3D.x() + recoPos3D.y() * recoPos3D.y());
111  const double Z = recoPos3D.z();
112  allHitsGraph->SetPoint(allHitsGraph->GetN(), R, Z);
113  }
114 
115  static int nevent(0);
116  TCanvas canv("trackCanvas", "CDC stereo hits in an event", 0, 0, 1600, 1200);
117  canv.cd();
118  allHitsGraph->Draw("APL*");
119  allHitsGraph->GetXaxis()->SetLimits(0, 120);
120  allHitsGraph->GetYaxis()->SetRangeUser(-180, 180);
121 
122  TGraph* foundHitsGraph = new TGraph();
123  foundHitsGraph->SetMarkerStyle(8);
124  foundHitsGraph->SetMarkerColor(2);
125 
126  for (const CDCRecoHit3D& recoHit3D : foundHits) {
127  const Vector3D& recoPos3D = recoHit3D.getRecoPos3D();
128  const double R = std::sqrt(recoPos3D.x() * recoPos3D.x() + recoPos3D.y() * recoPos3D.y());
129  const double Z = recoPos3D.z();
130  foundHitsGraph->SetPoint(foundHitsGraph->GetN(), R, Z);
131  }
132  foundHitsGraph->Draw("P");
133 
134  const double centerX = (AInBoxAlgorithm::BoxAlgorithm::centerX(node));
135  const double deltaX = (AInBoxAlgorithm::BoxAlgorithm::deltaX(node));
136  const double centerY = (AInBoxAlgorithm::BoxAlgorithm::centerY(node));
137  const double centerZ = (AInBoxAlgorithm::BoxAlgorithm::centerZ(node));
138 
139  TF1* candidateL = new TF1("candL", AInBoxAlgorithm::BoxAlgorithm::debugLine(), 0, 120);
140  TF1* candidateH = new TF1("candH", AInBoxAlgorithm::BoxAlgorithm::debugLine(), 0, 120);
141  TF1* candidateMean = new TF1("candMean", AInBoxAlgorithm::BoxAlgorithm::debugLine(), 0, 120);
142 
143  candidateL->SetParameters(centerX - deltaX, centerY, centerZ - 100.0 * deltaX);
144  candidateH->SetParameters(centerX + deltaX, centerY, centerZ + 100.0 * deltaX);
145  candidateMean->SetParameters(centerX, centerY, centerZ);
146 
147  candidateL->SetLineColor(9);
148  candidateH->SetLineColor(41);
149  candidateMean->SetLineColor(2);
150 
151  candidateL->Draw("same");
152  candidateH->Draw("same");
153  candidateMean->Draw("same");
154  canv.SaveAs(Form("CDCRLHits_%i.png", nevent));
155  nevent++;
156  }
157 
159  float getMaximumX() const
160  {
161  return m_maximumX;
162  }
163 
165  float getMaximumY() const
166  {
167  return m_maximumY;
168  }
169 
171  float getMaximumZ() const
172  {
173  return m_maximumZ;
174  }
175 
177  Width<0> getOverlapX() const
178  {
179  return m_overlapX;
180  }
181 
183  Width<1> getOverlapY() const
184  {
185  return m_overlapY;
186  }
187 
189  Width<2> getOverlapZ() const
190  {
191  return m_overlapZ;
192  }
193 
194  private:
196  float m_maximumX = 0;
197 
199  float m_maximumY = 0;
200 
202  float m_maximumZ = 0;
203 
205  Width<0> m_overlapX = 0;
206 
208  Width<1> m_overlapY = 0;
209 
211  Width<2> m_overlapZ = 0;
212  };
213  }
215 }
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::findSingleBest
std::vector< std::pair< HoughBox, std::vector< AHitPtr > > > findSingleBest(const Weight &minWeight)
Find only the leaf with the highest weight (~= number of items)
Definition: SimpleBoxDivisionHoughTree3D.h:81
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::writeDebugInfoToFile
void writeDebugInfoToFile(const std::string &filename __attribute__((unused)))
Write debug information into a ROOT file; not implemented.
Definition: SimpleBoxDivisionHoughTree3D.h:98
Belle2::TrackFindingCDC::CDCRecoHit3D
Class representing a three dimensional reconstructed hit.
Definition: CDCRecoHit3D.h:62
Belle2::TrackFindingCDC::BoxDivisionHoughTree< AHitPointerType, AHitDecisionAlgorithm ::HoughBox, divisionX, divisionY, divisionZ >::Width
typename HoughBox::template Width< I > Width
Type of the width in coordinate I.
Definition: BoxDivisionHoughTree.h:70
Belle2::TrackFindingCDC::WeightedFastHoughTree::findHeaviestLeafSingle
std::unique_ptr< std::pair< ADomain, std::vector< T > > > findHeaviestLeafSingle(AItemInDomainMeasure &weightItemInDomain, int maxLevel, ASkipNodePredicate &skipNode)
Go through all children until the maxLevel is reached and find the leaf with the highest weight.
Definition: WeightedFastHoughTree.h:188
Belle2::TrackFindingCDC::BoxDivisionHoughTree< AHitPtr, AInBoxAlgorithm::HoughBox, divisionX, divisionY, divisionZ >::getTree
HoughTree * getTree() const
Getter for the tree used in the search in the hough plane.
Definition: BoxDivisionHoughTree.h:204
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::getMaximumX
float getMaximumX() const
Return the maximum value in x direction.
Definition: SimpleBoxDivisionHoughTree3D.h:167
Belle2::TrackFindingCDC::Vector3D::x
double x() const
Getter for the x coordinate.
Definition: Vector3D.h:464
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::initialize
void initialize()
Initialize the tree with the given values.
Definition: SimpleBoxDivisionHoughTree3D.h:70
Belle2::TrackFindingCDC::BoxDivisionHoughTree
A fast hough algorithm with rectangular boxes, which are split linearly by a fixed number of division...
Definition: BoxDivisionHoughTree.h:44
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::m_overlapZ
Width< 2 > m_overlapZ
The overlap in Y direction.
Definition: SimpleBoxDivisionHoughTree3D.h:219
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::getOverlapX
Width< 0 > getOverlapX() const
Return the overlap in x direction.
Definition: SimpleBoxDivisionHoughTree3D.h:185
Belle2::TrackFindingCDC::BoxDivisionHoughTree::Node
typename HoughTree::Node Node
Type of the nodes used in the tree for the search.
Definition: BoxDivisionHoughTree.h:73
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::m_maximumY
float m_maximumY
The maximum value in y direction.
Definition: SimpleBoxDivisionHoughTree3D.h:207
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::getMaximumY
float getMaximumY() const
Return the maximum value in y direction.
Definition: SimpleBoxDivisionHoughTree3D.h:173
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::m_maximumX
float m_maximumX
The maximum value in X direction.
Definition: SimpleBoxDivisionHoughTree3D.h:204
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::HoughBox
typename AInBoxAlgorithm::HoughBox HoughBox
The HoughBox we use.
Definition: SimpleBoxDivisionHoughTree3D.h:45
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackFindingCDC::BoxDivisionHoughTree::initialize
void initialize()
Initialise the algorithm by constructing the hough tree from the parameters.
Definition: BoxDivisionHoughTree.h:174
Belle2::TrackFindingCDC::BoxDivisionHoughTree< AHitPtr, AInBoxAlgorithm::HoughBox, divisionX, divisionY, divisionZ >::getMaxLevel
int getMaxLevel() const
Getter for the currently set maximal level.
Definition: BoxDivisionHoughTree.h:210
Belle2::TrackFindingCDC::Vector3D
A three dimensional vector.
Definition: Vector3D.h:34
Belle2::TrackFindingCDC::Vector3D::y
double y() const
Getter for the y coordinate.
Definition: Vector3D.h:476
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::getOverlapY
Width< 1 > getOverlapY() const
Return the overlap in y direction.
Definition: SimpleBoxDivisionHoughTree3D.h:191
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::Width
typename HoughBox::template Width< I > Width
Type of the width in coordinate I.
Definition: SimpleBoxDivisionHoughTree3D.h:49
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::m_maximumZ
float m_maximumZ
The maximum value in z direction.
Definition: SimpleBoxDivisionHoughTree3D.h:210
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::Super
BoxDivisionHoughTree< AHitPtr, typename AInBoxAlgorithm::HoughBox, divisionX, divisionY, divisionZ > Super
The Super class.
Definition: SimpleBoxDivisionHoughTree3D.h:42
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::getMaximumZ
float getMaximumZ() const
Return the maximum value in Z direction.
Definition: SimpleBoxDivisionHoughTree3D.h:179
Belle2::TrackFindingCDC::Vector3D::z
double z() const
Getter for the z coordinate.
Definition: Vector3D.h:488
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::getOverlapZ
Width< 2 > getOverlapZ() const
Return the overlap in y direction.
Definition: SimpleBoxDivisionHoughTree3D.h:197
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::SimpleBoxDivisionHoughTree3D
SimpleBoxDivisionHoughTree3D(float maximumX, float maximumY, float maximumZ, Width< 0 > overlapX=0, Width< 1 > overlapY=0, Width< 2 > overlapZ=0)
Constructor using the given maximal level.
Definition: SimpleBoxDivisionHoughTree3D.h:53
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::m_overlapY
Width< 1 > m_overlapY
The overlap in Y direction.
Definition: SimpleBoxDivisionHoughTree3D.h:216
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::drawDebugPlot
void drawDebugPlot(const std::vector< CDCRecoHit3D > &allHits, const std::vector< CDCRecoHit3D > &foundHits, const typename AInBoxAlgorithm::HoughBox &node)
Draws found hits and node boundaries FIXME this is a copy-paste from DebugableSimpleBoxDivisionHoughT...
Definition: SimpleBoxDivisionHoughTree3D.h:108
Belle2::TrackFindingCDC::SimpleBoxDivisionHoughTree3D::m_overlapX
Width< 0 > m_overlapX
The overlap in X direction.
Definition: SimpleBoxDivisionHoughTree3D.h:213