Belle II Software  release-05-02-19
TOPalign.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Marko Staric *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #pragma once
12 
13 #include <top/reconstruction/TOPtrack.h>
14 #include <framework/gearbox/Const.h>
15 #include <framework/logging/Logger.h>
16 #include <vector>
17 #include <string>
18 #include <algorithm>
19 
20 namespace Belle2 {
25  namespace TOP {
26 
30  class TOPalign {
31  public:
32 
36  TOPalign();
37 
41  static void clearData();
42 
51  static int addData(int moduleID, int pixelID, double time, double timeError);
52 
58  static void setPhotonYields(double bkgPerModule, double scaleN0 = 1);
59 
64  void setModuleID(int moduleID) {m_moduleID = moduleID;}
65 
73  void setSteps(double position, double angle, double time, double refind);
74 
80  void setGrid(int NP, int NC)
81  {
82  m_NP = NP;
83  m_NC = NC;
84  }
85 
91  void setParameters(const std::vector<double>& parInit)
92  {
93  for (size_t i = 0; i < std::min(parInit.size(), m_parInit.size()); i++) {
94  m_parInit[i] = parInit[i];
95  }
96  m_par = m_parInit;
97  }
98 
103  void fixParameter(const std::string& name)
104  {
105  for (unsigned i = 0; i < m_parNames.size(); i++) {
106  if (name == m_parNames[i]) {
107  m_fixed[i] = true;
108  return;
109  }
110  }
111  B2ERROR("TOPalign::fixParameter: invalid parameter name '" << name << "'");
112  }
113 
118  void unfixParameter(const std::string& name)
119  {
120  for (unsigned i = 0; i < m_parNames.size(); i++) {
121  if (name == m_parNames[i]) {
122  m_fixed[i] = false;
123  return;
124  }
125  }
126  B2ERROR("TOPalign::unfixParameter: invalid parameter name '" << name << "'");
127  }
128 
132  void unfixAll()
133  {
134  for (unsigned i = 0; i < m_fixed.size(); i++) m_fixed[i] = false;
135  }
136 
143  int iterate(const TOPtrack& track, const Const::ChargedStable& hypothesis);
144 
148  void reset();
149 
154  int getModuleID() const {return m_moduleID;}
155 
161  const std::vector<float>& getParameters() const {return m_par;}
162 
167  const std::vector<std::string>& getParameterNames() const {return m_parNames;}
168 
174  std::vector<float> getErrors() const;
175 
180  const std::vector<float>& getErrorMatrix() const {return m_COV;}
181 
186  int getNumTracks() const {return m_numTracks;}
187 
192  int getNumUsedTracks() const {return m_numUsedTracks;}
193 
198  bool isValid() const {return m_valid;}
199 
204  int getNumOfPhotons() const {return m_numPhotons;}
205 
206  private:
207 
208  int m_moduleID = 0;
209  int m_opt = 0;
210  int m_NP = 0;
211  int m_NC = 0;
213  std::vector<std::string> m_parNames;
214  std::vector<float> m_parInit;
215  std::vector<float> m_par;
216  std::vector<float> m_steps;
217  std::vector<float> m_maxDpar;
218  std::vector<bool> m_fixed;
219  std::vector<float> m_COV;
220  int m_numTracks = 0;
221  int m_numUsedTracks = 0;
222  bool m_valid = false;
223  int m_numPhotons = 0;
225  std::vector<double> m_U;
227  };
228 
229 
230  } // TOP namespace
232 } // Belle2 namespace
233 
Belle2::TOP::TOPalign::getErrorMatrix
const std::vector< float > & getErrorMatrix() const
Returns error matrix of alignment parameters.
Definition: TOPalign.h:188
Belle2::TOP::TOPalign::m_fixed
std::vector< bool > m_fixed
true if parameter is fixed
Definition: TOPalign.h:226
Belle2::TOP::TOPalign::setPhotonYields
static void setPhotonYields(double bkgPerModule, double scaleN0=1)
Sets expected photon yields.
Definition: TOPalign.cc:106
Belle2::TOP::TOPalign::getParameters
const std::vector< float > & getParameters() const
Returns alignment parameters.
Definition: TOPalign.h:169
Belle2::TOP::TOPtrack
Class to hold reconstructed track, interface to fortran.
Definition: TOPtrack.h:42
Belle2::TOP::TOPalign::setSteps
void setSteps(double position, double angle, double time, double refind)
Sets steps for numerical calculation of derivatives.
Definition: TOPalign.cc:59
Belle2::TOP::TOPalign::getNumTracks
int getNumTracks() const
Returns track counter.
Definition: TOPalign.h:194
Belle2::TOP::TOPalign::setGrid
void setGrid(int NP, int NC)
Sets grid for averaging of time-of-propagation in analytic PDF.
Definition: TOPalign.h:88
Belle2::TOP::TOPalign::m_parInit
std::vector< float > m_parInit
initial parameter values
Definition: TOPalign.h:222
Belle2::TOP::TOPalign::getNumUsedTracks
int getNumUsedTracks() const
Returns number of tracks used in current result.
Definition: TOPalign.h:200
Belle2::TOP::TOPalign::m_NC
int m_NC
grid for averaging: number of Cerenkov angles
Definition: TOPalign.h:219
Belle2::TOP::TOPalign::m_numUsedTracks
int m_numUsedTracks
number of tracks used
Definition: TOPalign.h:229
Belle2::TOP::TOPalign::m_numPhotons
int m_numPhotons
number of photons used for log likelihood calculation
Definition: TOPalign.h:231
Belle2::TOP::TOPalign::clearData
static void clearData()
Clear data list.
Definition: TOPalign.cc:73
Belle2::TOP::TOPalign::TOPalign
TOPalign()
Constructor.
Definition: TOPalign.cc:38
Belle2::TOP::TOPalign::setParameters
void setParameters(const std::vector< double > &parInit)
Sets initial values of parameters (overwrites current parameters!) Order is: translations in x,...
Definition: TOPalign.h:99
Belle2::TOP::TOPalign::m_numTracks
int m_numTracks
track counter
Definition: TOPalign.h:228
Belle2::TOP::TOPalign::m_maxDpar
std::vector< float > m_maxDpar
maximal parameter changes in one iteration
Definition: TOPalign.h:225
Belle2::TOP::TOPalign::m_COV
std::vector< float > m_COV
covariance matrix
Definition: TOPalign.h:227
Belle2::TOP::TOPalign::getErrors
std::vector< float > getErrors() const
Returns errors on alignment parameters.
Definition: TOPalign.cc:183
Belle2::TOP::TOPalign::fixParameter
void fixParameter(const std::string &name)
Fixes parameter with its name given as argument.
Definition: TOPalign.h:111
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TOP::TOPalign::addData
static int addData(int moduleID, int pixelID, double time, double timeError)
Add data.
Definition: TOPalign.cc:78
Belle2::TOP::TOPalign::getModuleID
int getModuleID() const
Returns module ID.
Definition: TOPalign.h:162
Belle2::TOP::TOPalign::m_NP
int m_NP
grid for averaging: number of emission points along track
Definition: TOPalign.h:218
Belle2::TOP::TOPalign::unfixParameter
void unfixParameter(const std::string &name)
Unfixes parameter with its name given as argument.
Definition: TOPalign.h:126
Belle2::TOP::TOPalign::m_parNames
std::vector< std::string > m_parNames
parameter names
Definition: TOPalign.h:221
Belle2::Const::ChargedStable
Provides a type-safe way to pass members of the chargedStableSet set.
Definition: Const.h:465
Belle2::TOP::TOPalign::m_moduleID
int m_moduleID
module ID
Definition: TOPalign.h:216
Belle2::TOP::TOPalign::reset
void reset()
Reset the object.
Definition: TOPalign.cc:173
Belle2::TOP::TOPalign::iterate
int iterate(const TOPtrack &track, const Const::ChargedStable &hypothesis)
Run a single iteration.
Definition: TOPalign.cc:114
Belle2::TOP::TOPalign::getNumOfPhotons
int getNumOfPhotons() const
Returns number of photons used for log likelihood calculation.
Definition: TOPalign.h:212
Belle2::TOP::TOPalign::unfixAll
void unfixAll()
Unfixes all parameters.
Definition: TOPalign.h:140
Belle2::TOP::TOPalign::setModuleID
void setModuleID(int moduleID)
Sets module ID.
Definition: TOPalign.h:72
Belle2::TOP::TOPalign::isValid
bool isValid() const
Checks if the results are valid.
Definition: TOPalign.h:206
Belle2::TOP::TOPalign::m_opt
int m_opt
PDF option (=rough)
Definition: TOPalign.h:217
Belle2::TOP::TOPalign::m_steps
std::vector< float > m_steps
step sizes
Definition: TOPalign.h:224
Belle2::TOP::TOPalign::m_par
std::vector< float > m_par
current parameter values
Definition: TOPalign.h:223
Belle2::TOP::TOPalign::getParameterNames
const std::vector< std::string > & getParameterNames() const
Returns alignment parameter names.
Definition: TOPalign.h:175
Belle2::TOP::TOPalign::m_U
std::vector< double > m_U
matrix (neg.
Definition: TOPalign.h:233
Belle2::TOP::TOPalign::m_valid
bool m_valid
validity of results
Definition: TOPalign.h:230