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