Belle II Software  release-05-01-25
CDCCalibrationCollector.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: CDC Group *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include "cdc/modules/cdcCalibrationCollector/CDCCalibrationCollector.h"
12 #include <cdc/translators/RealisticTDCCountTranslator.h>
13 #include <framework/datastore/StoreObjPtr.h>
14 #include <framework/datastore/StoreArray.h>
15 #include <framework/datastore/RelationArray.h>
16 
17 #include <mdst/dataobjects/TrackFitResult.h>
18 #include <mdst/dataobjects/Track.h>
19 #include <tracking/dataobjects/RecoTrack.h>
20 #include <tracking/trackFindingCDC/rootification/StoreWrappedObjPtr.h>
21 #include <tracking/trackFindingCDC/eventdata/hits/CDCWireHit.h>
22 #include <tracking/trackFindingCDC/eventdata/tracks/CDCTrack.h>
23 #include <tracking/trackFindingCDC/topology/CDCWireTopology.h>
24 
25 #include <genfit/TrackPoint.h>
26 #include <genfit/KalmanFitterInfo.h>
27 #include <genfit/MeasurementOnPlane.h>
28 #include <genfit/MeasuredStateOnPlane.h>
29 
30 #include <Math/ProbFuncMathCore.h>
31 
32 #include <cdc/dataobjects/WireID.h>
33 #include <cdc/geometry/CDCGeometryPar.h>
34 
35 #include <TH1F.h>
36 
37 //using namespace std;
38 using namespace Belle2;
39 using namespace CDC;
40 using namespace genfit;
41 using namespace TrackFindingCDC;
42 
43 
44 REG_MODULE(CDCCalibrationCollector)
45 
46 
48 {
49  setDescription("Collector module for cdc calibration");
50  setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
51  addParam("recoTracksColName", m_recoTrackArrayName, "Name of collection hold genfit::Track", std::string(""));
52  addParam("bField", m_bField, "If true -> #Params ==5 else #params ==4 for calculate P-Val", false);
53  addParam("calExpectedDriftTime", m_calExpectedDriftTime, "if true module will calculate expected drift time, it take a time",
54  true);
55  addParam("storeTrackParams", m_storeTrackParams, "Store Track Parameter or not, it will be multicount for each hit", false);
56  addParam("eventT0Extraction", m_eventT0Extraction, "use event t0 extract t0 or not", true);
57  addParam("minimumPt", m_minimumPt, "Tracks with tranverse momentum smaller than this value will not used", 0.15);
58  addParam("isCosmic", m_isCosmic, "True when we process cosmic events, else False (collision)", m_isCosmic);
59  addParam("effStudy", m_effStudy, "When true module collects info only necessary for wire eff study", false);
60 }
61 
63 {
64 }
65 
67 {
68  StoreArray<Belle2::Track> storeTrack(m_trackArrayName);
69  StoreArray<RecoTrack> recoTracks(m_recoTrackArrayName);
70  StoreArray<Belle2::TrackFitResult> storeTrackFitResults(m_trackFitResultArrayName);
71  StoreArray<Belle2::CDCHit> cdcHits(m_cdcHitArrayName);
72  StoreWrappedObjPtr<std::vector<CDCTrack>> cdcTracks(m_cdcTrackVectorName);
73  RelationArray relRecoTrackTrack(recoTracks, storeTrack, m_relRecoTrackTrackName);
74  //Store names to speed up creation later
75  m_recoTrackArrayName = recoTracks.getName();
76  m_trackFitResultArrayName = storeTrackFitResults.getName();
77  m_relRecoTrackTrackName = relRecoTrackTrack.getName();
78 
79  if (!m_effStudy) { // by default collects calibration data
80  auto m_tree = new TTree(m_treeName.c_str(), "tree for cdc calibration");
81  m_tree->Branch<Float_t>("x_mea", &x_mea);
82  m_tree->Branch<Float_t>("x_u", &x_u);
83  m_tree->Branch<Float_t>("x_b", &x_b);
84  m_tree->Branch<Float_t>("alpha", &alpha);
85  m_tree->Branch<Float_t>("theta", &theta);
86  m_tree->Branch<Float_t>("t", &t);
87  m_tree->Branch<UShort_t>("adc", &adc);
88  // m_tree->Branch<int>("boardID", &boardID);
89  m_tree->Branch<UChar_t>("lay", &lay);
90  m_tree->Branch<Float_t>("weight", &weight);
91  m_tree->Branch<UShort_t>("IWire", &IWire);
92  m_tree->Branch<Float_t>("Pval", &Pval);
93  m_tree->Branch<Float_t>("ndf", &ndf);
94  if (m_storeTrackParams) {
95  m_tree->Branch<Float_t>("d0", &d0);
96  m_tree->Branch<Float_t>("z0", &z0);
97  m_tree->Branch<Float_t>("phi0", &phi0);
98  m_tree->Branch<Float_t>("tanL", &tanL);
99  m_tree->Branch<Float_t>("omega", &omega);
100  }
101 
102  if (m_calExpectedDriftTime) { // expected drift time, calculated form xfit
103  m_tree->Branch<Float_t>("t_fit", &t_fit);
104  }
105 
106  registerObject<TTree>("tree", m_tree);
107  }
108  if (m_effStudy) { //if m_effStudy is changed to true prepares to only run wire efficiency study
109  auto m_efftree = new TTree(m_effTreeName.c_str(), "tree for wire efficiency");
110  m_efftree->Branch<unsigned short>("layerID", &layerID);
111  m_efftree->Branch<unsigned short>("wireID", &wireID);
112  m_efftree->Branch<float>("z", &z);
113  m_efftree->Branch<bool>("isFound", &isFound);
114 
115  registerObject<TTree>("efftree", m_efftree);
116  }
117 
118  auto m_hNDF = new TH1F("hNDF", "NDF of fitted track;NDF;Tracks", 71, -1, 70);
119  auto m_hPval = new TH1F("hPval", "p-values of tracks;pVal;Tracks", 1000, 0, 1);
120  auto m_hEventT0 = new TH1F("hEventT0", "Event T0", 1000, -100, 100);
121  auto m_hNTracks = new TH1F("hNTracks", "Number of tracks", 50, 0, 10);
122  auto m_hOccupancy = new TH1F("hOccupancy", "occupancy", 100, 0, 1.0);
123 
124  registerObject<TH1F>("hNDF", m_hNDF);
125  registerObject<TH1F>("hPval", m_hPval);
126  registerObject<TH1F>("hEventT0", m_hEventT0);
127  registerObject<TH1F>("hNTracks", m_hNTracks);
128  registerObject<TH1F>("hOccupancy", m_hOccupancy);
129 }
130 
132 {
133  const StoreArray<Belle2::Track> storeTrack(m_trackArrayName);
134  const StoreArray<Belle2::TrackFitResult> storeTrackFitResults(m_trackFitResultArrayName);
135  const StoreArray<Belle2::CDCHit> cdcHits(m_cdcHitArrayName);
136  const StoreWrappedObjPtr<std::vector<CDCTrack>> cdcTracks(m_cdcTrackVectorName);
137  const StoreArray<Belle2::RecoTrack> recoTracks(m_recoTrackArrayName);
138  const RelationArray relTrackTrack(recoTracks, storeTrack, m_relRecoTrackTrackName);
139 
140  /* CDCHit distribution */
141  // make evt t0 incase we dont use evt t0
142  evtT0 = 0;
143 
144  //reject events don't have eventT0
145  if (m_eventT0Extraction) {
146  // event with is fail to extract t0 will be exclude from analysis
147  if (m_eventTimeStoreObject.isValid() && m_eventTimeStoreObject->hasEventT0()) {
148  evtT0 = m_eventTimeStoreObject->getEventT0();
149  getObjectPtr<TH1F>("hEventT0")->Fill(evtT0);
150  } else {
151  return;
152  }
153  }
154  // Collects the WireID and Layer of every hit in this event
155  // Used in wire efficiency building
156  std::vector<unsigned short> wiresInCDCTrack;
157 
158  for (CDCTrack& cdcTrack : *cdcTracks) {
159  for (CDCRecoHit3D& cdcHit : cdcTrack) {
160  unsigned short eWireID = cdcHit.getWire().getEWire();
161  wiresInCDCTrack.push_back(eWireID);
162  }
163  }
164  // WireID collection finished
165 
166  const int nTr = recoTracks.getEntries();
167 
168  // Skip events which have number of charged tracks <= 1.
169  int nCTracks = 0;
170  for (int i = 0; i < nTr; ++i) {
171  RecoTrack* track = recoTracks[i];
172  if (track->getDirtyFlag()) continue;
173  if (!track->getTrackFitStatus()->isFitted()) continue;
174  const genfit::FitStatus* fs = track->getTrackFitStatus();
175  if (!fs || !fs->isFitConverged()) continue; //not fully convergence
176 
177  const Belle2::Track* b2track = track->getRelatedFrom<Belle2::Track>();
178  if (!b2track) {B2DEBUG(99, "No relation found"); continue;}
179  const Belle2::TrackFitResult* fitresult = b2track->getTrackFitResult(Const::muon);
180  if (!fitresult) {
181  B2WARNING("track was fitted but Relation not found");
182  continue;
183  }
184 
185  short charge = fitresult->getChargeSign();
186  if (fabs(charge) > 0) {
187  nCTracks++;
188  }
189  }
190  if (nCTracks <= 1) {
191  return ;
192  } else {
193  getObjectPtr<TH1F>("hNTracks")->Fill(nCTracks);
194  }
195 
196 
197  const int nHits = cdcHits.getEntries();
198  const int nWires = 14336;
199  float oc = static_cast<float>(nHits) / static_cast<float>(nWires);
200  getObjectPtr<TH1F>("hOccupancy")->Fill(oc);
201 
202  for (int i = 0; i < nTr; ++i) {
203  RecoTrack* track = recoTracks[i];
204  if (track->getDirtyFlag()) continue;
205  const genfit::FitStatus* fs = track->getTrackFitStatus();
206  if (!fs || !fs->isFitted()) continue;
207  if (!fs->isFitConverged()) continue; //not fully convergence
208 
209  const Belle2::Track* b2track = track->getRelatedFrom<Belle2::Track>();
210  if (!b2track) {B2DEBUG(99, "No relation found"); continue;}
211  const Belle2::TrackFitResult* fitresult = b2track->getTrackFitResult(Const::muon);
212 
213  if (!fitresult) {
214  B2WARNING("track was fitted but Relation not found");
215  continue;
216  }
217  ndf = fs->getNdf();
218  if (!m_bField) {
219  ndf += 1;
220  }
221 
222  getObjectPtr<TH1F>("hPval")->Fill(Pval);
223  getObjectPtr<TH1F>("hNDF")->Fill(ndf);
224  B2DEBUG(99, "ndf = " << ndf);
225  B2DEBUG(99, "Pval = " << Pval);
226 
227  if (ndf < 15) continue; //hard code,
228  double Chi2 = fs->getChi2();
229  Pval = std::max(0., ROOT::Math::chisquared_cdf_c(Chi2, ndf));
230  //store track parameters
231 
232  d0 = fitresult->getD0();
233  z0 = fitresult->getZ0();
234  tanL = fitresult->getTanLambda();
235  omega = fitresult->getOmega();
236  phi0 = fitresult->getPhi0() * 180 / M_PI;
237 
238  // Rejection of suspicious cosmic tracks.
239  // phi0 of cosmic track must be negative in our definition!
240  if (m_isCosmic == true && phi0 > 0.0) continue;
241 
242  //cut at Pt
243  if (fitresult->getTransverseMomentum() < m_minimumPt) continue;
244  if (!m_effStudy) { // all harvest to fill the tree if collecting calibration info
245  try {
246  harvest(track);
247  } catch (...) {
248  }
249  }
250  if (m_effStudy) { // call buildEfficiencies for efficiency study
251  // Request tracks coming from IP
252  if (fitresult->getD0() > 2 || fitresult->getZ0() > 5) continue;
253  const Helix helixFit = fitresult->getHelix();
254  buildEfficiencies(wiresInCDCTrack, helixFit);
255  }
256  }
257 
258 }
259 
261 {
262 }
263 
265 {
266  B2DEBUG(99, "start collect hit");
267  static CDCGeometryPar& cdcgeo = CDCGeometryPar::Instance();
269 
270  for (const RecoHitInformation::UsedCDCHit* hit : track->getCDCHitList()) {
271  const genfit::TrackPoint* tp = track->getCreatedTrackPoint(track->getRecoHitInformation(hit));
272  if (!tp) continue;
273  lay = hit->getICLayer();
274  IWire = hit->getIWire();
275  adc = hit->getADCCount();
276  unsigned short tdc = hit->getTDCCount();
277  WireID wireid(lay, IWire);
279  if (!kfi) {B2DEBUG(199, "No Fitter Info: Layer " << hit->getICLayer()); continue;}
280  for (unsigned int iMeas = 0; iMeas < kfi->getNumMeasurements(); ++iMeas) {
281  if ((kfi->getWeights().at(iMeas)) > 0.5) {
282  // int boardID = cdcgeo.getBoardID(WireID(lay,IWire));
283  const genfit::MeasuredStateOnPlane& mop = kfi->getFittedState();
284  const TVector3 pocaOnWire = mop.getPlane()->getO();//Local wire position
285  const TVector3 pocaOnTrack = mop.getPlane()->getU();//residual direction
286  const TVector3 pocaMom = mop.getMom();
287  alpha = cdcgeo.getAlpha(pocaOnWire, pocaMom) ;
288  theta = cdcgeo.getTheta(pocaMom);
289  x_mea = kfi->getMeasurementOnPlane(iMeas)->getState()(0);
290  x_b = kfi->getFittedState(true).getState()(3);// x fit biased
291  x_u = kfi->getFittedState(false).getState()(3);//x fit unbiased
292  weight = kfi->getWeights().at(iMeas);
293 
294  int lr;
295  if (x_u > 0) lr = 1;
296  else lr = 0;
297 
298  //Convert to outgoing
299  if (fabs(alpha) > M_PI / 2) {
300  x_b *= -1;
301  x_u *= -1;
302  }
303  x_mea = copysign(x_mea, x_b);
304  lr = cdcgeo.getOutgoingLR(lr, alpha);
305  theta = cdcgeo.getOutgoingTheta(alpha, theta);
306  alpha = cdcgeo.getOutgoingAlpha(alpha);
307 
308  B2DEBUG(99, "x_unbiased " << x_u << " |left_right " << lr);
309  if (m_calExpectedDriftTime) { t_fit = cdcgeo.getDriftTime(abs(x_u), lay, lr, alpha , theta);}
310  alpha *= 180 / M_PI;
311  theta *= 180 / M_PI;
312  //estimate drift time
313  t = tdcTrans->getDriftTime(tdc, wireid, mop.getTime(), pocaOnWire.Z(), adc);
314  getObjectPtr<TTree>("tree")->Fill();
315  } //NDF
316  // }//end of if isU
317  }//end of for
318  }//end of for tp
319 }//end of func
320 
321 const CDCWire& CDCCalibrationCollectorModule::getIntersectingWire(const TVector3& xyz, const CDCWireLayer& layer,
322  const Helix& helixFit) const
323 {
324  Vector3D crosspoint;
325  if (layer.isAxial())
326  crosspoint = Vector3D(xyz);
327  else {
328  const CDCWire& oneWire = layer.getWire(1);
329  double newR = oneWire.getWirePos2DAtZ(xyz.z()).norm();
330  double arcLength = helixFit.getArcLength2DAtCylindricalR(newR);
331  TVector3 xyzOnWire = helixFit.getPositionAtArcLength2D(arcLength);
332  crosspoint = Vector3D(xyzOnWire);
333  }
334 
335  const CDCWire& wire = layer.getClosestWire(crosspoint);
336 
337  return wire;
338 }
339 
340 void CDCCalibrationCollectorModule::buildEfficiencies(std::vector<unsigned short> wireHits, const Helix helixFit)
341 {
342  static const TrackFindingCDC::CDCWireTopology& wireTopology = CDCWireTopology::getInstance();
343  for (const CDCWireLayer& wireLayer : wireTopology.getWireLayers()) {
344  const double radiusofLayer = wireLayer.getRefCylindricalR();
345  //simple extrapolation of fit
346  const double arcLength = helixFit.getArcLength2DAtCylindricalR(radiusofLayer);
347  const TVector3 xyz = helixFit.getPositionAtArcLength2D(arcLength);
348  if (!xyz.x()) continue;
349  const CDCWire& wireIntersected = getIntersectingWire(xyz, wireLayer, helixFit);
350  unsigned short crossedWire = wireIntersected.getEWire();
351  unsigned short crossedCWire = wireIntersected.getNeighborCW()->getEWire();
352  unsigned short crossedCCWire = wireIntersected.getNeighborCCW()->getEWire();
353 
354  if (find(wireHits.begin(), wireHits.end(), crossedWire) != wireHits.end()
355  || find(wireHits.begin(), wireHits.end(), crossedCWire) != wireHits.end()
356  || find(wireHits.begin(), wireHits.end(), crossedCCWire) != wireHits.end())
357  isFound = true;
358  else
359  isFound = false;
360 
361  wireID = wireIntersected.getIWire();
362  layerID = wireIntersected.getICLayer();
363  z = xyz.z();
364  getObjectPtr<TTree>("efftree")->Fill();
365  }
366 }
367 
368 
369 
Belle2::CalibrationCollectorModule
Calibration collector module base class.
Definition: CalibrationCollectorModule.h:44
genfit::TrackPoint
Object containing AbsMeasurement and AbsFitterInfo objects.
Definition: TrackPoint.h:46
Belle2::RelationArray
Low-level class to create/modify relations between StoreArrays.
Definition: RelationArray.h:72
Belle2::WireID
Class to identify a wire inside the CDC.
Definition: WireID.h:44
genfit::FitStatus::getChi2
double getChi2() const
Get chi^2 of the fit.
Definition: FitStatus.h:120
Belle2::Vector3D
HepGeom::Vector3D< double > Vector3D
3D Vector
Definition: Cell.h:35
Belle2::CDC::CDCGeometryPar::getOutgoingTheta
double getOutgoingTheta(const double alpha, const double theta) const
Converts incoming- to outgoing-theta.
Definition: CDCGeometryPar.cc:2700
Belle2::CDC::CDCCalibrationCollectorModule::harvest
void harvest(Belle2::RecoTrack *track)
collect hit information of fitted track.
Definition: CDCCalibrationCollector.cc:264
genfit::MeasuredStateOnPlane
#StateOnPlane with additional covariance matrix.
Definition: MeasuredStateOnPlane.h:39
Belle2::Track::getTrackFitResult
const TrackFitResult * getTrackFitResult(const Const::ChargedStable &chargedStable) const
Access to TrackFitResults.
Definition: Track.cc:19
Belle2::CDC::RealisticTDCCountTranslator
Translator mirroring the realistic Digitization.
Definition: RealisticTDCCountTranslator.h:36
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::CDC::CDCCalibrationCollectorModule::finish
void finish() override
Termination action.
Definition: CDCCalibrationCollector.cc:260
Belle2::CDCHit
Class containing the result of the unpacker in raw data and the result of the digitizer in simulation...
Definition: CDCHit.h:51
Belle2::TrackFitResult::getOmega
double getOmega() const
Getter for omega.
Definition: TrackFitResult.h:194
genfit
Defines for I/O streams used for error and debug printing.
Definition: AlignablePXDRecoHit.h:19
genfit::TrackPoint::getKalmanFitterInfo
KalmanFitterInfo * getKalmanFitterInfo(const AbsTrackRep *rep=nullptr) const
Helper to avoid casting.
Definition: TrackPoint.cc:180
Belle2::StoreAccessorBase::getName
const std::string & getName() const
Return name under which the object is saved in the DataStore.
Definition: StoreAccessorBase.h:130
Belle2::TrackFindingCDC::StoreWrappedObjPtr
This class is for convenience access and registration of objects, that are stored inside the StoreWra...
Definition: StoreWrappedObjPtr.h:40
Belle2::TrackFitResult::getTanLambda
double getTanLambda() const
Getter for tanLambda.
Definition: TrackFitResult.h:206
Belle2::TrackFindingCDC::CDCWireTopology::getWireLayers
const std::vector< Belle2::TrackFindingCDC::CDCWireLayer > & getWireLayers() const
Getter for the underlying storing layer vector.
Definition: CDCWireTopology.h:159
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::TrackFitResult::getZ0
double getZ0() const
Getter for z0.
Definition: TrackFitResult.h:200
Belle2::RecoTrack
This is the Reconstruction Event-Data Model Track.
Definition: RecoTrack.h:78
Belle2::CDC::CDCGeometryPar
The Class for CDC Geometry Parameters.
Definition: CDCGeometryPar.h:75
Belle2::CDC::CDCCalibrationCollectorModule::~CDCCalibrationCollectorModule
virtual ~CDCCalibrationCollectorModule()
Destructor.
Definition: CDCCalibrationCollector.cc:62
Belle2::CDC::CDCGeometryPar::getOutgoingLR
unsigned short getOutgoingLR(const unsigned short lr, const double alpha) const
Converts incoming-lr to outgoing-lr.
Definition: CDCGeometryPar.cc:2680
genfit::KalmanFitterInfo
Collects information needed and produced by a AbsKalmanFitter implementations and is specific to one ...
Definition: KalmanFitterInfo.h:44
Belle2::CDC::CDCGeometryPar::Instance
static CDCGeometryPar & Instance(const CDCGeometry *=nullptr)
Static method to get a reference to the CDCGeometryPar instance.
Definition: CDCGeometryPar.cc:41
Belle2::CDC::CDCGeometryPar::getTheta
double getTheta(const TVector3 &momentum) const
Returns track incident angle (theta in rad.).
Definition: CDCGeometryPar.cc:2674
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
genfit::FitStatus::isFitConverged
bool isFitConverged(bool inAllPoints=true) const
Did the fit converge (in all Points or only partially)?
Definition: FitStatus.h:105
Belle2::CDC::CDCCalibrationCollectorModule
Collect hit information for cdc calibration with CAF.
Definition: CDCCalibrationCollector.h:44
Belle2::CDC::CDCGeometryPar::getDriftTime
double getDriftTime(double dist, unsigned short layer, unsigned short lr, double alpha, double theta) const
Return the drift time to the sense wire.
Definition: CDCGeometryPar.cc:2457
Belle2::TrackFitResult::getTransverseMomentum
double getTransverseMomentum() const
Getter for the absolute value of the transverse momentum at the perigee.
Definition: TrackFitResult.h:140
Belle2::TrackFitResult::getPhi0
double getPhi0() const
Getter for phi0.
Definition: TrackFitResult.h:184
genfit::FitStatus::getNdf
double getNdf() const
Get the degrees of freedom of the fit.
Definition: FitStatus.h:122
Belle2::CDC::Helix
Helix parameter class.
Definition: Helix.h:51
Belle2::CDC::CDCCalibrationCollectorModule::prepare
void prepare() override
Initializes the Module.
Definition: CDCCalibrationCollector.cc:66
Belle2::CDC::CDCGeometryPar::getOutgoingAlpha
double getOutgoingAlpha(const double alpha) const
Converts incoming- to outgoing-alpha.
Definition: CDCGeometryPar.cc:2687
Belle2::CDC::RealisticTDCCountTranslator::getDriftTime
double getDriftTime(unsigned short tdcCount, const WireID &wireID, double timeOfFlightEstimator, double z, unsigned short adcCount) override
Get Drift time.
Definition: RealisticTDCCountTranslator.cc:53
genfit::KalmanFitterInfo::getFittedState
const MeasuredStateOnPlane & getFittedState(bool biased=true) const override
Get unbiased or biased (default) smoothed state.
Definition: KalmanFitterInfo.cc:179
Belle2::CDC::CDCCalibrationCollectorModule::getIntersectingWire
const TrackFindingCDC::CDCWire & getIntersectingWire(const TVector3 &xyz, const TrackFindingCDC::CDCWireLayer &layer, const Helix &helixFit) const
extrapolates the helix fit to a given layer and finds the wire which it would be hitting
Definition: CDCCalibrationCollector.cc:321
Belle2::Const::muon
static const ChargedStable muon
muon particle
Definition: Const.h:534
Belle2::Track
Class that bundles various TrackFitResults.
Definition: Track.h:35
genfit::KalmanFitterInfo::getWeights
std::vector< double > getWeights() const
Get weights of measurements.
Definition: KalmanFitterInfo.cc:168
Belle2::StoreArray< Belle2::Track >
Belle2::CDC::CDCCalibrationCollectorModule::collect
void collect() override
Event action, collect information for calibration.
Definition: CDCCalibrationCollector.cc:131
Belle2::TrackFitResult::getHelix
Helix getHelix() const
Conversion to framework Helix (without covariance).
Definition: TrackFitResult.h:230
genfit::FitStatus
Class where important numbers and properties of a fit can be stored.
Definition: FitStatus.h:80
genfit::FitStatus::isFitted
bool isFitted() const
Has the track been fitted?
Definition: FitStatus.h:94
Belle2::CDC::CDCCalibrationCollectorModule::buildEfficiencies
void buildEfficiencies(std::vector< unsigned short > wireHits, const Helix helixFit)
fills efficiency objects
Definition: CDCCalibrationCollector.cc:340
Belle2::CDC::CDCGeometryPar::getAlpha
double getAlpha(const TVector3 &posOnWire, const TVector3 &momentum) const
Returns track incident angle in rphi plane (alpha in rad.).
Definition: CDCGeometryPar.cc:2661
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::TrackFindingCDC::CDCWireTopology
Class representating the sense wire arrangement in the whole of the central drift chamber.
Definition: CDCWireTopology.h:54
Belle2::TrackFitResult::getChargeSign
short getChargeSign() const
Return track charge (1 or -1).
Definition: TrackFitResult.h:160
Belle2::TrackFitResult::getD0
double getD0() const
Getter for d0.
Definition: TrackFitResult.h:178