CbmRoot
CbmMvdReadoutSimple.cxx
Go to the documentation of this file.
1 // -------------------------------------------------------------------------
2 // ----- CbmMvdReadoutSimple source file -----
3 // ----- Created 17/10/16 by P. Sitzmann -----
4 // -------------------------------------------------------------------------
5 
6 #include "CbmMvdReadoutSimple.h"
7 
8 #include "CbmMCTrack.h"
9 #include "CbmMvdPoint.h"
10 
11 #include "CbmMvdDetector.h"
12 #include "CbmMvdSensor.h"
13 
14 #include "tools/CbmMvdGeoHandler.h"
15 
16 
17 //-- Include from Fair --//
18 #include "FairLogger.h"
19 
20 #include "FairTrackParam.h"
21 
22 
23 //-- Include from Root --//
24 #include "TCanvas.h"
25 #include "TF1.h"
26 #include "TMath.h"
27 #include <TFile.h>
28 
29 #include <iostream>
30 
31 using std::cout;
32 using std::endl;
33 
34 // ----- Default constructor -------------------------------------------
36  : CbmMvdReadoutSimple::CbmMvdReadoutSimple("MvdReadoutSimple", 0) {}
37 // -------------------------------------------------------------------------
38 
39 // ----- Standard constructor ------------------------------------------
40 CbmMvdReadoutSimple::CbmMvdReadoutSimple(const char* name, Int_t iVerbose)
41  : FairTask(name, iVerbose)
42  , foutFile(nullptr)
43  , fshow(kFALSE)
44  , fMvdMCBank()
45  , fMvdMCHitsStations()
46  , fWordsPerRegion(nullptr)
47  , fWordsPerRegion2(nullptr)
48  , fWordsPerWorstRegion(nullptr)
49  , fWordsPerSuperRegion(nullptr)
50  , fWorstSuperPerEvent(nullptr)
51  , fMvdBankDist(nullptr)
52  , fMvdMCWorst(nullptr)
53  , fMvdMCWorstDelta(nullptr)
54  , fMvdDataLoadPerSensor(nullptr)
55  , fMvdDataLoadHotSensor(nullptr)
56  , fMvdDataPerRegion()
57  , fMvdDataPerSuperRegion()
58  , fMcPoints()
59  , fListMCTracks()
60  , fEventNumber(0) {}
61 // -------------------------------------------------------------------------
62 
63 // ----- Destructor ----------------------------------------------------
65 // -------------------------------------------------------------------------
66 
67 // -------------------------------------------------------------------------
69  FairRootManager* ioman = FairRootManager::Instance();
70  if (!ioman) { LOG(fatal) << "RootManager not instantised!"; }
71 
72  fMcPoints = (TClonesArray*) ioman->GetObject("MvdPoint"); // PileUp Mc points
73  fListMCTracks = (TClonesArray*) ioman->GetObject("MCTrack");
74 
75  if (!fMcPoints) LOG(fatal) << "Mvd Pile Up Mc array missing";
76 
78 
79  return kSUCCESS;
80 }
81 // -------------------------------------------------------------------------
82 
83 // -------------------------------------------------------------------------
85  for (Int_t i = 0; i < 63; i++) {
86  fMvdMCBank[i] =
87  new TH2F(Form("fMvdMCBank%d", i),
88  "Mvd mc distribution worst spot only delta electrons",
89  300,
90  -2,
91  0,
92  1500,
93  3.5,
94  0);
95  }
96 
97  fMvdMCHitsStations[0] = new TH2F(
98  "fMvdMCStation0", "Mvd mc distribution", 2, -2.5, -0.5, 3, -0.5, 2.5);
99  fMvdMCHitsStations[1] = new TH2F(
100  "fMvdMCStation1", "Mvd mc distribution", 4, -4.5, -0.5, 6, -0.5, 5.5);
101  fMvdMCHitsStations[2] = new TH2F(
102  "fMvdMCStation2", "Mvd mc distribution", 6, -7.5, -1.5, 9, -1.5, 7.5);
103  fMvdMCHitsStations[3] = new TH2F(
104  "fMvdMCStation3", "Mvd mc distribution", 8, -9.5, -1.5, 12, -1.5, 10.5);
105 
107  new TH1F("fWordsPerRegion", "Words send to a region", 65, 0, 64);
108  fWordsPerRegion->GetXaxis()->SetTitle("Regionnumber");
109  fWordsPerRegion->GetYaxis()->SetTitle("Average Entries per Region");
110 
112  new TH2F("fWordsPerRegion2",
113  "Words send to a region, errors sigma of gauss fit",
114  64,
115  0,
116  63,
117  150,
118  0,
119  15);
120  fWordsPerRegion2->GetXaxis()->SetTitle("Regionnumber");
121  fWordsPerRegion2->GetYaxis()->SetTitle("Average Entries per Region");
122 
123  fWordsPerWorstRegion = new TH1F(
124  "fWordsPerWorstRegion", "Most worst send to a region per Event", 65, 0, 64);
125  fWordsPerWorstRegion->GetXaxis()->SetTitle("words");
126  fWordsPerWorstRegion->GetYaxis()->SetTitle("Entries");
127 
128  fWordsPerSuperRegion = new TH1F(
129  "fWordsPerSuperRegion", "Words send to a super region", 1000, 0, 400);
130  fWordsPerSuperRegion->GetXaxis()->SetTitle("words");
131  fWordsPerSuperRegion->GetYaxis()->SetTitle("Entries");
132 
133  fWorstSuperPerEvent = new TH1F("fWorstSuperRegion",
134  "Most words send to super region per Event",
135  1000,
136  0,
137  400);
138  fWorstSuperPerEvent->GetXaxis()->SetTitle("words");
139  fWorstSuperPerEvent->GetYaxis()->SetTitle("Entries");
140 
141  fMvdMCWorst = new TH2F(
142  "fMvdMCWorst", "Mvd mc distribution worst spot", 300, -2, 0, 1500, 3.5, 0);
143  fMvdMCWorst->GetYaxis()->SetTitle("y [cm]");
144  fMvdMCWorst->GetXaxis()->SetTitle("x [cm]");
145 
147  new TH2F("fMvdMCWorstDelta",
148  "Mvd mc distribution worst spot only delta electrons",
149  300,
150  -2,
151  0,
152  1500,
153  3.5,
154  0);
155  fMvdMCWorstDelta->GetYaxis()->SetTitle("y [cm]");
156  fMvdMCWorstDelta->GetXaxis()->SetTitle("x [cm]");
157 
159  new TH1I("fMvdDataLoadPerSensor", "Mvd Dataload per Sensor", 300, 0, 300);
160  fMvdDataLoadPerSensor->GetXaxis()->SetTitle("Sensor number");
161  fMvdDataLoadPerSensor->GetYaxis()->SetTitle("Entries");
162 
163  fMvdDataLoadHotSensor = new TH1I(
164  "fMvdDataLoadHotSensor", "Mvd Dataload in worst Sensor", 2000, 0, 2000);
165  fMvdDataLoadHotSensor->GetXaxis()->SetTitle("number of words");
166  fMvdDataLoadHotSensor->GetYaxis()->SetTitle("Entries");
167 
168  fMvdBankDist =
169  new TH2I("fMvdBankDist", "Avarage Hits per Region", 63, 0, 63, 50, 0, 50);
170  fMvdBankDist->GetXaxis()->SetTitle("Region number");
171  fMvdBankDist->GetYaxis()->SetTitle("Entries");
172 
173  for (Int_t i = 0; i < 64; i++) {
174  fMvdDataPerRegion[i] = new TH1F(Form("fMvdDataPerRegion[%d]", i),
175  Form("Words send to region %d", i),
176  100,
177  0,
178  100);
179  fMvdDataPerRegion[i]->GetXaxis()->SetTitle("Words");
180  fMvdDataPerRegion[i]->GetYaxis()->SetTitle("Entries");
181  }
182 
183 
184  for (Int_t i = 0; i < 16; i++) {
186  new TH1F(Form("fMvdDataPerSuperRegion[%d]", i),
187  Form("Words send to superregion %d", i),
188  400,
189  0,
190  400);
191  fMvdDataPerSuperRegion[i]->GetXaxis()->SetTitle("Words");
192  fMvdDataPerSuperRegion[i]->GetYaxis()->SetTitle("Entries");
193  }
194 }
195 // -------------------------------------------------------------------------
196 
197 
198 // -------------------------------------------------------------------------
199 void CbmMvdReadoutSimple::Exec(Option_t* /*opt*/) {
200  fEventNumber++;
201 
202  Float_t wordsPerRegion[64] = {0};
203  Float_t wordsPerSuper[16] = {0};
204 
205  Float_t yPosMin = -0.73;
206  Float_t yPosMax = 2.5;
207 
208  for (Int_t k = 0; k < fMcPoints->GetEntriesFast(); k++) {
209  CbmMvdPoint* curPoint = (CbmMvdPoint*) fMcPoints->At(k);
210  fMvdDataLoadPerSensor->Fill(curPoint->GetStationNr(), 1.5);
211 
212  if (curPoint->GetZ() < 8) {
213  if (curPoint->GetX() > -1.93 && curPoint->GetX() <= -0.55
214  && curPoint->GetY() >= yPosMin && curPoint->GetY() <= yPosMax) {
215  fMvdMCWorst->Fill(curPoint->GetX(), curPoint->GetY());
216  if (curPoint->GetTrackID() == -3)
217  fMvdMCWorstDelta->Fill(curPoint->GetX(), curPoint->GetY());
218  for (Int_t nRegion = 0; nRegion < 64; nRegion++) {
219  if (curPoint->GetY() >= (yPosMin + (nRegion * 0.05))
220  && curPoint->GetY() < (yPosMin + ((nRegion + 1) * 0.05))) {
221  fMvdMCBank[nRegion]->Fill(curPoint->GetX(), curPoint->GetY());
222  wordsPerRegion[nRegion] = wordsPerRegion[nRegion] + 1.5;
223  if (wordsPerRegion[nRegion] > 100)
224  LOG(info) << " Region " << nRegion << " has an overflow";
225  fWordsPerRegion->Fill(nRegion, 1.5);
226  break;
227  }
228  }
229  }
230  }
231  if (curPoint->GetZ() < 8) {
232  if (curPoint->GetX() < -0.5 && curPoint->GetY() > -0.5)
233  fMvdMCHitsStations[0]->Fill(curPoint->GetX(), curPoint->GetY());
234  } else if (curPoint->GetZ() < 13) {
235  if (curPoint->GetX() < -0.5 && curPoint->GetY() > -1.5)
236  fMvdMCHitsStations[1]->Fill(curPoint->GetX(), curPoint->GetY());
237  } else if (curPoint->GetZ() < 18) {
238  if (curPoint->GetX() < -1.5 && curPoint->GetY() > -1.5)
239  fMvdMCHitsStations[2]->Fill(curPoint->GetX(), curPoint->GetY());
240  } else {
241  if (curPoint->GetX() < -1.5 && curPoint->GetY() > -1.5)
242  fMvdMCHitsStations[3]->Fill(curPoint->GetX(), curPoint->GetY());
243  }
244  }
245  LOG(info) << "//--------------- New Event " << fEventNumber
246  << " -----------------------\\";
247 
248  Int_t i = 0;
249  Int_t wordsInWorst = 0;
250  Int_t wordsInWorstReg = 0;
251  Int_t wordsTotal = 0;
252 
253  for (Int_t supReg = 0; supReg < 16; supReg++) {
254  for (Int_t k = 0; k < 4; k++) {
255  wordsPerSuper[supReg] += wordsPerRegion[i];
256  fMvdDataPerRegion[i]->Fill(wordsPerRegion[i]);
257  if (wordsPerRegion[i] > wordsInWorstReg)
258  wordsInWorstReg = wordsPerRegion[i];
259  LOG(debug) << "Words in Region " << i << ": " << wordsPerRegion[i];
260  i++;
261  }
262 
263  LOG(debug) << " Words in super region " << supReg << ": "
264  << wordsPerSuper[supReg];
265  if (wordsPerSuper[supReg] > 400)
266  LOG(info) << "SuperRegion " << supReg << " has an overflow";
267  fWordsPerSuperRegion->Fill(wordsPerSuper[supReg]);
268  fMvdDataPerSuperRegion[supReg]->Fill(wordsPerSuper[supReg]);
269 
270  if (wordsPerSuper[supReg] > wordsInWorst)
271  wordsInWorst = wordsPerSuper[supReg];
272 
273  wordsTotal += wordsPerSuper[supReg];
274  }
275 
276  fWorstSuperPerEvent->Fill(wordsInWorst);
277  fWordsPerWorstRegion->Fill(wordsInWorstReg);
278  fMvdDataLoadHotSensor->Fill(wordsTotal);
279 
280  LOG(info) << "//--------------- End Event -----------------------\\";
281 }
282 // -------------------------------------------------------------------------
283 
284 // -------------------------------------------------------------------------
286  foutFile->cd();
287 
288  Float_t scale = 1. / (Float_t) fEventNumber;
289  for (Int_t iPad = 0; iPad < 4; iPad++) {
290  fMvdMCHitsStations[iPad]->Scale(scale);
291  }
292  fMvdDataLoadPerSensor->Scale(scale);
293  fWordsPerRegion->Scale(scale);
294 
295 
296  if (fshow) DrawHistograms();
297  WriteHistograms();
298 }
299 // -------------------------------------------------------------------------
300 
301 // -------------------------------------------------------------------------
303 
304  for (Int_t iPad = 0; iPad < 4; iPad++) {
305  fMvdMCHitsStations[iPad]->Write();
306  }
307 
308  fMvdDataLoadPerSensor->Write();
309 
310  fMvdDataLoadHotSensor->Write();
311 
312  fWordsPerSuperRegion->Write();
313 
314  fWorstSuperPerEvent->Write();
315 
316  fWordsPerRegion->Write();
317 
318  fWordsPerWorstRegion->Write();
319 
320  for (Int_t i = 0; i < 64; i++) {
321  fMvdDataPerRegion[i]->Fit("gaus");
322  TF1* func = fMvdDataPerRegion[i]->GetFunction("gaus");
323  Double_t param0 = func->GetParameter(0);
324  Double_t param1 = func->GetParameter(1);
325  Double_t param2 = func->GetParameter(2);
326  Double_t chi2 = func->GetChisquare();
327  cout << " // - " << i << " -- param 0 = " << param0
328  << " -- param 1 = " << param1 << " -- param 2 = " << param2
329  << " -- chi2 = " << chi2 << endl;
330  fWordsPerRegion2->Fill(i, param1);
331  fWordsPerRegion2->SetBinError(i, 0, param1, param2);
332  fMvdDataPerRegion[i]->Write();
333  }
334  fWordsPerRegion2->Write();
335 
336  for (Int_t i = 0; i < 16; i++)
337  fMvdDataPerSuperRegion[i]->Write();
338 }
339 // -------------------------------------------------------------------------
340 
341 // -------------------------------------------------------------------------
343  TCanvas* mcCanvas1 = new TCanvas();
344  mcCanvas1->Divide(2, 2);
345  for (Int_t iPad = 0; iPad < 4; iPad++) {
346  mcCanvas1->cd(iPad + 1);
347  fMvdMCHitsStations[iPad]->Draw("COLZ");
348  }
349 
350  TCanvas* dataLoad = new TCanvas();
351  dataLoad->Divide(2, 1);
352  dataLoad->cd(1);
353  fMvdDataLoadPerSensor->Draw("BAR");
354  dataLoad->cd(2);
355  fMvdDataLoadHotSensor->Draw("BAR");
356 
357  TCanvas* mcCanvas2 = new TCanvas();
358  mcCanvas2->Divide(1, 2);
359  mcCanvas2->cd(1);
360  fWordsPerSuperRegion->Draw();
361  mcCanvas2->cd(2);
362  fWorstSuperPerEvent->Draw();
363 
364  /*
365 TCanvas* mcCanvas3 = new TCanvas();
366 fMvdMCHitsStations[0]->Draw("COLZ");
367 
368 TCanvas* mcCanvas4 = new TCanvas();
369 mcCanvas4->Divide(8,8);
370 for(Int_t pad = 0; pad < 63; pad++)
371  {
372  mcCanvas4->cd(pad+1);
373  fMvdMCBank[pad]->Draw("COL");
374  fMvdMCBank[pad]->Write();
375  cout << "Bank " << pad << " avarage Entries " << fMvdMCBank[pad]->GetEntries()/fEventNumber << endl;
376  fMvdBankDist->Fill(pad, fMvdMCBank[pad]->GetEntries()/fEventNumber);
377  }
378 TCanvas* mcCanvas5 = new TCanvas();
379 fMvdBankDist->Draw();
380 
381 TCanvas* mcCanvas6 = new TCanvas();
382 fMvdMCWorst->Draw("COL");
383 fMvdMCWorst->Write();
384 */
385 
386  TCanvas* regionCanvas = new TCanvas();
387  regionCanvas->Divide(1, 2);
388  regionCanvas->cd(1);
389  fWordsPerRegion->Draw();
390  regionCanvas->cd(2);
391  fWordsPerRegion2->Draw();
392 
393  TCanvas* regionsCanvas[64];
394  TCanvas* supregionsCanvas[16];
395 
396  for (Int_t i = 0; i < 64; i++) {
397  regionsCanvas[i] = new TCanvas();
398  regionsCanvas[i]->cd();
399  fMvdDataPerRegion[i]->Draw();
400  }
401 
402  for (Int_t i = 0; i < 16; i++) {
403  supregionsCanvas[i] = new TCanvas();
404  supregionsCanvas[i]->cd();
405  fMvdDataPerSuperRegion[i]->Draw();
406  }
407 }
408 // -------------------------------------------------------------------------
409 
CbmMvdReadoutSimple::fMvdMCHitsStations
TH2F * fMvdMCHitsStations[4]
Definition: CbmMvdReadoutSimple.h:45
CbmMvdDetector.h
CbmMvdReadoutSimple::fWordsPerWorstRegion
TH1F * fWordsPerWorstRegion
Definition: CbmMvdReadoutSimple.h:48
CbmMvdPoint::GetStationNr
Int_t GetStationNr() const
Definition: CbmMvdPoint.h:78
ClassImp
ClassImp(CbmMvdReadoutSimple)
CbmMvdReadoutSimple::fListMCTracks
TClonesArray * fListMCTracks
Definition: CbmMvdReadoutSimple.h:60
CbmMvdReadoutSimple::WriteHistograms
void WriteHistograms()
Definition: CbmMvdReadoutSimple.cxx:302
CbmMvdReadoutSimple::Init
InitStatus Init()
Definition: CbmMvdReadoutSimple.cxx:68
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmMvdSensor.h
CbmMvdReadoutSimple::fMvdDataLoadHotSensor
TH1I * fMvdDataLoadHotSensor
Definition: CbmMvdReadoutSimple.h:55
CbmMvdReadoutSimple::fWordsPerSuperRegion
TH1F * fWordsPerSuperRegion
Definition: CbmMvdReadoutSimple.h:49
CbmMvdReadoutSimple::Finish
void Finish()
Definition: CbmMvdReadoutSimple.cxx:285
CbmMvdReadoutSimple
Definition: CbmMvdReadoutSimple.h:21
CbmMvdReadoutSimple::fEventNumber
Int_t fEventNumber
Definition: CbmMvdReadoutSimple.h:62
CbmMvdReadoutSimple::fshow
Bool_t fshow
Definition: CbmMvdReadoutSimple.h:43
CbmMvdReadoutSimple::Exec
void Exec(Option_t *opt)
Definition: CbmMvdReadoutSimple.cxx:199
CbmMvdReadoutSimple::fMvdBankDist
TH2I * fMvdBankDist
Definition: CbmMvdReadoutSimple.h:51
CbmMvdReadoutSimple::SetupHistograms
void SetupHistograms()
Definition: CbmMvdReadoutSimple.cxx:84
CbmMvdReadoutSimple::fWordsPerRegion
TH1F * fWordsPerRegion
Definition: CbmMvdReadoutSimple.h:46
CbmMvdGeoHandler.h
Helper class to extract information from the GeoManager. Addapted from TrdGeoHandler byFlorian Uhlig ...
CbmMvdPoint.h
CbmMvdPoint
Definition: CbmMvdPoint.h:28
CbmMvdReadoutSimple::fMvdMCWorstDelta
TH2F * fMvdMCWorstDelta
Definition: CbmMvdReadoutSimple.h:53
CbmMvdReadoutSimple.h
CbmMvdReadoutSimple::fMcPoints
TClonesArray * fMcPoints
Definition: CbmMvdReadoutSimple.h:59
CbmMvdReadoutSimple::fMvdDataPerRegion
TH1F * fMvdDataPerRegion[64]
Definition: CbmMvdReadoutSimple.h:56
CbmMvdReadoutSimple::DrawHistograms
void DrawHistograms()
Definition: CbmMvdReadoutSimple.cxx:342
CbmMvdReadoutSimple::foutFile
TFile * foutFile
Definition: CbmMvdReadoutSimple.h:41
CbmMvdReadoutSimple::fWordsPerRegion2
TH2F * fWordsPerRegion2
Definition: CbmMvdReadoutSimple.h:47
CbmMvdReadoutSimple::fMvdDataLoadPerSensor
TH1I * fMvdDataLoadPerSensor
Definition: CbmMvdReadoutSimple.h:54
CbmMvdReadoutSimple::fMvdMCBank
TH2F * fMvdMCBank[63]
Definition: CbmMvdReadoutSimple.h:44
CbmMCTrack.h
CbmMvdReadoutSimple::CbmMvdReadoutSimple
CbmMvdReadoutSimple()
Definition: CbmMvdReadoutSimple.cxx:35
CbmMvdReadoutSimple::fMvdDataPerSuperRegion
TH1F * fMvdDataPerSuperRegion[16]
Definition: CbmMvdReadoutSimple.h:57
CbmMvdReadoutSimple::fMvdMCWorst
TH2F * fMvdMCWorst
Definition: CbmMvdReadoutSimple.h:52
CbmMvdReadoutSimple::fWorstSuperPerEvent
TH1F * fWorstSuperPerEvent
Definition: CbmMvdReadoutSimple.h:50
CbmMvdReadoutSimple::~CbmMvdReadoutSimple
~CbmMvdReadoutSimple()
Definition: CbmMvdReadoutSimple.cxx:64