CbmRoot
CbmRecoQa.cxx
Go to the documentation of this file.
1 
6 #include <iostream>
7 #include <string>
8 #include <vector>
9 
10 #include "TClonesArray.h"
11 #include "TFile.h"
12 #include "TH1.h"
13 #include "TH1F.h"
14 #include "TList.h"
15 #include "TROOT.h"
16 #include "TString.h"
17 #include "TTask.h"
18 #include "TVector3.h"
19 
20 #include "CbmMuchPixelHit.h"
21 #include "CbmMuchPoint.h"
22 
23 #include "CbmTrdHit.h"
24 #include "CbmTrdPoint.h"
25 
26 #include "CbmTofHit.h"
27 #include "CbmTofPoint.h"
28 
29 #include "FairDetector.h"
30 #include "FairLogger.h"
31 #include "FairMCPoint.h"
32 #include "FairTask.h"
33 
34 #include "CbmStsDigi.h"
35 #include "CbmStsHit.h"
36 #include "CbmStsPoint.h"
37 
38 #include "CbmMvdHit.h"
39 #include "CbmMvdPoint.h"
40 
41 #include "CbmHit.h"
42 #include "CbmLink.h"
43 #include "CbmMCDataArray.h"
44 #include "CbmMatch.h"
45 #include "CbmRecoQa.h"
46 
48 
56  std::vector<std::pair<std::string, std::array<int, 4>>> decNames,
57  std::string outName,
58  int verbose_l)
59  : FairTask("CbmRecoQa")
60  , pullresfile(nullptr)
61  , verbosity(verbose_l)
62  , detectors(decNames)
63  , hists(
64  std::vector<std::vector<TH1F*>>(detectors.size(), std::vector<TH1F*>(6)))
65  , eventList(nullptr)
66  , fManager(nullptr)
67  , mcManager(nullptr)
68  , outname(outName) {
69  if (!instance) { instance = this; };
70 };
71 
72 // --- Destructor
74 
75 // --- Init and Re-Init
76 InitStatus CbmRecoQa::ReInit() { return kSUCCESS; };
77 
78 InitStatus CbmRecoQa::Init() {
79 
80  fManager = FairRootManager::Instance();
81  mcManager = (CbmMCDataManager*) fManager->GetObject("MCDataManager");
82  eventList = (CbmMCEventList*) fManager->GetObject("MCEventList.");
83 
84  if (!mcManager) {
85  LOG(error) << "Could not initialise MCDataManager" << std::endl;
86  return kERROR;
87  }
88  if (!fManager) {
89  LOG(warn) << "Could not initialise FairRootManager" << std::endl;
90  return kERROR;
91  }
92 
93  for (unsigned int i = 0; i < detectors.size(); i++) {
94 
95  // Histogram Setup
96  std::string decName = detectors[i].first;
97  std::string px = "Px_" + decName;
98  std::string py = "Py_" + decName;
99  std::string pt = "Pt_" + decName;
100  std::string rx = "x_" + decName;
101  std::string ry = "y_" + decName;
102  std::string rt = "t_" + decName;
103  std::string pullx = decName + " Pull x";
104  std::string pully = decName + " Pull y";
105  std::string pullt = decName + " Pull t";
106  std::string resx = decName + " Residual x";
107  std::string resy = decName + " Residual y";
108  std::string rest = decName + " Residual t";
109 
110  TH1F* pullxHist = new TH1F(px.c_str(),
111  pullx.c_str(),
112  100,
113  -1 * detectors[i].second[0],
114  detectors[i].second[0]);
115  TH1F* pullyHist = new TH1F(py.c_str(),
116  pully.c_str(),
117  100,
118  -1 * detectors[i].second[0],
119  detectors[i].second[0]);
120  TH1F* pulltHist = new TH1F(pt.c_str(),
121  pullt.c_str(),
122  100,
123  -1 * detectors[i].second[0],
124  detectors[i].second[0]);
125 
126  TH1F* residualxHist = new TH1F(rx.c_str(),
127  resx.c_str(),
128  100,
129  -1 * detectors[i].second[1],
130  detectors[i].second[1]);
131  TH1F* residualyHist = new TH1F(ry.c_str(),
132  resy.c_str(),
133  100,
134  -1 * detectors[i].second[2],
135  detectors[i].second[2]);
136  TH1F* residualtHist = new TH1F(rt.c_str(),
137  rest.c_str(),
138  100,
139  -1 * detectors[i].second[3],
140  detectors[i].second[3]);
141 
142  // pullxHist->SetCanExtend(TH1::kAllAxes);
143  // pullyHist->SetCanExtend(TH1::kAllAxes);
144  // pulltHist->SetCanExtend(TH1::kAllAxes);
145  // residualxHist->SetCanExtend(TH1::kAllAxes);
146  // residualyHist->SetCanExtend(TH1::kAllAxes);
147  // residualtHist->SetCanExtend(TH1::kAllAxes);
148 
149  pullxHist->GetXaxis()->SetTitle("Pull x");
150  pullyHist->GetXaxis()->SetTitle("Pull y");
151  pulltHist->GetXaxis()->SetTitle("Pull t");
152  residualxHist->GetXaxis()->SetTitle("Residual x");
153  residualyHist->GetXaxis()->SetTitle("Residual y");
154  residualtHist->GetXaxis()->SetTitle("Residual t");
155 
156  hists[i][0] = pullxHist;
157  hists[i][1] = pullyHist;
158  hists[i][2] = pulltHist;
159  hists[i][3] = residualxHist;
160  hists[i][4] = residualyHist;
161  hists[i][5] = residualtHist;
162 
163  if (verbosity > 1)
164  LOG(info) << "CbmRecoQa Success Initiliasing Histograms for: " << decName;
165  }
166 
167  return kSUCCESS;
168 }
169 
170 // --- Finish Event
171 // Record all Data for Detectors defined
173 
174  static int event = 0;
175  LOG(info) << "CbmRecoQa for Event " << event++;
176  for (unsigned int k = 0; k < detectors.size(); k++) {
177  instance->record(detectors[k].first, k);
178  }
179 }
180 
181 // --- Finish Task
182 // Save Data in File
184 
185  std::string filename = outname + ".qa.hists.root";
186  pullresfile = new TFile(filename.c_str(), "Recreate");
187 
188  for (unsigned int i = 0; i < hists.size(); i++) {
189  gDirectory->mkdir(detectors[i].first.c_str());
190  gDirectory->cd(detectors[i].first.c_str());
191  for (unsigned int k = 0; k < hists[i].size(); k++) {
192  if (verbosity > 2)
193  LOG(info) << "CbmRecoQa Histogramm " << hists[i][k]->GetName()
194  << " Entries " << hists[i][k]->GetEntries();
195  if (hists[i][k]->GetEntries() > 0) hists[i][k]->Write();
196  }
197  gDirectory->cd("..");
198  }
199 }
200 
201 // --- Write Data in Historgrams, depending on detectorw
202 void CbmRecoQa::record(std::string decName, int decNum) {
203 
204  if (verbosity > 1) LOG(info) << "CbmRecoQa Record called for: " << decName;
205 
206  // Data aquasition
207 
208  TClonesArray* listHits = nullptr;
209  TClonesArray* listHitMatches = nullptr;
210  // TClonesArray* listDigiMatches = nullptr;
211  CbmMCDataArray* listPoints = nullptr;
212  // TClonesArray* listClusters = nullptr; (FU) unused
213  TClonesArray* listClusterMatches = nullptr;
214 
215  if (decName == "sts") {
216  listHits = (TClonesArray*) (fManager->GetObject("StsHit"));
217  listHitMatches = (TClonesArray*) (fManager->GetObject("StsHitMatch"));
218  // listClusters= (TClonesArray*)(fManager->GetObject("StsCluster")); (FU) unused
219  listClusterMatches =
220  (TClonesArray*) (fManager->GetObject("StsClusterMatch"));
221  listPoints = mcManager->InitBranch("StsPoint");
222  if (listPoints == nullptr) { LOG(warn) << "No StsPoint data!"; }
223 
224 
225  if (listHits && listHitMatches) {
226 
227  TVector3 hitPos(.0, .0, .0);
228  TVector3 mcPos(.0, .0, .0);
229  TVector3 hitErr(.0, .0, .0);
230 
231  int nEnt = listHits->GetEntries();
232  if (verbosity > 2)
233  LOG(info) << "CbmRecoQa for " << decName << " found " << nEnt
234  << " Hit Entries";
235  for (int j = 0; j < nEnt; j++) {
236 
237  float bestWeight = 0.f;
238  float totalWeight = 0.f;
239  // int iMCPoint = -1; (FU) unused
240  CbmLink link;
241 
242  CbmStsHit* curr_hit = dynamic_cast<CbmStsHit*>(listHits->At(j));
243  CbmStsPoint* curr_point = 0;
244 
245  if (listClusterMatches) {
246 
247  const CbmMatch* front_match = static_cast<const CbmMatch*>(
248  listClusterMatches->At(curr_hit->GetFrontClusterId()));
249  const CbmMatch* back_match = static_cast<const CbmMatch*>(
250  listClusterMatches->At(curr_hit->GetBackClusterId()));
251  CbmMatch curr_match;
252  if (verbosity > 3)
253  LOG(info) << "Front Match: " << front_match->ToString()
254  << " Back Match: " << back_match->ToString();
255 
256  for (int frontlink_c = 0; frontlink_c < front_match->GetNofLinks();
257  frontlink_c++) {
258 
259  const CbmLink& frontLink = front_match->GetLink(frontlink_c);
260 
261  for (int backlink_c = 0; backlink_c < back_match->GetNofLinks();
262  backlink_c++) {
263 
264  const CbmLink& backLink = back_match->GetLink(backlink_c);
265  if (verbosity > 3)
266  LOG(info) << "FrontLink: " << frontLink.ToString()
267  << " BackLink: " << backLink.ToString();
268  if (frontLink == backLink) {
269  curr_match.AddLink(frontLink);
270  curr_match.AddLink(backLink);
271  }
272  }
273  }
274 
275  //LOG(info) << curr_match.ToString();
276 
277  for (int iLink = 0; iLink < curr_match.GetNofLinks(); iLink++) {
278  float tmpweight = curr_match.GetLink(iLink).GetWeight();
279  totalWeight += tmpweight;
280  //LOG(info) << " Match Link Nr.: " << iLink;
281  if (tmpweight > bestWeight) {
282  bestWeight = tmpweight;
283  // iMCPoint = curr_match.GetLink(iLink).GetIndex(); (FU) unused
284  link = curr_match.GetLink(iLink);
285  if (verbosity > 3) LOG(info) << "Found Link for current HIT";
286  }
287  }
288 
289  curr_point = (CbmStsPoint*) (listPoints->Get(
290  link.GetFile(), link.GetEntry(), link.GetIndex()));
291 
292  if (curr_point == 0) continue;
293  double mcTime =
294  curr_point->GetTime()
295  + eventList->GetEventTime(link.GetEntry(), link.GetFile());
296 
297  curr_hit->Position(hitPos);
298  curr_hit->PositionError(hitErr);
299  if (verbosity > 3) LOG(info) << "Calculated Position Error";
300 
301  mcPos.SetX(curr_point->GetX(hitPos.Z()));
302  mcPos.SetY(curr_point->GetY(hitPos.Z()));
303  mcPos.SetZ(hitPos.Z());
304  if (verbosity > 3) LOG(info) << "Calculated MCPos";
305 
306 
307  if (hitErr.X() != 0)
308  hists[decNum][0]->Fill((hitPos.X() - mcPos.X())
309  / curr_hit->GetDx());
310  if (hitErr.Y() != 0)
311  hists[decNum][1]->Fill((hitPos.Y() - mcPos.Y())
312  / curr_hit->GetDy());
313  hists[decNum][2]->Fill((curr_hit->GetTime() - mcTime)
314  / curr_hit->GetTimeError());
315 
316  hists[decNum][3]->Fill((hitPos.X() - mcPos.X()) * 10 * 1000);
317  hists[decNum][4]->Fill((hitPos.Y() - mcPos.Y()) * 10 * 1000);
318  hists[decNum][5]->Fill(curr_hit->GetTime() - mcTime);
319  } else {
320  LOG(warn) << "CBMRECOQA WARNING :-- No Sts Cluster Matches found!"
321  << std::endl;
322  }
323  }
324  } else {
325  LOG(warn) << "CBMRECOQA WARNING :-- NO Data for Reco QA found! "
326  << "Detector: " << decName << std::endl;
327  }
328  } else if (decName == "mvd") {
329  listHits = (TClonesArray*) (fManager->GetObject("MvdHit"));
330  listHitMatches = (TClonesArray*) (fManager->GetObject("MvdHitMatch"));
331  //listClusters= (TClonesArray*)(fManager->GetObject("MvdCluster"));
332  //listClusterMatches = (TClonesArray*)(fManager->GetObject("MvdClusterMatch"));
333  TClonesArray* listMvdPoints =
334  (TClonesArray*) (fManager->GetObject("MvdPoint"));
335 
336  if (listHits && listHitMatches) {
337 
338  TVector3 hitPos(.0, .0, .0);
339  TVector3 mcPos(.0, .0, .0);
340  TVector3 hitErr(.0, .0, .0);
341 
342  int nEnt = listHits->GetEntries();
343  if (verbosity > 2)
344  LOG(info) << "CbmRecoQa for " << decName << " found " << nEnt
345  << " Hit Entries";
346  for (int j = 0; j < nEnt; j++) {
347 
348  // float bestWeight = 0.f;
349  // float totalWeight = 0.f;
350  int iMC = -1;
351 
352  CbmMvdHit* curr_hit = dynamic_cast<CbmMvdHit*>(listHits->At(j));
353  CbmMatch* curr_match = dynamic_cast<CbmMatch*>(listHitMatches->At(j));
354  CbmMvdPoint* curr_point = 0;
355 
356  if (curr_match->GetNofLinks() > 0)
357  iMC = curr_match->GetLink(0).GetIndex();
358  if (iMC < 0) continue;
359 
360  curr_point = (CbmMvdPoint*) (listMvdPoints->At(iMC));
361 
362  if (curr_point == 0) continue;
363 
364  curr_hit->Position(hitPos);
365  curr_hit->PositionError(hitErr);
366  if (verbosity > 3) LOG(info) << "Calculated Position Error";
367 
368  mcPos.SetX((curr_point->GetX() + curr_point->GetXOut()) / 2.);
369  mcPos.SetY((curr_point->GetY() + curr_point->GetYOut()) / 2.);
370  mcPos.SetZ(hitPos.Z());
371  if (verbosity > 3) LOG(info) << "Calculated MCPos";
372 
373 
374  if (hitErr.X() != 0)
375  hists[decNum][0]->Fill((hitPos.X() - mcPos.X()) / curr_hit->GetDx());
376  if (hitErr.Y() != 0)
377  hists[decNum][1]->Fill((hitPos.Y() - mcPos.Y()) / curr_hit->GetDy());
378 
379  hists[decNum][3]->Fill((hitPos.X() - mcPos.X()) * 10 * 1000);
380  hists[decNum][4]->Fill((hitPos.Y() - mcPos.Y()) * 10 * 1000);
381  }
382  } else {
383  LOG(warn) << "CBMRECOQA WARNING :-- NO Data for Reco QA found! "
384  << "Detector: " << decName << std::endl;
385  }
386  } else if (decName == "much") {
387 
388  listHits = (TClonesArray*) (fManager->GetObject("MuchPixelHit"));
389  listHitMatches = (TClonesArray*) (fManager->GetObject("MuchPixelHitMatch"));
390  listPoints = mcManager->InitBranch("MuchPoint");
391  if (listPoints == nullptr) { LOG(warn) << "No MuchPoint data!"; }
392 
393  if (listHits && listHitMatches) {
394 
395  TVector3 hitPos(.0, .0, .0);
396  TVector3 mcPos(.0, .0, .0);
397  TVector3 hitErr(.0, .0, .0);
398 
399  int nEnt = listHits->GetEntries();
400  if (verbosity > 2)
401  LOG(info) << "CbmRecoQa for " << decName << " found " << nEnt
402  << " Hit Entries";
403  for (int j = 0; j < nEnt; j++) {
404 
405  float bestWeight = 0.f;
406  float totalWeight = 0.f;
407  int iMC = -1;
408  CbmLink link;
409 
410  CbmMuchPixelHit* curr_hit =
411  dynamic_cast<CbmMuchPixelHit*>(listHits->At(j));
412  CbmMatch* curr_match = dynamic_cast<CbmMatch*>(listHitMatches->At(j));
413 
414  for (int iLink = 0; iLink < curr_match->GetNofLinks(); iLink++) {
415  float tmpweight = curr_match->GetLink(iLink).GetWeight();
416  totalWeight += tmpweight;
417  //LOG(info) << " Match Link Nr.: " << iLink;
418  if (tmpweight > bestWeight) {
419  bestWeight = tmpweight;
420  iMC = curr_match->GetLink(iLink).GetIndex();
421  link = curr_match->GetLink(iLink);
422  if (verbosity > 3) LOG(info) << "Found Link for current HIT";
423  }
424  }
425 
426  if (iMC < 0) continue;
427 
428  CbmMuchPoint* curr_point = (CbmMuchPoint*) (listPoints->Get(
429  link.GetFile(), link.GetEntry(), link.GetEntry()));
430  double mcTime =
431  curr_point->GetTime()
432  + eventList->GetEventTime(link.GetEntry(), link.GetFile());
433  //LOG(info) << "Much Hit " << j << " Time: " << mcTime << " " << curr_hit->GetTime() << " " << curr_hit->GetTimeError();
434 
435  if (curr_point == 0) continue;
436 
437  curr_hit->Position(hitPos);
438  curr_hit->PositionError(hitErr);
439  if (verbosity > 3) LOG(info) << "Calculated Position Error";
440 
441  mcPos.SetX((curr_point->GetXIn() + curr_point->GetXOut()) / 2.);
442  mcPos.SetY((curr_point->GetYIn() + curr_point->GetYOut()) / 2.);
443  mcPos.SetZ(hitPos.Z());
444  if (verbosity > 3) LOG(info) << "Calculated MCPos";
445 
446  curr_hit->Position(hitPos);
447 
448  if (hitErr.X() != 0)
449  hists[decNum][0]->Fill((hitPos.X() - mcPos.X()) / curr_hit->GetDx());
450  if (hitErr.Y() != 0)
451  hists[decNum][1]->Fill((hitPos.Y() - mcPos.Y()) / curr_hit->GetDy());
452  if (hitErr.Y() != 0)
453  hists[decNum][2]->Fill(
454  ((curr_hit->GetTime() - mcTime) / curr_hit->GetTimeError()));
455 
456  hists[decNum][3]->Fill((hitPos.X() - mcPos.X()));
457  hists[decNum][4]->Fill((hitPos.Y() - mcPos.Y()));
458  hists[decNum][5]->Fill(curr_hit->GetTime() - mcTime);
459  }
460  } else {
461  LOG(warn) << "CBMRECOQA WARNING :-- NO Data for Reco QA found! "
462  << "Detector: " << decName << std::endl;
463  }
464 
465  } else {
466  LOG(warn) << "CBMRECOQA WARNING :-- NO matching Detector found ! "
467  << std::endl;
468  }
469 }
CbmMuchPoint::GetXOut
Double_t GetXOut() const
Definition: CbmMuchPoint.h:73
CbmMatch
Definition: CbmMatch.h:22
CbmRecoQa::fManager
FairRootManager * fManager
Definition: CbmRecoQa.h:30
CbmRecoQa::pullresfile
TFile * pullresfile
Definition: CbmRecoQa.h:25
CbmMuchPoint
Definition: CbmMuchPoint.h:21
CbmPixelHit::Position
void Position(TVector3 &pos) const
Copies hit position to pos.
Definition: CbmPixelHit.cxx:64
CbmMatch::GetLink
const CbmLink & GetLink(Int_t i) const
Definition: CbmMatch.h:35
CbmMatch::GetNofLinks
Int_t GetNofLinks() const
Definition: CbmMatch.h:38
CbmRecoQa::verbosity
int verbosity
Definition: CbmRecoQa.h:26
CbmPixelHit::GetDx
Double_t GetDx() const
Definition: CbmPixelHit.h:85
CbmRecoQa::Init
virtual InitStatus Init()
Definition: CbmRecoQa.cxx:78
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmRecoQa::ReInit
virtual InitStatus ReInit()
Definition: CbmRecoQa.cxx:76
CbmRecoQa::mcManager
CbmMCDataManager * mcManager
Definition: CbmRecoQa.h:31
CbmMCDataManager::InitBranch
CbmMCDataArray * InitBranch(const char *name)
Definition: CbmMCDataManager.cxx:106
CbmRecoQa::FinishEvent
virtual void FinishEvent()
Definition: CbmRecoQa.cxx:172
CbmMCDataArray.h
CbmMCDataArray
Access to a MC data branch for time-based analysis.
Definition: CbmMCDataArray.h:35
CbmStsHit::GetFrontClusterId
Int_t GetFrontClusterId() const
Definition: CbmStsHit.h:93
CbmRecoQa::eventList
CbmMCEventList * eventList
Definition: CbmRecoQa.h:29
CbmRecoQa
Definition: CbmRecoQa.h:22
CbmStsPoint
Definition: CbmStsPoint.h:27
CbmMatch.h
CbmPixelHit::GetDy
Double_t GetDy() const
Definition: CbmPixelHit.h:86
CbmHit::GetTimeError
Double_t GetTimeError() const
Definition: CbmHit.h:76
CbmMvdHit
Definition: CbmMvdHit.h:29
CbmMuchPoint::GetYOut
Double_t GetYOut() const
Definition: CbmMuchPoint.h:74
CbmStsHit
data class for a reconstructed 3-d hit in the STS
Definition: CbmStsHit.h:31
CbmMvdPoint.h
CbmStsDigi.h
CbmHit.h
CbmMvdPoint
Definition: CbmMvdPoint.h:28
CbmMuchPoint.h
CbmMCDataArray::Get
TObject * Get(const CbmLink *lnk)
Definition: CbmMCDataArray.h:47
CbmHit::GetTime
Double_t GetTime() const
Definition: CbmHit.h:75
CbmMCEventList::GetEventTime
Double_t GetEventTime(UInt_t event, UInt_t file)
Event start time.
Definition: CbmMCEventList.cxx:86
CbmMatch::AddLink
void AddLink(const CbmLink &newLink)
Definition: CbmMatch.cxx:42
CbmRecoQa::hists
std::vector< std::vector< TH1F * > > hists
Definition: CbmRecoQa.h:28
CbmTrdHit.h
Class for hits in TRD detector.
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmMCEventList
Container class for MC events with number, file and start time.
Definition: CbmMCEventList.h:38
first
bool first
Definition: LKFMinuit.cxx:143
CbmRecoQa.h
CbmMCDataManager
Task class creating and managing CbmMCDataArray objects.
Definition: CbmMCDataManager.h:27
CbmRecoQa::detectors
std::vector< std::pair< std::string, std::array< int, 4 > > > detectors
Definition: CbmRecoQa.h:27
CbmRecoQa::instance
static CbmRecoQa * instance
Definition: CbmRecoQa.h:39
CbmStsPoint.h
CbmStsHit::GetBackClusterId
Int_t GetBackClusterId() const
Definition: CbmStsHit.h:69
CbmMuchPoint::GetYIn
Double_t GetYIn() const
Definition: CbmMuchPoint.h:71
CbmMvdPoint::GetXOut
Double_t GetXOut() const
Definition: CbmMvdPoint.h:70
CbmMvdHit.h
CbmStsPoint::GetX
Double_t GetX(Double_t z) const
Definition: CbmStsPoint.cxx:96
CbmTofPoint.h
CbmMuchPixelHit.h
Class for pixel hits in MUCH detector.
CbmTrdPoint.h
CbmRecoQa::record
void record(std::string decName, int i)
Definition: CbmRecoQa.cxx:202
CbmMvdPoint::GetYOut
Double_t GetYOut() const
Definition: CbmMvdPoint.h:71
CbmMuchPoint::GetXIn
Double_t GetXIn() const
Definition: CbmMuchPoint.h:70
CbmRecoQa::outname
std::string outname
Definition: CbmRecoQa.h:32
CbmRecoQa::FinishTask
virtual void FinishTask()
Definition: CbmRecoQa.cxx:183
CbmMuchPixelHit
Definition: CbmMuchPixelHit.h:17
CbmRecoQa::CbmRecoQa
CbmRecoQa(std::vector< std::pair< std::string, std::array< int, 4 >>> decNames, std::string outName="test", int verbose_l=0)
Definition: CbmRecoQa.cxx:55
CbmPixelHit::PositionError
void PositionError(TVector3 &dpos) const
Copies hit position error to pos.
Definition: CbmPixelHit.cxx:68
CbmRecoQa::~CbmRecoQa
~CbmRecoQa()
Definition: CbmRecoQa.cxx:73
CbmStsPoint::GetY
Double_t GetY(Double_t z) const
Definition: CbmStsPoint.cxx:106
CbmMatch::ToString
virtual std::string ToString() const
Return string representation of the object.
Definition: CbmMatch.cxx:21
CbmStsHit.h
Data class for a reconstructed hit in the STS.