CbmRoot
CbmTrdOccupancyQa.cxx
Go to the documentation of this file.
1 #include "CbmTrdOccupancyQa.h"
2 
3 #include "CbmTrdAddress.h"
4 #include "CbmTrdCluster.h"
5 #include "CbmTrdDigi.h"
6 #include "CbmTrdGeoHandler.h"
7 #include "CbmTrdParModDigi.h"
8 #include "CbmTrdParModGeo.h"
9 #include "CbmTrdParSetDigi.h"
10 #include "CbmTrdParSetGeo.h"
11 
12 #include <FairRootManager.h>
13 #include <FairRunAna.h>
14 #include <FairRuntimeDb.h>
15 
16 #include "CbmDigiManager.h"
17 #include <TBox.h>
18 #include <TCanvas.h>
19 #include <TClonesArray.h>
20 #include <TColor.h>
21 #include <TF1.h>
22 #include <TFile.h>
23 #include <TGeoManager.h>
24 #include <TH1.h>
25 #include <TH2.h>
26 #include <TMath.h>
27 #include <TPaveText.h>
28 #include <TProfile.h>
29 #include <TProfile2D.h>
30 #include <TStopwatch.h>
31 
32 
33 #include <cmath>
34 #include <iostream>
35 
36 using std::cin;
37 using std::cout;
38 using std::endl;
39 using std::fabs;
40 using std::pair;
41 
42 CbmTrdOccupancyQa::CbmTrdOccupancyQa()
43  : CbmTrdOccupancyQa("TrdOccupancy", "", "", 1e-6, kFALSE) {}
44 
45 CbmTrdOccupancyQa::CbmTrdOccupancyQa(const char* name,
46  const char*,
47  const char* geo,
48  Double_t triggerThreshold,
49  Bool_t plotMergedResults)
50  : FairTask(name)
51  , fClusters(NULL)
52  , fDigiPar(NULL)
53  , fGeoPar(NULL)
54  , fGeoHandler(new CbmTrdGeoHandler())
55  , fDigiChargeSpectrum(
56  new TH1I("DigiChargeSpectrum", "DigiChargeSpectrum", 1e6, 0, 1.0e-3))
57  , fLayerDummy(new TH2I("LayerDummy", "", 1200, -400, 400, 1000, -300, 300))
58  , fmin(0)
59  , fmax(13)
60  , fModuleMap()
61  , fModuleMapIt()
62  , fModuleOccupancyMap()
63  , fModuleOccupancyMapIt()
64  , fModuleOccupancyMemoryMap()
65  , fModuleOccupancyMemoryMapIt()
66  , fLayerOccupancyMap()
67  , fLayerOccupancyMapIt()
68  , fLayerAverageOccupancyMap()
69  , fLayerAverageOccupancyMapIt()
70  , fTriggerThreshold(triggerThreshold)
71  , fNeigbourReadout(true)
72  , fPlotMergedResults(plotMergedResults)
73  , fGeo(geo) {
74  fLayerDummy->SetXTitle("x-coordinate [cm]");
75  fLayerDummy->SetYTitle("y-coordinate [cm]");
76  fLayerDummy->GetXaxis()->SetLabelSize(0.02);
77  fLayerDummy->GetYaxis()->SetLabelSize(0.02);
78  fLayerDummy->GetZaxis()->SetLabelSize(0.02);
79  fLayerDummy->GetXaxis()->SetTitleSize(0.02);
80  fLayerDummy->GetXaxis()->SetTitleOffset(1.5);
81  fLayerDummy->GetYaxis()->SetTitleSize(0.02);
82  fLayerDummy->GetYaxis()->SetTitleOffset(2);
83  fLayerDummy->GetZaxis()->SetTitleSize(0.02);
84  fLayerDummy->GetZaxis()->SetTitleOffset(-2);
85  fLayerDummy->SetContour(99);
86  fLayerDummy->GetZaxis()->SetRangeUser(fmin, fmax);
87  for (fModuleOccupancyMapIt = fModuleOccupancyMap.begin();
88  fModuleOccupancyMapIt != fModuleOccupancyMap.end();
89  ++fModuleOccupancyMapIt) {
90  //delete fModuleOccupancyMapIt->second;
91  }
92 }
93 
94 CbmTrdOccupancyQa::~CbmTrdOccupancyQa() {
95 
96  FairRootManager* ioman = FairRootManager::Instance();
97  ioman->Write();
98  //CreateLayerView();
99  //SaveHistos2File();
100  fClusters->Delete();
101  delete fClusters;
102 
103  if (fDigiPar) delete fDigiPar;
104  // if(fModuleInfo)
105  // delete fModuleInfo;
106  for (fModuleMapIt = fModuleMap.begin(); fModuleMapIt != fModuleMap.end();
107  ++fModuleMapIt) {
108  delete fModuleMapIt->second;
109  }
110  fModuleMap.clear();
111  for (fLayerOccupancyMapIt = fLayerOccupancyMap.begin();
112  fLayerOccupancyMapIt != fLayerOccupancyMap.end();
113  ++fLayerOccupancyMapIt) {
114  delete fLayerOccupancyMapIt->second;
115  }
116  fLayerOccupancyMap.clear();
117  for (fModuleOccupancyMapIt = fModuleOccupancyMap.begin();
118  fModuleOccupancyMapIt != fModuleOccupancyMap.end();
119  ++fModuleOccupancyMapIt) {
120  delete fModuleOccupancyMapIt->second;
121  }
122  fModuleOccupancyMap.clear();
123  for (fModuleOccupancyMemoryMapIt = fModuleOccupancyMemoryMap.begin();
124  fModuleOccupancyMemoryMapIt != fModuleOccupancyMemoryMap.end();
125  ++fModuleOccupancyMemoryMapIt) {
126  delete fModuleOccupancyMemoryMapIt->second;
127  }
128  fModuleOccupancyMemoryMap.clear();
129  delete fDigiChargeSpectrum;
130  delete fLayerDummy;
131 }
132 void CbmTrdOccupancyQa::SetParContainers() {
133  cout << " * CbmTrdOccupancyQa * :: SetParContainers() " << endl;
134  // Get Base Container
135  FairRunAna* ana = FairRunAna::Instance();
136  FairRuntimeDb* rtdb = ana->GetRuntimeDb();
137  fDigiPar = (CbmTrdParSetDigi*) (rtdb->getContainer("CbmTrdParSetDigi"));
138  fGeoPar = (CbmTrdParSetGeo*) (rtdb->getContainer("CbmTrdParSetGeo"));
139 }
140 // ---- ReInit -------------------------------------------------------
141 InitStatus CbmTrdOccupancyQa::ReInit() {
142  cout << " * CbmTrdClusterizer * :: ReInit() " << endl;
143  FairRunAna* ana = FairRunAna::Instance();
144  FairRuntimeDb* rtdb = ana->GetRuntimeDb();
145  fDigiPar = (CbmTrdParSetDigi*) (rtdb->getContainer("CbmTrdParSetDigi"));
146  fGeoPar = (CbmTrdParSetGeo*) (rtdb->getContainer("CbmTrdParSetGeo"));
147  return kSUCCESS;
148 }
149 // --------------------------------------------------------------------
150 
151 // ---- Init ----------------------------------------------------------
152 InitStatus CbmTrdOccupancyQa::Init() {
153  cout << " * CbmTrdOccupancyQa * :: Init() " << endl;
154  FairRootManager* ioman = FairRootManager::Instance();
155 
156  CbmDigiManager::Instance()->Init();
157  if (!CbmDigiManager::Instance()->IsPresent(ECbmModuleId::kTrd)) LOG(fatal);
158 
159  fClusters = (TClonesArray*) ioman->GetObject("TrdCluster");
160  if (!fClusters) {
161  cout << "-W CbmTrdHitDensityQa::Init: No TrdCluster array!" << endl;
162  cout << " Task will be inactive" << endl;
163  return kERROR;
164  }
165  fGeoHandler->Init();
166 
167  return kSUCCESS;
168 }
169 // --------------------------------------------------------------------
170 void CbmTrdOccupancyQa::SetNeighbourTrigger(Bool_t trigger) {
171  fNeigbourReadout = trigger;
172 }
173 
174 // ---- Exec ----------------------------------------------------------
175 void CbmTrdOccupancyQa::Exec(Option_t*) {
176  printf("================CbmTrdOccupancyQa=====================\n");
177  //fTriggerThreshold = CbmTrdClusterFinderFast::GetTriggerThreshold();
178  printf("TriggerThreshold: %.2E\n", fTriggerThreshold);
179  printf("NeigbourReadout:%i\n", Int_t(fNeigbourReadout));
180  // Bool_t debug = false;
181  //TFile *outFile = new TFile("data/CbmTrdOccupancyQa.root","UPDATE","output of CbmTrdOccupancyQa");
182  TStopwatch timer;
183  timer.Start();
184  Int_t nEntries = CbmDigiManager::Instance()->GetNofDigis(ECbmModuleId::kTrd);
185  Int_t digiCounter(0), digiTriggerCounter(0);
186  TString title;
187  /*
188  fLayerDummy = new TH2I("LayerDummy","",1200,-600,600,1000,-500,500);
189  fLayerDummy->SetXTitle("x-coordinate [cm]");
190  fLayerDummy->SetYTitle("y-coordinate [cm]");
191  //fDigiChargeSpectrum = (TH1I*)outFile->Get("DigiChargeSpectrum");
192  //if (!fDigiChargeSpectrum)
193  fDigiChargeSpectrum = new TH1I("DigiChargeSpectrum","DigiChargeSpectrum",1e6,0,1e-3);
194  */
195  const CbmTrdDigi* digi = NULL;
196  CbmTrdCluster* cluster = NULL;
197  if (fNeigbourReadout == true) {
198  Int_t nCluster = fClusters->GetEntries();
199  for (Int_t iCluster = 0; iCluster < nCluster; iCluster++) {
200  //cout << iCluster << endl;
201  cluster = (CbmTrdCluster*) fClusters->At(
202  iCluster); //pointer to the acvit cluster
203  Int_t nDigisInCluster = cluster->GetNofDigis();
204  for (Int_t iDigi = 0; iDigi < nDigisInCluster; iDigi++) {
205  digi =
206  CbmDigiManager::Instance()->Get<CbmTrdDigi>(cluster->GetDigi(iDigi));
207  //digi = (CbmTrdDigi*)fDigis->At(cluster->GetDigi(iDigi));
208  if (digi->GetCharge() > fTriggerThreshold)
209  digiTriggerCounter++;
210  else
211  digiCounter++;
212  Int_t digiAddress = digi->GetAddress();
213  Int_t moduleAddress = CbmTrdAddress::GetModuleAddress(digiAddress);
214  // Int_t moduleId = CbmTrdAddress::GetModuleId(moduleAddress);
215  Int_t Station = CbmTrdAddress::GetLayerId(moduleAddress) / 4
216  + 1; //fGeoHandler->GetStation(moduleId);
217  Int_t Layer = CbmTrdAddress::GetLayerId(moduleAddress) % 4
218  + 1; //fGeoHandler->GetLayer(moduleId);
219  Int_t combiId = 10 * Station + Layer;
220  CbmTrdParModDigi* fModuleInfo =
221  (CbmTrdParModDigi*) fDigiPar->GetModulePar(moduleAddress);
222  CbmTrdParModGeo* fModuleGeo =
223  (CbmTrdParModGeo*) fGeoPar->GetModulePar(moduleAddress);
224  Int_t nRows = fModuleInfo->GetNofRows();
225  Int_t nCols = fModuleInfo->GetNofColumns();
226  if (fModuleOccupancyMap.find(moduleAddress)
227  == fModuleOccupancyMap.end()) {
228  title.Form("Module_%i", moduleAddress);
229  //fModuleOccupancyMap[moduleAddress] = (TH2F*)outFile->Get(title);
230  //if (!fModuleOccupancyMap[moduleAddress])
231  fModuleOccupancyMap[moduleAddress] = new TH2I(
232  title, title, nCols, -0.5, nCols - 0.5, nRows, -0.5, nRows - 0.5);
233  //fModuleOccupancyMap[moduleAddress]->Reset();
234  title.Form("M_%i", moduleAddress);
235  fModuleOccupancyMemoryMap[moduleAddress] =
236  new TH1F(title, title, 10000, 0, 100);
237  fModuleMap[moduleAddress] = new OccupancyModule();
238  fModuleMap[moduleAddress]->Station = Station;
239  fModuleMap[moduleAddress]->Layer = Layer;
240  fModuleMap[moduleAddress]->ModuleSizeX = fModuleInfo->GetSizeX();
241  fModuleMap[moduleAddress]->ModuleSizeY = fModuleInfo->GetSizeY();
242  fModuleMap[moduleAddress]->ModulePositionX = fModuleGeo->GetX();
243  fModuleMap[moduleAddress]->ModulePositionY = fModuleGeo->GetY();
244  fModuleMap[moduleAddress]->moduleAddress = moduleAddress;
245  }
246 
247  if (fLayerOccupancyMap.find(combiId) == fLayerOccupancyMap.end()) {
248  title.Form("S%i_L%i", Station, Layer);
249  printf(" new %s\n", title.Data());
250  //fLayerOccupancyMap[combiId] = new TH2F(title,title,1200,-600,600,1000,-500,500);
251  fLayerOccupancyMap[combiId] = new TCanvas(title, title, 1200, 1000);
252  fLayerOccupancyMap[combiId]->cd();
253  fLayerDummy->DrawCopy("colz");
254  title.Form("Station%i_Layer%i", Station, Layer);
255  printf(" new %s\n", title.Data());
256  fLayerAverageOccupancyMap[combiId] =
257  new TProfile(title, title, 1e6, 0, 1.0e-3);
258  fLayerAverageOccupancyMap[combiId]->SetYTitle(
259  "Average layer occupancy [%]");
260  fLayerAverageOccupancyMap[combiId]->SetXTitle("Trigger threshold");
261  }
262 
263  Int_t iCol(CbmTrdAddress::GetColumnId(digiAddress)),
264  iRow(CbmTrdAddress::GetRowId(digiAddress)),
265  iSec(CbmTrdAddress::GetSectorId(digiAddress));
266  iRow = fModuleInfo->GetModuleRow(iSec, iRow);
267  Int_t ixBin(iCol + 1), iyBin(iRow + 1);
268  if (fModuleOccupancyMap[moduleAddress]->GetBinContent(ixBin, iyBin)
269  == 0)
270  fModuleOccupancyMap[moduleAddress]->Fill(iCol, iRow);
271  }
272  }
273  } else {
274  for (Int_t iDigi = 0; iDigi < nEntries; iDigi++) {
275  digi = CbmDigiManager::Instance()->Get<CbmTrdDigi>(iDigi);
276  // digi = (CbmTrdDigi*) fDigis->At(iDigi);
277 
278  Int_t moduleAddress = CbmTrdAddress::GetModuleAddress(
279  digi->GetAddress()); //digi->GetDetId();
280  Int_t Station = CbmTrdAddress::GetLayerId(moduleAddress) / 4
281  + 1; //fGeoHandler->GetStation(moduleId);
282  Int_t Layer = CbmTrdAddress::GetLayerId(moduleAddress) % 4
283  + 1; //fGeoHandler->GetLayer(moduleId);
284  Int_t combiId = 10 * Station + Layer;
285  //printf("Station %i Layer %i combiId %i\n", Station, Layer, combiId);
286  CbmTrdParModDigi* fModuleInfo =
287  (CbmTrdParModDigi*) fDigiPar->GetModulePar(moduleAddress);
288  CbmTrdParModGeo* fModuleGeo =
289  (CbmTrdParModGeo*) fGeoPar->GetModulePar(moduleAddress);
290  Int_t nRows = fModuleInfo->GetNofRows();
291  Int_t nCols = fModuleInfo->GetNofColumns();
292  fDigiChargeSpectrum->Fill(digi->GetCharge());
293  if (digi->GetCharge() > fTriggerThreshold) {
294  digiTriggerCounter++;
295  digiCounter++;
296  if (fModuleOccupancyMap.find(moduleAddress)
297  == fModuleOccupancyMap.end()) {
298  title.Form("Module_%i", moduleAddress);
299  //fModuleOccupancyMap[moduleAddress] = (TH2F*)outFile->Get(title);
300  //if (!fModuleOccupancyMap[moduleAddress])
301  fModuleOccupancyMap[moduleAddress] = new TH2I(
302  title, title, nCols, -0.5, nCols - 0.5, nRows, -0.5, nRows - 0.5);
303  //fModuleOccupancyMap[moduleAddress]->Reset();
304  title.Form("M_%i", moduleAddress);
305  fModuleOccupancyMemoryMap[moduleAddress] =
306  new TH1F(title, title, 10000, 0, 100);
307  fModuleMap[moduleAddress] = new OccupancyModule();
308  fModuleMap[moduleAddress]->Station = Station;
309  fModuleMap[moduleAddress]->Layer = Layer;
310  fModuleMap[moduleAddress]->ModuleSizeX = fModuleInfo->GetSizeX();
311  fModuleMap[moduleAddress]->ModuleSizeY = fModuleInfo->GetSizeY();
312  fModuleMap[moduleAddress]->ModulePositionX = fModuleGeo->GetX();
313  fModuleMap[moduleAddress]->ModulePositionY = fModuleGeo->GetY();
314  fModuleMap[moduleAddress]->moduleAddress = moduleAddress;
315  }
316  if (fLayerOccupancyMap.find(combiId) == fLayerOccupancyMap.end()) {
317  title.Form("S%i_L%i", Station, Layer);
318  printf(" new %s\n", title.Data());
319  //fLayerOccupancyMap[combiId] = new TH2F(title,title,1200,-600,600,1000,-500,500);
320  fLayerOccupancyMap[combiId] = new TCanvas(title, title, 1200, 1000);
321  fLayerOccupancyMap[combiId]->cd();
322  fLayerDummy->DrawCopy("colz");
323  title.Form("Station%i_Layer%i", Station, Layer);
324  printf(" new %s\n", title.Data());
325  fLayerAverageOccupancyMap[combiId] =
326  new TProfile(title, title, 1e6, 0, 1.0e-3);
327  fLayerAverageOccupancyMap[combiId]->SetYTitle(
328  "Average layer occupancy [%]");
329  fLayerAverageOccupancyMap[combiId]->SetXTitle("Trigger threshold");
330  }
331 
332  Int_t digiAddress(digi->GetAddress()),
333  iCol(CbmTrdAddress::GetColumnId(digiAddress)),
334  iRow(CbmTrdAddress::GetRowId(digiAddress)),
335  iSec(CbmTrdAddress::GetSectorId(digiAddress));
336  iRow = fModuleInfo->GetModuleRow(iSec, iRow);
337  Int_t ixBin(iCol + 1), iyBin(iRow + 1);
338  if (fModuleOccupancyMap[moduleAddress]->GetBinContent(ixBin, iyBin)
339  == 0)
340  fModuleOccupancyMap[moduleAddress]->Fill(iCol, iRow);
341  /*
342  if(fNeigbourReadout) {
343  if (nCols > nRows) { // Neigbour readout does not depend on a signal on the neighbouring channel
344  if (ixBin > 1 && fModuleOccupancyMap[moduleAddress]->GetBinContent(ixBin-1,iyBin) == 0){//left
345  fModuleOccupancyMap[moduleAddress]->Fill(iCol-1,iRow);
346  digiCounter++;
347  }
348  if (ixBin < nCols && fModuleOccupancyMap[moduleAddress]->GetBinContent(ixBin+1,iyBin) == 0){//right
349  fModuleOccupancyMap[moduleAddress]->Fill(iCol+1,iRow);
350  digiCounter++;
351  }
352  }
353  else {
354  if (iyBin < nRows && fModuleOccupancyMap[moduleAddress]->GetBinContent(ixBin,iyBin+1) == 0){//upper
355  fModuleOccupancyMap[moduleAddress]->Fill(iCol,iRow+1);
356  digiCounter++;
357  }
358  if (iyBin > 1 && fModuleOccupancyMap[moduleAddress]->GetBinContent(ixBin,iyBin-1) == 0){//lower
359  fModuleOccupancyMap[moduleAddress]->Fill(iCol,iRow-1);
360  digiCounter++;
361  }
362  }
363  }
364  */
365  //fModuleOccupancyMap[moduleAddress]->SetBinContent(digi->GetCol()+1, digi->GetRow()+1, 1);
366  }
367  }
368  }
369  printf(
370  "%6i primary triggered digis\n%6i including neigbouring triggered digis\n",
371  digiTriggerCounter,
372  digiCounter);
373  CopyEvent2MemoryMap();
374  //CreateLayerView();
375  //SaveHistos2File();
376  timer.Stop();
377  Double_t rtime = timer.RealTime();
378  Double_t ctime = timer.CpuTime();
379 
380  printf("\n\n******************** Reading Test **********************\n");
381  printf(" RealTime=%f seconds, CpuTime=%f seconds\n", rtime, ctime);
382  printf("*********************************************************\n\n");
383 }
384 // ---- FinishTask-----------------------------------------------------
385 void CbmTrdOccupancyQa::FinishEvent() {
386  if (fClusters) {
387  fClusters->Clear("C");
388  fClusters->Delete();
389  }
390  for (fModuleOccupancyMapIt = fModuleOccupancyMap.begin();
391  fModuleOccupancyMapIt != fModuleOccupancyMap.end();
392  ++fModuleOccupancyMapIt) {
393  fModuleOccupancyMapIt->second->Reset();
394  }
395 }
396 void CbmTrdOccupancyQa::FinishTask() {
397  CreateLayerView();
398  SaveHistos2File();
399 }
400 // ---- Register ------------------------------------------------------
401 void CbmTrdOccupancyQa::Register() {
402  //FairRootManager::Instance()->Register("TrdDigi","Trd Digi", fDigiCollection, IsOutputBranchPersistent("TrdDigi"));
403  //FairRootManager::Instance()->Register("TrdDigiMatch","Trd Digi Match", fDigiMatchCollection, IsOutputBranchPersistent("TrdDigiMatch"));
404 }
405 // --------------------------------------------------------------------
406 void CbmTrdOccupancyQa::CopyEvent2MemoryMap() {
407  for (fModuleOccupancyMapIt = fModuleOccupancyMap.begin();
408  fModuleOccupancyMapIt != fModuleOccupancyMap.end();
409  ++fModuleOccupancyMapIt) {
410  fModuleOccupancyMemoryMap[fModuleOccupancyMapIt->first]->Fill(
411  Float_t(fModuleOccupancyMapIt->second->Integral())
412  / Float_t(fModuleOccupancyMapIt->second->GetNbinsX()
413  * fModuleOccupancyMapIt->second->GetNbinsY())
414  * 100.);
415  fLayerAverageOccupancyMap
416  [10 * (fModuleMap[fModuleOccupancyMapIt->first]->Station)
417  + (fModuleMap[fModuleOccupancyMapIt->first]->Layer)]
418  ->Fill(fTriggerThreshold,
419  Float_t(fModuleOccupancyMapIt->second->Integral())
420  / Float_t(fModuleOccupancyMapIt->second->GetNbinsX()
421  * fModuleOccupancyMapIt->second->GetNbinsY())
422  * 100.);
423  }
424 }
425 void CbmTrdOccupancyQa::SwitchToMergedFile() {
426  TString filename = "data/result.root";
427  TFile* Target = new TFile(filename, "READ");
428  if (Target) {
429  gDirectory = Target->CurrentDirectory();
430  for (fModuleOccupancyMapIt = fModuleOccupancyMap.begin();
431  fModuleOccupancyMapIt != fModuleOccupancyMap.end();
432  ++fModuleOccupancyMapIt) {
433  TString histName = fModuleOccupancyMapIt->second->GetName();
434  std::cout << histName << std::endl;
435  ;
436  fModuleOccupancyMap[fModuleOccupancyMapIt->first] =
437  (TH2I*) Target->Get("TrdOccupancy/Module/" + histName)
438  ->Clone(histName + "_result");
439  }
440  } else {
441  cout << "No merged data/result.root found" << endl;
442  }
443 }
444 void CbmTrdOccupancyQa::CreateLayerView() {
445  if (fPlotMergedResults) SwitchToMergedFile();
446  TString origpath = gDirectory->GetPath();
447  printf("\n%s\n", origpath.Data());
448  TString newpath = origpath;
449  newpath.ReplaceAll("eds", "oc_qa");
450  newpath.ReplaceAll(":/", "");
451 
452  TFile* tempFile = new TFile(newpath, "recreate");
453  gDirectory = tempFile->CurrentDirectory();
454  gDirectory->pwd();
455  cout << "CreateLayerView" << endl;
456  Bool_t debug = false;
457  TString title;
458  for (fModuleOccupancyMemoryMapIt = fModuleOccupancyMemoryMap.begin();
459  fModuleOccupancyMemoryMapIt != fModuleOccupancyMemoryMap.end();
460  ++fModuleOccupancyMemoryMapIt) {
461  /*
462  TBox *b = new TBox(
463  fModuleMap[fModuleOccupancyMapIt->first]->ModulePositionX - fModuleMap[fModuleOccupancyMapIt->first]->ModuleSizeX,
464  fModuleMap[fModuleOccupancyMapIt->first]->ModulePositionY - fModuleMap[fModuleOccupancyMapIt->first]->ModuleSizeY,
465  fModuleMap[fModuleOccupancyMapIt->first]->ModulePositionX + fModuleMap[fModuleOccupancyMapIt->first]->ModuleSizeX,
466  fModuleMap[fModuleOccupancyMapIt->first]->ModulePositionY + fModuleMap[fModuleOccupancyMapIt->first]->ModuleSizeY
467  );
468  */
469  /*
470  Float_t occupancy = Float_t(fModuleOccupancyMapIt->second->Integral()) / Float_t(fModuleOccupancyMapIt->second->GetNbinsX() * fModuleOccupancyMapIt->second->GetNbinsY()) * 100.;
471  if (debug)
472  printf("M: %i O:%6.2f ",fModuleOccupancyMapIt->first,occupancy);
473  fModuleOccupancyMemoryMap[fModuleOccupancyMapIt->first]->Fill(occupancy);
474  */
475  Float_t occupancy = fModuleOccupancyMemoryMapIt->second->GetMean(1);
476  Float_t occupancyE = fModuleOccupancyMemoryMapIt->second->GetRMS(1);
477  if (debug) printf(" <O>:%6.2f ", occupancy);
478  title.Form("%.1f#pm%.1f%%", occupancy, occupancyE);
479  if (debug) printf("%s\n", title.Data());
480  TPaveText* text = new TPaveText(
481  fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModulePositionX
482  - fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModuleSizeX,
483  fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModulePositionY
484  - fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModuleSizeY,
485  fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModulePositionX
486  + fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModuleSizeX,
487  fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModulePositionY
488  + fModuleMap[fModuleOccupancyMemoryMapIt->first]->ModuleSizeY);
489  text->SetFillStyle(1001);
490  text->SetLineColor(1);
491 
492  //Double_t min(0), max(12);
493  std::vector<Int_t> fColors;
494  std::vector<Double_t> fZLevel;
495  for (Int_t i = 0; i < TColor::GetNumberOfColors(); i++) {
496  fColors.push_back(TColor::GetColorPalette(i));
497  //fZLevel.push_back(min + TMath::Power(10, TMath::Log10(max) / TColor::GetNumberOfColors() * i));// log scale
498  fZLevel.push_back(
499  fmin + (fmax / TColor::GetNumberOfColors() * i)); // lin scale
500  }
501  // Int_t color(0), j(0);
502  Int_t j(0);
503  while ((occupancy > fZLevel[j]) && (j < (Int_t) fZLevel.size())) {
504  //printf (" %i<%i %i %E <= %E\n",j,(Int_t)fZLevel.size(),fColors[j], rate, fZLevel[j]);
505  j++;
506  }
507  text->SetFillColor(fColors[j]);
508  if (j >= (Int_t) fZLevel.size()) text->SetFillColor(2);
509 
510  if (fColors[j] < 65) text->SetTextColor(kWhite);
511  /*
512  Int_t maxCol = 4;
513  if (occupancy > 0 && occupancy <= 2.5)
514  text->SetFillColor(kViolet
515  //+Int_t((occupancy-0)/6.*maxCol)
516  );
517  if (occupancy > 2.5 && occupancy <= 5){
518  text->SetFillColor(kAzure
519  //+Int_t((occupancy-6)/6.*maxCol)
520  );
521  text->SetTextColor(kWhite);
522  }
523  if (occupancy > 5 && occupancy <= 7.5)
524  text->SetFillColor(kTeal
525  //+Int_t((occupancy-12)/6.*maxCol)
526  );
527  if (occupancy > 7.5 && occupancy <= 10)
528  text->SetFillColor(kSpring
529  //+Int_t((occupancy-18)/6.*maxCol)
530  );
531  if (occupancy > 10 && occupancy <= 12.5)
532  text->SetFillColor(kYellow
533  //+Int_t((occupancy-24)/6.*maxCol)
534  );
535  if (occupancy > 12.5 && occupancy <= 15)
536  text->SetFillColor(kOrange
537  //+Int_t((occupancy-24)/6.*maxCol)
538  );
539  if (occupancy > 15)
540  text->SetFillColor(kRed);
541  */
542  /*
543  if (occupancy > 0 && occupancy <= 5)
544  text->SetFillColor(kViolet
545  //+Int_t((occupancy-0)/6.*maxCol)
546  );
547  if (occupancy > 5 && occupancy <= 10){
548  text->SetFillColor(kAzure
549  //+Int_t((occupancy-6)/6.*maxCol)
550  );
551  text->SetTextColor(kWhite);
552  }
553  if (occupancy > 10 && occupancy <= 15)
554  text->SetFillColor(kTeal
555  //+Int_t((occupancy-12)/6.*maxCol)
556  );
557  if (occupancy > 15 && occupancy <= 20)
558  text->SetFillColor(kSpring
559  //+Int_t((occupancy-18)/6.*maxCol)
560  );
561  if (occupancy > 20 && occupancy <= 25)
562  text->SetFillColor(kYellow
563  //+Int_t((occupancy-24)/6.*maxCol)
564  );
565  if (occupancy > 25 && occupancy <= 30)
566  text->SetFillColor(kOrange
567  //+Int_t((occupancy-24)/6.*maxCol)
568  );
569  if (occupancy > 30)
570  text->SetFillColor(kRed);
571  */
572  fLayerOccupancyMap[fModuleMap[fModuleOccupancyMemoryMapIt->first]->Station
573  * 10
574  + fModuleMap[fModuleOccupancyMemoryMapIt->first]->Layer]
575  ->cd();
576  //b->Draw("same");
577  //text->SetFillStyle(0);
578  text->AddText(title);
579  text->Draw("same");
580  }
581  gDirectory->pwd();
582 
583  TDirectoryFile* TopLayerDir;
584  gDirectory->GetObject("Occupancy", TopLayerDir);
585  if (!TopLayerDir) {
586  std::cout << "unknown Directory Occupancy" << std::endl;
587  gDirectory->mkdir("Occupancy");
588  gDirectory->GetObject("Occupancy", TopLayerDir);
589  }
590  TopLayerDir->cd();
591  gDirectory->pwd();
592  TDirectoryFile* Triggerdir;
593  title.Form("TH%.2E", fTriggerThreshold);
594  gDirectory->GetObject(title, Triggerdir);
595  if (!Triggerdir) {
596  std::cout << "unknown Directory Trigger Dir" << std::endl;
597  gDirectory->mkdir(title);
598  gDirectory->GetObject(title, Triggerdir);
599  }
600  Triggerdir->cd();
601  gDirectory->pwd();
602  for (fLayerOccupancyMapIt = fLayerOccupancyMap.begin();
603  fLayerOccupancyMapIt != fLayerOccupancyMap.end();
604  ++fLayerOccupancyMapIt) {
605  fLayerOccupancyMapIt->second->Write("", TObject::kOverwrite);
606  title.Form("pics/Occupancy_%.2E_%s_%s.pdf",
607  fTriggerThreshold,
608  TString(fLayerOccupancyMapIt->second->GetTitle()).Data(),
609  fGeo.Data());
610  fLayerOccupancyMapIt->second->SaveAs(title);
611  title.Form("pics/Occupancy_%.2E_%s_%s.png",
612  fTriggerThreshold,
613  TString(fLayerOccupancyMapIt->second->GetTitle()).Data(),
614  fGeo.Data());
615  fLayerOccupancyMapIt->second->SaveAs(title);
616  }
617  gDirectory->Cd("..");
618  TCanvas* c = new TCanvas("c", "c", 800, 600);
619  Int_t counter = 0;
620  for (fLayerAverageOccupancyMapIt = fLayerAverageOccupancyMap.begin();
621  fLayerAverageOccupancyMapIt != fLayerAverageOccupancyMap.end();
622  ++fLayerAverageOccupancyMapIt) {
623  fLayerAverageOccupancyMapIt->second->Write("", TObject::kOverwrite);
624 
625  counter++;
626  fLayerAverageOccupancyMapIt->second->SetLineColor(counter);
627  fLayerAverageOccupancyMapIt->second->SetMarkerColor(counter);
628  fLayerAverageOccupancyMapIt->second->SetMarkerStyle(20);
629  c->cd()->SetLogy(1);
630  c->cd()->SetLogx(1);
631  if (fLayerAverageOccupancyMapIt == fLayerAverageOccupancyMap.begin()) {
632  fLayerAverageOccupancyMapIt->second->GetXaxis()->SetRangeUser(1.0e-10,
633  1.0e-3);
634  fLayerAverageOccupancyMapIt->second->GetYaxis()->SetRangeUser(0.1, 100);
635  fLayerAverageOccupancyMapIt->second->DrawCopy();
636  } else
637  fLayerAverageOccupancyMapIt->second->DrawCopy("same");
638  }
639  title.Form("pics/Occupancy_%.2E_%s.pdf", fTriggerThreshold, fGeo.Data());
640  c->SaveAs(title);
641  title.Form("pics/Occupancy_%.2E_%s.png", fTriggerThreshold, fGeo.Data());
642  c->SaveAs(title);
643  gDirectory->Cd("..");
644  c->Close();
645 
646  tempFile->Close();
647  gDirectory->Cd(origpath);
648  gDirectory->pwd();
649 }
650 void CbmTrdOccupancyQa::SetNeighbourReadout(Bool_t neighbourReadout) {
651  fNeigbourReadout = neighbourReadout;
652 }
653 void CbmTrdOccupancyQa::SetTriggerThreshold(Double_t triggerthreshold) {
654  fTriggerThreshold = triggerthreshold;
655 }
656 void CbmTrdOccupancyQa::SaveHistos2File() {
657  TString origpath = gDirectory->GetPath();
658  printf("\n%s\n", origpath.Data());
659  TString newpath = origpath;
660  newpath.ReplaceAll("eds", "oc_qa");
661  newpath.ReplaceAll(":/", "");
662 
663  TFile* tempFile = new TFile(newpath, "Update");
664  gDirectory = tempFile->CurrentDirectory();
665  gDirectory->pwd();
666 
667  cout << "SaveHistos2File" << endl;
668  TString title;
669  /*
670  TFile *outFile = new TFile("data/CbmTrdOccupancyQa.root","UPDATE","output of CbmTrdOccupancyQa");
671  outFile->cd();
672  */
673  gDirectory->pwd();
674  TDirectoryFile* Basedir;
675  gDirectory->GetObject("TrdOccupancy", Basedir);
676  if (!Basedir) {
677  gDirectory->mkdir("TrdOccupancy");
678  std::cout << "unknown Directory Basedir" << std::endl;
679  gDirectory->GetObject("TrdOccupancy", Basedir);
680  }
681  Basedir->cd();
682 
683  fDigiChargeSpectrum->Write("", TObject::kOverwrite);
684  /*
685  for (fLayerOccupancyMapIt = fLayerOccupancyMap.begin();
686  fLayerOccupancyMapIt != fLayerOccupancyMap.end(); ++fLayerOccupancyMapIt) {
687  fLayerOccupancyMapIt->second->Write("", TObject::kOverwrite);
688  }
689  */
690  /*
691  if (!gDirectory->Cd("Module2D"))
692  gDirectory->mkdir("Module2D");
693  gDirectory->Cd("Module2D");
694  for (fModuleOccupancyMapIt = fModuleOccupancyMap.begin();
695  fModuleOccupancyMapIt != fModuleOccupancyMap.end(); ++fModuleOccupancyMapIt) {
696  title.Form("Occupancy: %.2f%%",Float_t(fModuleOccupancyMapIt->second->Integral()) / Float_t(fModuleOccupancyMapIt->second->GetNbinsX() * fModuleOccupancyMapIt->second->GetNbinsY()) * 100.);
697  fModuleOccupancyMapIt->second->SetTitle(title);
698  fModuleOccupancyMapIt->second->Write("", TObject::kOverwrite);
699  }
700  gDirectory->Cd("..");
701  */
702  TDirectoryFile* Moduledir;
703  gDirectory->GetObject("Module1D", Moduledir);
704  if (!Moduledir) {
705  gDirectory->mkdir("Module1D");
706  std::cout << "unknown Directory Module1D" << std::endl;
707  gDirectory->GetObject("Module1D", Moduledir);
708  }
709  Moduledir->cd();
710  gDirectory->pwd();
711 
712  for (fModuleOccupancyMemoryMapIt = fModuleOccupancyMemoryMap.begin();
713  fModuleOccupancyMemoryMapIt != fModuleOccupancyMemoryMap.end();
714  ++fModuleOccupancyMemoryMapIt) {
715  fModuleOccupancyMemoryMapIt->second->Write("", TObject::kOverwrite);
716  }
717  gDirectory->Cd("..");
718  gDirectory->Cd("..");
719  //outFile->Close();
720 
721  tempFile->Close();
722  gDirectory->Cd(origpath);
723  gDirectory->pwd();
724 }
725 ClassImp(CbmTrdOccupancyQa)
CbmTrdAddress::GetModuleAddress
static UInt_t GetModuleAddress(UInt_t address)
Return unique module ID from address.
Definition: CbmTrdAddress.h:117
CbmTrdDigi::GetAddress
Int_t GetAddress() const
Address getter for module in the format defined by CbmTrdDigi (format of CbmTrdAddress can be accesse...
Definition: CbmTrdDigi.h:92
CbmTrdOccupancyQa
Definition: CbmTrdOccupancyQa.h:76
CbmTrdOccupancyQa::FinishTask
virtual void FinishTask()
Definition: CbmTrdOccupancyQa.cxx:396
CbmTrdAddress.h
Helper class to convert unique channel ID back and forth.
CbmTrdOccupancyQa::fGeoHandler
CbmTrdGeoHandler * fGeoHandler
Definition: CbmTrdOccupancyQa.h:108
CbmTrdOccupancyQa::fTriggerThreshold
Double_t fTriggerThreshold
Definition: CbmTrdOccupancyQa.h:123
CbmTrdOccupancyQa::fNeigbourReadout
Bool_t fNeigbourReadout
Definition: CbmTrdOccupancyQa.h:124
CbmTrdParSetGeo.h
CbmTrdParModDigi::GetNofRows
Int_t GetNofRows() const
Definition: CbmTrdParModDigi.cxx:340
CbmDigiManager::Init
InitStatus Init()
Initialisation.
Definition: CbmDigiManager.cxx:71
CbmTrdOccupancyQa::Register
void Register()
Definition: CbmTrdOccupancyQa.cxx:401
CbmTrdOccupancyQa::ReInit
virtual InitStatus ReInit()
Definition: CbmTrdOccupancyQa.cxx:141
CbmTrdAddress::GetSectorId
static UInt_t GetSectorId(UInt_t address)
Return sector ID from address.
Definition: CbmTrdAddress.h:88
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmTrdOccupancyQa::FinishEvent
virtual void FinishEvent()
Definition: CbmTrdOccupancyQa.cxx:385
CbmTrdOccupancyQa::CreateLayerView
void CreateLayerView()
Definition: CbmTrdOccupancyQa.cxx:444
CbmTrdOccupancyQa::fModuleOccupancyMemoryMap
std::map< Int_t, TH1F * > fModuleOccupancyMemoryMap
Definition: CbmTrdOccupancyQa.h:117
CbmTrdOccupancyQa::fModuleMap
std::map< Int_t, OccupancyModule * > fModuleMap
Definition: CbmTrdOccupancyQa.h:113
CbmTrdAddress::GetLayerId
static UInt_t GetLayerId(UInt_t address)
Return layer ID from address.
Definition: CbmTrdAddress.h:69
CbmTrdOccupancyQa::fLayerOccupancyMap
std::map< Int_t, TCanvas * > fLayerOccupancyMap
Definition: CbmTrdOccupancyQa.h:119
CbmDigiManager::GetNofDigis
static Int_t GetNofDigis(ECbmModuleId systemId)
Definition: CbmDigiManager.cxx:62
CbmTrdOccupancyQa::SwitchToMergedFile
void SwitchToMergedFile()
Definition: CbmTrdOccupancyQa.cxx:425
CbmTrdParModDigi::GetSizeX
Double_t GetSizeX() const
Definition: CbmTrdParModDigi.h:109
CbmTrdParModGeo.h
OccupancyModule
struct OccupancyModule OccupancyModule
CbmTrdOccupancyQa::SetNeighbourReadout
void SetNeighbourReadout(Bool_t neighbourReadout)
Definition: CbmTrdOccupancyQa.cxx:650
CbmTrdGeoHandler
Definition: CbmTrdGeoHandler.h:29
CbmDigiManager::Instance
static CbmDigiManager * Instance()
Static instance.
Definition: CbmDigiManager.h:93
CbmCluster::GetNofDigis
Int_t GetNofDigis() const
Number of digis in cluster.
Definition: CbmCluster.h:69
CbmTrdCluster
Data Container for TRD clusters.
Definition: CbmTrdCluster.h:23
CbmTrdDigi.h
CbmTrdOccupancyQa::fDigiPar
CbmTrdParSetDigi * fDigiPar
Definition: CbmTrdOccupancyQa.h:106
CbmTrdOccupancyQa::fDigiChargeSpectrum
TH1I * fDigiChargeSpectrum
Definition: CbmTrdOccupancyQa.h:109
CbmTrdOccupancyQa::fmax
Double_t fmax
Definition: CbmTrdOccupancyQa.h:112
CbmTrdParModGeo
Definition of geometry for one TRD module.
Definition: CbmTrdParModGeo.h:12
CbmTrdGeoHandler.h
Helper class to extract information from the GeoManager.
CbmDigiManager::Get
const Digi * Get(Int_t index) const
Get a digi object.
Definition: CbmDigiManager.h:52
CbmTrdOccupancyQa::SaveHistos2File
void SaveHistos2File()
Definition: CbmTrdOccupancyQa.cxx:656
CbmTrdParModDigi
Definition of chamber gain conversion for one TRD module.
Definition: CbmTrdParModDigi.h:14
CbmTrdGeoHandler::Init
void Init(Bool_t isSimulation=kFALSE)
Definition: CbmTrdGeoHandler.cxx:45
CbmTrdOccupancyQa::fModuleMapIt
std::map< Int_t, OccupancyModule * >::iterator fModuleMapIt
Definition: CbmTrdOccupancyQa.h:114
CbmTrdOccupancyQa::fLayerDummy
TH2I * fLayerDummy
Definition: CbmTrdOccupancyQa.h:110
CbmTrdOccupancyQa.h
CbmTrdOccupancyQa::CbmTrdOccupancyQa
CbmTrdOccupancyQa()
Definition: CbmTrdOccupancyQa.cxx:42
CbmTrdParModDigi.h
CbmTrdOccupancyQa::fLayerAverageOccupancyMap
std::map< Int_t, TProfile * > fLayerAverageOccupancyMap
Definition: CbmTrdOccupancyQa.h:121
CbmTrdOccupancyQa::fGeo
TString fGeo
Definition: CbmTrdOccupancyQa.h:126
CbmTrdParSet::GetModulePar
virtual const CbmTrdParMod * GetModulePar(Int_t detId) const
Definition: CbmTrdParSet.cxx:45
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmTrdParSetGeo
Definition: CbmTrdParSetGeo.h:9
CbmTrdOccupancyQa::fClusters
TClonesArray * fClusters
Definition: CbmTrdOccupancyQa.h:105
CbmTrdParSetDigi.h
ECbmModuleId::kTrd
@ kTrd
Transition Radiation Detector.
CbmTrdOccupancyQa::fModuleOccupancyMemoryMapIt
std::map< Int_t, TH1F * >::iterator fModuleOccupancyMemoryMapIt
Definition: CbmTrdOccupancyQa.h:118
counter
int counter
Definition: CbmMvdSensorDigiToHitTask.cxx:72
CbmDigiManager.h
CbmTrdOccupancyQa::SetTriggerThreshold
void SetTriggerThreshold(Double_t triggerthreshold)
Definition: CbmTrdOccupancyQa.cxx:653
CbmTrdOccupancyQa::fPlotMergedResults
Bool_t fPlotMergedResults
Definition: CbmTrdOccupancyQa.h:125
CbmTrdParModDigi::GetSizeY
Double_t GetSizeY() const
Definition: CbmTrdParModDigi.h:110
CbmTrdOccupancyQa::fModuleOccupancyMapIt
std::map< Int_t, TH2I * >::iterator fModuleOccupancyMapIt
Definition: CbmTrdOccupancyQa.h:116
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:60
CbmTrdOccupancyQa::fmin
Double_t fmin
Definition: CbmTrdOccupancyQa.h:111
CbmTrdDigi
Definition: CbmTrdDigi.h:14
CbmTrdParSetDigi
Definition: CbmTrdParSetDigi.h:15
CbmTrdOccupancyQa::SetParContainers
virtual void SetParContainers()
Definition: CbmTrdOccupancyQa.cxx:132
CbmTrdOccupancyQa::Init
virtual InitStatus Init()
Definition: CbmTrdOccupancyQa.cxx:152
CbmTrdCluster.h
Data Container for TRD clusters.
CbmTrdParModGeo::GetX
virtual Double_t GetX() const
Definition: CbmTrdParModGeo.h:22
CbmTrdParModGeo::GetY
virtual Double_t GetY() const
Definition: CbmTrdParModGeo.h:27
CbmTrdOccupancyQa::CopyEvent2MemoryMap
void CopyEvent2MemoryMap()
Definition: CbmTrdOccupancyQa.cxx:406
CbmTrdParModDigi::GetModuleRow
Int_t GetModuleRow(Int_t &sectorId, Int_t &rowId) const
Definition: CbmTrdParModDigi.cxx:373
CbmTrdOccupancyQa::fLayerAverageOccupancyMapIt
std::map< Int_t, TProfile * >::iterator fLayerAverageOccupancyMapIt
Definition: CbmTrdOccupancyQa.h:122
CbmTrdDigi::GetCharge
Double_t GetCharge() const
Charge getter for SPADIC.
Definition: CbmTrdDigi.cxx:133
CbmTrdOccupancyQa::SetNeighbourTrigger
void SetNeighbourTrigger(Bool_t trigger)
Definition: CbmTrdOccupancyQa.cxx:170
CbmTrdOccupancyQa::Exec
virtual void Exec(Option_t *option)
Definition: CbmTrdOccupancyQa.cxx:175
CbmCluster::GetDigi
Int_t GetDigi(Int_t index) const
Get digi at position index.
Definition: CbmCluster.h:76
CbmTrdOccupancyQa::fGeoPar
CbmTrdParSetGeo * fGeoPar
Definition: CbmTrdOccupancyQa.h:107
CbmTrdOccupancyQa::fLayerOccupancyMapIt
std::map< Int_t, TCanvas * >::iterator fLayerOccupancyMapIt
Definition: CbmTrdOccupancyQa.h:120
CbmTrdOccupancyQa::~CbmTrdOccupancyQa
virtual ~CbmTrdOccupancyQa()
Definition: CbmTrdOccupancyQa.cxx:94
CbmTrdParModDigi::GetNofColumns
Int_t GetNofColumns() const
Definition: CbmTrdParModDigi.cxx:321
CbmTrdAddress::GetRowId
static UInt_t GetRowId(UInt_t address)
Return row ID from address.
Definition: CbmTrdAddress.h:98
CbmTrdAddress::GetColumnId
static UInt_t GetColumnId(UInt_t address)
Return column ID from address.
Definition: CbmTrdAddress.h:107
CbmTrdOccupancyQa::fModuleOccupancyMap
std::map< Int_t, TH2I * > fModuleOccupancyMap
Definition: CbmTrdOccupancyQa.h:115