FairRoot/PandaRoot
runOnlineDisplayMCCheck.C
Go to the documentation of this file.
1 #include "drawhelpers.hxx"
2 
3 void SaveAndUpdateHisto(TH1* currenthisto, TFile& storagefile)
4 {
5  if (storagefile.Get(currenthisto->GetName()) != 0) {
6  currenthisto->Add((TH1*)storagefile.Get(currenthisto->GetName()));
7  }
8  cout << currenthisto->GetName() << ": " << currenthisto->GetEntries() << endl;
9  currenthisto->Write(currenthisto->GetName(), TObject::kWriteDelete);
10 }
11 
12 int runOnlineDisplayMCCheck(Int_t maximumTime = 2050) {
13  gROOT->Reset();
14  // Verbosity level (0=quiet, 1=event level, 2=track level, 3=debug)
15  Int_t iVerbose = 0;
16 
17  enum DetectorId {
19 
20  TStopwatch timer;
21  timer.Start();
22  gROOT->LoadMacro("$VMCWORKDIR/gconfig/rootlogon.C");
23  rootlogon();
24 
25  TDatabasePDG *pdg = new TDatabasePDG();
26 
27  // Reco Parameters
28  Double_t stthitlifetime = 245;
29 
30  TString simFileName = "Sim_Dpm_500.root";
31  TString parFileName = "Sim_Dpm_500_params.root";
32  TString digiFileName = "Sim_Dpm_500_digi.root";
33  TString recoFileName = "Sim_Dpm_500_reco.root";
34  TString outFileName = "Sim_Dpm_500_streamdisplay.root";
35  TString histofilename = "onlinetrkmccheck.root";
36 
37  TFile filedigi(digiFileName.Data());
38  TFile filereco(recoFileName.Data());
39  TFile filerecopixel(recoFileName.Data());
40 
41  TTree *treedigi = (TTree*) filedigi.Get("pndsim");
42  TTree *recotree = (TTree*) filereco.Get("pndsim");
43 
44  TFile simfile(simFileName.Data());
45  TTree* mctree = (TTree*) simfile.Get("pndsim");
46  TClonesArray* mcarray = new TClonesArray("PndMCTrack");
47  mctree->SetBranchAddress("MCTrack", &mcarray);
48 
49  FairRunAna *fRun = new FairRunAna();
50  fRun->SetInputFile(simFileName.Data());
51  fRun->AddFriend(recoFileName.Data());
52  fRun->AddFriend(digiFileName.Data());
53  fRun->SetOutputFile(outFileName.Data());
54  fRun->Init();
55  FairRuntimeDb* rtdb = fRun->GetRuntimeDb();
56 
57  PndOnlineGeometryManager *online_geometry = new PndOnlineGeometryManager(rtdb, parFileName.Data());
58 
59  PndOnlineManager *online = new PndOnlineManager();
60  online->AddGeometryManager(online_geometry);
61 
62  online->LoadStream(PndOnline::kSTT,"STTSortedHits",stthitlifetime);
63  online->LoadStream(PndOnline::kMVDPixel,"MVDHitsPixel",stthitlifetime);
64  online->LoadStream(PndOnline::kMVDStrip,"MVDHitsStrip",stthitlifetime);
65  //online->SetHESRRevolutionDuration(3000);
66 
67  TClonesArray* tubearray = online_geometry->GetDetectorGeometry(PndOnline::kSTT);
68 
69  // load task
70  PndOnlineSttTripletFinder *triplet_finder = new PndOnlineSttTripletFinder(online, tubearray);
71  online->AddTask((FairTask*)triplet_finder);
72 
73  // initialize display
74  TCanvas* c1 = new TCanvas("c1");
75  c1->Range(-42,-42,42,42);
76  c1->SetCanvasSize(1200, 1200);
77  TText* mytext = new TText();
78 
79  // event loop!
80  int current_time = 0;
81  int drawtime = 0;
82  int delta_t = 2050;
83  int drawdeltat = 5;
84  int final_time = maximumTime;
85 
86  online->SetHESRRevolutionDuration(delta_t);
87 
88  online->Init();
89 
90  TObjArray* online_fairhits_stt = 0;
91  TObjArray* online_tracks = 0;
92  TObjArray* online_tracks2 = 0;
93  TObjArray* online_fairhits_mvd = 0;
94 
95  std::vector<double> mctracktime;
96  std::vector<double> mctrackpt;
97  std::vector<double> mctrackphi;
98  std::vector<int> mctrackstatus;
99  std::vector<double> recotracktime;
100  std::vector<double> recotrackpt;
101  std::vector<double> recotrackphi;
102  std::vector<int> recotrackstatus;
103  std::vector<int> recotrackmcid;
104  int totalmctracks = 0;
105  int goodmctracks = 0;
106  int totalrecotracks = 0;
107  int goodrecotracks = 0;
108  int norecotracks = 0;
109  int deltaphifailed = 0;
110  int deltaptfailed = 0;
111  int deltabothfailed = 0;
112  int secondaries = 0;
113  int globalmaxmcid = 0;
114 
115  //std::map<int, std::map<int, int> > foundmap90;
116  //std::map<int, std::map<int, int> > foundmap80;
117 
118  int mcidgood90 = 0;
119  int mcidgood80 = 0;
120  int mcidbad = 0;
121 
122  //output histos
123  TH2D* hmcall = new TH2D("hmcall", "hmcall", 1500, -15, 15, 75, 0, 1.5);
124  TH2D* hmcprime = new TH2D("hmcprime", "hmcprime", 1500, -15, 15, 75, 0, 1.5);
125  TH2D* hmcrecoable = new TH2D("hmcrecoable", "hmcrecoable", 1500, -15, 15, 75, 0, 1.5);
126  TH2D* hdafuq = new TH2D("hdafuq", "hdafuq", 1500, -15, 15, 75, 0, 1.5);
127  TH2D* hmcrecoed = new TH2D("hmcrecoed", "hmcrecoed", 1500, -15, 15, 75, 0, 1.5);
128  TH2D* hmcnotrecoed = new TH2D("hmcnotrecoed", "hmcnotrecoed", 1500, -15, 15, 75, 0, 1.5);
129 
130  TH1D* hphires = new TH1D("hphires", "hphires", 3000, -15, 15);
131  TH2D* hptres = new TH2D("hptres", "hptres", 300, -1500, 1500, 200, -1, 1);
132  TH1D* hphiresmc = new TH1D("hphiresmc", "hphiresmc", 3000, -15, 15);
133  TH2D* hptresmc = new TH2D("hptresmc", "hptresmc", 300, -1500, 1500, 200, -1, 1);
134  TH1D* hphiresmc90 = new TH1D("hphiresmc90", "hphiresmc90", 3000, -15, 15);
135  TH2D* hptresmc90 = new TH2D("hptresmc90", "hptresmc90", 300, -1500, 1500, 200, -1, 1);
136  TH1D* hphiresmc80 = new TH1D("hphiresmc80", "hphiresmc80", 3000, -15, 15);
137  TH2D* hptresmc80 = new TH2D("hptresmc80", "hptresmc80", 300, -1500, 1500, 200, -1, 1);
138  TH1D* hphiresmcbad = new TH1D("hphiresmcbad", "hphiresmcbad", 3000, -15, 15);
139  TH2D* hptresmcbad = new TH2D("hptresmcbad", "hptresmcbad", 300, -1500, 1500, 200, -1, 1);
140 
141 
142  for (current_time = 0; current_time < final_time; current_time += delta_t) {
143  online->Clear();
144  //online->LoadHits( PndOnlineManager::kHESRRevolution );
145  online->ClearTracks();
146  online->LoadHits( delta_t );
147  //online->LoadHits( current_time );
148  online->Process();
149  //online->Clear();
150 
151  online_tracks = online->GetTrackObjectList();
152  int maxmcid = 0;
153  int minmcid = 0;
154  int firsti = 0;
155  if (online_tracks->GetEntriesFast() > 0) {
156  PndOnlineTrack* firstTrack = (PndOnlineTrack*) online_tracks->At(firsti);
157  while (firstTrack->DrawMode() != PndOnlineTrack::kVertexMomentumAsCircleHack) {
158  firsti++;
159  if (firsti > online_tracks->GetEntriesFast()) {
160  cout << "No drawable online Track in Array. We have a severe problem." << endl;
161  break;
162  }
163  firstTrack = (PndOnlineTrack*) online_tracks->At(firsti);
164  }
165  maxmcid = firstTrack->MCEntryNumber();
166  minmcid = firstTrack->MCEntryNumber();
167  }
168  for (int i = 0; i < online_tracks->GetEntriesFast(); i++) {
169  PndOnlineTrack* myTrack = (PndOnlineTrack*) online_tracks->At(i);
170  if (myTrack->DrawMode() == PndOnlineTrack::kVertexMomentumAsCircleHack) {
171  TVector2 innerhit;
172  innerhit.Set(myTrack->Vertex().x(), myTrack->Vertex().y());
173  TVector2 outerhit;
174  outerhit.Set(myTrack->Momentum().x(), myTrack->Momentum().y());
175  TVector2 circleorigin;
176  Double_t circleradius;
177  Int_t clockwise;
178  CalculateCircle(&innerhit, &outerhit, &circleorigin, &circleradius, &clockwise);
179  Double_t trackpt = circleradius * 2.0 * 3.0;
180  Double_t orthophi = circleorigin.Phi();
181  Double_t trackphi = orthophi + 0.5*TMath::Pi()*clockwise;
182  if (trackphi < 0) {
183  trackphi += 2*TMath::Pi();
184  }
185  if (trackphi > 2*TMath::Pi()) {
186  trackphi -= 2*TMath::Pi();
187  }
188 
189  cout << "Reco Track MCEntryNumber: " << myTrack->MCEntryNumber() << endl;
190  cout << "Reco Track MCEntryMajority: " << myTrack->MCEntryMajority() << endl;
191  cout << "Reco Track MCTrackNumber: " << myTrack->MCTrackNumber() << endl;
192  cout << "Reco Track MCTrackMajority: " << myTrack->MCTrackMajority() << endl;
193  cout << "Reco Track HitCount: " << myTrack->HitCount() << endl;
194 
195  mctree->GetEntry(myTrack->MCEntryNumber());
196  cout << "MC Entry count: " << mcarray->GetEntriesFast() << endl;
197  PndMCTrack* mctrack = (PndMCTrack*)mcarray->At(myTrack->MCTrackNumber());
198  if (mctrack != 0) {
199  if (mctrack->GetMotherID() == -1) {
200  Double_t mcphi = mctrack->GetMomentum().Phi();
201  Double_t mcpt = mctrack->GetMomentum().Pt()*1000;
202  if (mcphi < 0) {
203  mcphi += 2*TMath::Pi();
204  }
205  Double_t calcphimc = mcphi;
206  if (calcphimc < 0) calcphimc += 2*TMath::Pi();
207  Double_t calcphireco = trackphi;
208  if (calcphireco < 0) calcphireco += 2*TMath::Pi();
209  Double_t deltaphi = calcphireco - calcphimc;
210  if (deltaphi < 0) {
211  deltaphi *= -1;
212  }
213  if (deltaphi > TMath::Pi()) {
214  deltaphi = 2*TMath::Pi() - deltaphi;
215  }
216  Double_t deltapt = trackpt - mcpt;
217  Int_t chargesign = 1;
218  if (pdg->GetParticle(mctrack->GetPdgCode())->Charge() < 0) {
219  chargesign = -1;
220  }
221  hphiresmc->Fill(deltaphi*chargesign);
222  hptresmc->Fill(mcpt, deltapt/mcpt);
223  if ( (myTrack->MCEntryMajority() > 0.9) && (myTrack->MCTrackMajority() > 0.9) ) {
224  hphiresmc90->Fill(deltaphi*chargesign);
225  hptresmc90->Fill(mcpt, deltapt/mcpt);
226  }
227  if ( (myTrack->MCEntryMajority() > 0.8) && (myTrack->MCTrackMajority() > 0.8) ) {
228  hphiresmc80->Fill(deltaphi*chargesign);
229  hptresmc80->Fill(mcpt, deltapt/mcpt);
230  } else {
231  hphiresmcbad->Fill(deltaphi*chargesign);
232  hptresmcbad->Fill(mcpt, deltapt/mcpt);
233  }
234  }
235  } else {
236  cout << "MC Track has NULL pointer. Oooops." << endl;
237  continue;
238  }
239 
240  //bool isrecoable = true;
241  if (mctrack->GetMomentum().Pt() < 0.052) isrecoable = false;
242  if ( (mctrack->GetMomentum().Pt() < 0.052) && (mctrack->GetNPoints(kSTT) > 0) ) {
243  cout << "Low pt hit" << endl;
244  }
245  if ( (myTrack->MCEntryMajority() > 0.9) && (myTrack->MCTrackMajority() > 0.9) ) {
246  //(foundmap90[myTrack->MCEntryNumber()][myTrack->MCTrackNumber()])++
247  mcidgood90++;
248  }
249  if ( (myTrack->MCEntryMajority() > 0.8) && (myTrack->MCTrackMajority() > 0.8) ) {
250  //(foundmap80[myTrack->MCEntryNumber()][myTrack->MCTrackNumber()])++
251  mcidgood80++;
252  } else {
253  mcidbad++;
254  }
255 
256  int mcentrynumber = myTrack->MCEntryNumber();
257  if (mcentrynumber > maxmcid) {
258  maxmcid = mcentrynumber;
259  }
260  if (mcentrynumber < minmcid) {
261  minmcid = mcentrynumber;
262  }
263 
264  //std::cout << "RecoTrack: " << totalrecotracks << ", MCEntryNumber: " << mcentrynumber << " DigiTime: " << myTrack->T0Min() << ", T0time: " << myTrack->T0() << ", Pt " << trackpt << ", Phi: " << trackphi << " TriPhi: " << innerhit.Phi() << ", Radius: " << circleradius << ", DeltaT0: " << myTrack->T0Max()-myTrack->T0Min() << endl;
265 
266  //std::cout << "RecoTrack: " << totalrecotracks << ", MCEntryNumber: " << mcentrynumber << ", Pt " << trackpt << ", Phi: " << trackphi << " TriPhi: " << innerhit.Phi() << ", Radius: " << circleradius << endl;
267 
268  //recotracktime.push_back(myTrack->T0Min());
269  recotrackpt.push_back(trackpt);
270  recotrackphi.push_back(trackphi);
271  recotrackstatus.push_back(0);
272  recotrackmcid.push_back(mcentrynumber);
273  totalrecotracks++;
274  }
275  }
276  if (maxmcid > globalmaxmcid) {
277  globalmaxmcid = maxmcid;
278  }
279  //cout << "MC ID Range: " << minmcid << ", " << globalmaxmcid << endl;
280  for (int j = minmcid; j < mctree->GetEntriesFast(); j++) {
281  cout << "MC outer loop: " << j << endl;
282  if (j > globalmaxmcid) {
283  break;
284  }
285  mctree->GetEntry(j);
286  cout << "MC Entry count: " << mcarray->GetEntriesFast() << endl;
287  for (int k = 0; k < mcarray->GetEntriesFast(); k++) {
288  PndMCTrack* mctrack = (PndMCTrack*)mcarray->At(k);
289  if (mctrack == 0) {
290  cout << "MC Track has NULL pointer. Oooops." << endl;
291  continue;
292  }
293  hmcall->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
294  bool isrecoable = false;
295  if (mctrack->GetMotherID() == -1) {
296  hmcprime->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
297  isrecoable = true;
298  if (mctrack->GetMomentum().Pz() < -0.4) {
299  cout << "Critical Track at Event: " << j << " Track " << k << " Particle: " << mctrack->GetPdgCode() << endl;
300  //cout << "Particle Mass: " << pdg->GetParticle(mctrack->GetPdgCode())->Mass() << " Charge: " << pdg->GetParticle(mctrack->GetPdgCode())->Charge() << endl;
301  //cout << "MCTrack STT Hits: " << mctrack->GetNPoints(kSTT) << ", Pz: " << mctrack->GetMomentum().Pz() << ", Pt: " << mctrack->GetMomentum().Pt() << ", Phi: " << mctrack->GetMomentum().Phi() << ", Radius: " << mctrack->GetMomentum().Pt() / (3.0*2.0) << endl;
302  cout << "MCVertex: " << mctrack->GetStartVertex().Mag() << " XYZ: " << mctrack->GetStartVertex().x() << ", " << mctrack->GetStartVertex().y() << ", " << mctrack->GetStartVertex().z() << endl;
303  }
304  }
305  //if (mctrack->GetMotherID() == -1 && mctrack->GetNPoints(kSTT) > 0) {
306  //cout << TRho::Instance()->GetPDG()->GetParticle(mctrack->GetPdgCode())->Mass() << endl;
307  //}
308 
309  //bool isrecoable = true;
310  if (mctrack->GetMomentum().Pt() < 0.052) isrecoable = false;
311  if ( (mctrack->GetMomentum().Pt() < 0.052) && (mctrack->GetNPoints(kSTT) > 0) ) {
312  cout << "Low pt hit" << endl;
313  }
314  if ( (mctrack->GetMomentum().Pt()/mctrack->GetMomentum().Pz() < 28.0/110.0) && (mctrack->GetMomentum().Pz() > 0) ) isrecoable = false;
315  if ( (-(mctrack->GetMomentum().Pt())/mctrack->GetMomentum().Pz() < 28.0/55.0) && (mctrack->GetMomentum().Pz() < 0) ) isrecoable = false;
316  if ( (mctrack->GetMotherID() == -1) && isrecoable && (mctrack->GetNPoints(kSTT) < 3) ) {
317  if (TMath::Abs(pdg->GetParticle(mctrack->GetPdgCode())->Charge()) < 1) isrecoable = false;
318  if (pdg->GetParticle(mctrack->GetPdgCode())->Mass() > 1) isrecoable = false;
319  }
320  if ( (mctrack->GetMotherID() == -1) && isrecoable && (mctrack->GetMomentum().Mag()/pdg->GetParticle(mctrack->GetPdgCode())->Mass() < 0.4) ) {
321  cout << "Low betagamma particle" << endl;
322  isrecoable = false;
323  }
324  if ( (mctrack->GetMotherID() == -1) && isrecoable && (mctrack->GetNPoints(kSTT) < 3) ) {
325  cout << "Meh, reco possibility mismatch." << endl;
326  hdafuq->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
327  cout << "Dafuq at Event: " << j << " Track " << k << " Particle: " << mctrack->GetPdgCode() << " Charge: " << pdg->GetParticle(mctrack->GetPdgCode())->Charge() << endl;
328  mctrack->GetMomentum().Print();
329  isrecoable = false;
330  }
331 
332  //if (mctrack->GetMotherID() == -1 && mctrack->GetNPoints(kSTT) > 0) {
333  if (mctrack->GetMotherID() == -1 && isrecoable) {
334  hmcrecoable->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
335  Double_t mcphi = mctrack->GetMomentum().Phi();
336  Double_t mcpt = mctrack->GetMomentum().Pt()*1000;
337  if (mcphi < 0) {
338  mcphi += 2*TMath::Pi();
339  }
340  //cout << "MCTrack: " << j << " " << k << ", Pt: " << mcpt << ", Phi: " << mcphi << ", Radius: " << mcpt / (3.0*2.0) << endl;
341  //cout << "MCVertex: " << mctrack->GetStartVertex().Mag() << " XYZ: " << mctrack->GetStartVertex().x() << ", " << mctrack->GetStartVertex().y() << ", " << mctrack->GetStartVertex().z() << endl;
342  int trackstatus = 0;
343  if (mctrack->GetStartVertex().Mag() < 2) {
344  for (int r = 0; r < recotrackmcid.size(); r++) {
345  if (recotrackmcid.at(r) == j) {
346  //Double_t deltaphi = recotrackphi.at(r) - mcphi;
347  Double_t calcphimc = mcphi;
348  if (calcphimc < 0) calcphimc += 2*TMath::Pi();
349  Double_t calcphireco = recotrackphi.at(r);
350  if (calcphireco < 0) calcphireco += 2*TMath::Pi();
351  Double_t deltaphi = calcphireco - calcphimc;
352  if (deltaphi < 0) {
353  deltaphi *= -1;
354  }
355  if (deltaphi > TMath::Pi()) {
356  deltaphi = 2*TMath::Pi() - deltaphi;
357  }
358  Double_t deltapt = recotrackpt.at(r) - mcpt;
359 
360  Int_t chargesign = 1;
361  if (pdg->GetParticle(mctrack->GetPdgCode())->Charge() < 0) {
362  chargesign = -1;
363  }
364 
365  hphires->Fill(deltaphi*chargesign);
366  hptres->Fill(mcpt, deltapt/mcpt);
367 
368  bool phiok = (deltaphi < 0.1);
369 
370  if (deltapt < 0) {
371  deltapt *= -1;
372  }
373  bool ptok = ( (deltapt < mcpt * 0.1) || ((mcpt > 450)&&(recotrackpt.at(r) > 450)) );
374  if (phiok && ptok) {
375  trackstatus++;
376  recotrackstatus.at(r) = 1;
377  //hmcrecoed->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
378  } else {
379  //hmcnotrecoed->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
380  }
381  if ((!phiok) && ptok) {
382  deltaphifailed++;
383  }
384  if (phiok && (!ptok)) {
385  deltaptfailed++;
386  }
387  if ((!phiok) && (!ptok)) {
388  deltabothfailed++;
389  }
390  }
391  }
392  //cout << "Found " << trackstatus << " matching reco tracks." << endl;
393  if (trackstatus == 0) {
394  hmcnotrecoed->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
395  norecotracks++;
396  } else {
397  hmcrecoed->Fill(mctrack->GetMomentum().Pz(), mctrack->GetMomentum().Pt());
398  }
399  mctrackstatus.push_back(trackstatus);
400  goodmctracks++;
401  } else {
402  secondaries++;
403  }
404  }
405  }
406  }
407  }
408 
409  for (int s = 0; s < mctrackstatus.size(); s++) {
410  if (mctrackstatus.at(s) > 0) {
411  goodrecotracks++;
412  }
413  }
414 
415  int recotrackstatusgoodcount = 0;
416  int recotrackstatusbadcount = 0;
417 
418  for (int s = 0; s < recotrackstatus.size(); s++) {
419  if (recotrackmcid.at(s) <= globalmaxmcid) {
420  if (recotrackstatus.at(s) > 0) {
421  recotrackstatusgoodcount++;
422  } else {
423  recotrackstatusbadcount++;
424  }
425  }
426  }
427 
428  //for (std::map<int, std::map<int, int> >::iterator iterevent = foundmap90.begin(); iterevent != foundmap90.end(); iterevent++) {
429  //for (itertrack = *iterevent.begin(); itertrack != *iterevent.end(); itertrack++
430  //}
431 
432  cout << "Total reconstructed track count: " << totalrecotracks << endl;
433  cout << "Reached MC ID: " << globalmaxmcid << endl;
434  cout << "Reconstructable MC track count: " << goodmctracks << endl;
435  cout << "MC Tracks reconstructed: " << goodrecotracks << endl;
436  cout << "MC Tracks NOT reconstructed: " << norecotracks << endl;
437  cout << "Reco Status Good: " << recotrackstatusgoodcount << endl;
438  cout << "Reco Status Bad: " << recotrackstatusbadcount << endl;
439  cout << "Delta phi failed: " << deltaphifailed << endl;
440  cout << "Delta pt failed: " << deltaptfailed << endl;
441  cout << "Delta both failed: " << deltabothfailed << endl;
442  cout << "Secondary Tracks: " << secondaries << endl;
443 
444  cout << "MC Match Good 90: " << mcidgood90 << endl;
445  cout << "MC Match Good 80: " << mcidgood80 << endl;
446  cout << "MC Match Bad: " << mcidbad << endl;
447 
448  //c1->Print("hitdisplay.gif++");
449 
450  TFile outputstorage(histofilename, "UPDATE");
451  outputstorage.cd();
452  SaveAndUpdateHisto(hmcall, outputstorage);
453  SaveAndUpdateHisto(hmcprime, outputstorage);
454  SaveAndUpdateHisto(hmcrecoable, outputstorage);
455  SaveAndUpdateHisto(hdafuq, outputstorage);
456  SaveAndUpdateHisto(hmcrecoed, outputstorage);
457  SaveAndUpdateHisto(hmcnotrecoed, outputstorage);
458  SaveAndUpdateHisto(hphires, outputstorage);
459  SaveAndUpdateHisto(hptres, outputstorage);
460  SaveAndUpdateHisto(hphiresmc, outputstorage);
461  SaveAndUpdateHisto(hptresmc, outputstorage);
462  SaveAndUpdateHisto(hphiresmc90, outputstorage);
463  SaveAndUpdateHisto(hptresmc90, outputstorage);
464  SaveAndUpdateHisto(hphiresmc80, outputstorage);
465  SaveAndUpdateHisto(hptresmc80, outputstorage);
466  SaveAndUpdateHisto(hphiresmcbad, outputstorage);
467  SaveAndUpdateHisto(hptresmcbad, outputstorage);
468 
469  // ----- Finish -------------------------------------------------------
470  timer.Stop();
471  Double_t rtime = timer.RealTime();
472  Double_t ctime = timer.CpuTime();
473  cout << endl << endl;
474  cout << "Macro finished succesfully." << endl;
475  cout << "Real time " << rtime << " s, CPU time " << ctime << " s" << endl;
476  cout << endl;
477  // ------------------------------------------------------------------------
478  return 0;
479 }
PndMCTrack * mctrack
TFile filereco("MvdStt_Test_reco.root")
double r
Definition: RiemannTest.C:14
Int_t i
Definition: run_full.C:25
Int_t GetNPoints(DetectorId detId) const
Definition: PndMCTrack.cxx:120
Int_t GetPdgCode() const
Definition: PndMCTrack.h:73
TLorentzVector s
Definition: Pnd2DStar.C:50
TVector3 GetMomentum() const
Definition: PndMCTrack.h:78
FairRunAna * fRun
Definition: hit_dirc.C:58
static T Abs(const T &x)
Definition: PndCAMath.h:39
Double_t
TTree * treedigi
TStopwatch timer
Definition: hit_dirc.C:51
TFile filedigi("testdigi.root")
FairRuntimeDb * rtdb
Definition: hit_dirc.C:66
c1
Definition: plot_dirc.C:35
Double_t ctime
Definition: hit_dirc.C:114
Int_t iVerbose
TVector3 GetStartVertex() const
Definition: PndMCTrack.h:76
Int_t GetMotherID() const
Definition: PndMCTrack.h:74
Double_t rtime
Definition: hit_dirc.C:113
Double_t Pi
int runOnlineDisplayMCCheck(Int_t maximumTime=2050)
void SaveAndUpdateHisto(TH1 *currenthisto, TFile &storagefile)