AliRoot Core  3abf5b4 (3abf5b4)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MUONBusPatchEvolution.C
Go to the documentation of this file.
1 
27 #include "AliCDBEntry.h"
28 #include "AliCDBManager.h"
29 #include "AliDAQ.h"
30 #include "AliMergeableCollection.h"
31 #include "AliMpBusPatch.h"
32 #include "AliMpCDB.h"
33 #include "AliMpDDL.h"
34 #include "AliMpDEIterator.h"
35 #include "AliMpDEManager.h"
36 #include "AliMpDetElement.h"
38 #include "Riostream.h"
39 #include "TCanvas.h"
40 #include "TFile.h"
41 #include "TH1.h"
42 #include "TLegend.h"
43 #include "TLegendEntry.h"
44 #include "TStyle.h"
45 #include "TSystem.h"
46 #include <AliMpDDLStore.h>
47 #include <cassert>
48 #include <limits>
49 #include <map>
50 #include <set>
51 
52 //_________________________________________________________________________________________________
53 void AssertMapping() {
54 
55  AliCDBManager* cdbm = AliCDBManager::Instance();
56 
57  if (!cdbm->IsDefaultStorageSet()) {
58  cdbm->SetDefaultStorage(
59  "local:///cvmfs/alice-ocdb.cern.ch/calibration/data/2015/OCDB");
60  }
61 
62  cdbm->SetRun(0);
63 
65 }
66 
67 //_________________________________________________________________________________________________
68 AliMergeableCollection* BPEVO(Int_t runNumber, const char* output = "mchbepevo.complete.root")
69 {
70  AliCDBManager* cdb = AliCDBManager::Instance();
71 
72 
73  if (!cdb->IsDefaultStorageSet()) {
74  cdb->SetDefaultStorage(
75  "local:///cvmfs/alice-ocdb.cern.ch/calibration/data/2015/OCDB");
76  }
77 
78  AssertMapping();
79 
80  cdb->SetRun(runNumber);
81 
82  AliCDBEntry* e = cdb->Get("MUON/Calib/BPEVO");
83 
84  AliMergeableCollection* hc = static_cast<AliMergeableCollection*>(e->GetObject()->Clone());
85 
86  if (!hc) return 0x0;
87 
88  AliMUONBusPatchEvolution bpe(*hc);
89 
90  bpe.Augment();
91 
92  if (output)
93  {
94  TFile* fout = TFile::Open(gSystem->ExpandPathName(output), "RECREATE");
95  hc->Write("bpevo");
96  delete fout;
97  }
98 
99  return hc;
100 }
101 
102 //_________________________________________________________________________________________________
103 AliMergeableCollection* BPEVO(const char* daoutput, const char* output = "mchbepevo.complete.root")
104 {
106 
107  TFile* f = TFile::Open(gSystem->ExpandPathName(daoutput));
108 
109  AliMergeableCollection* hc = static_cast<AliMergeableCollection*>(f->Get(
110  "bpevo"))->Clone();
111 
112  delete f;
113 
114  AliMUONBusPatchEvolution bpe(*hc);
115 
116  bpe.Augment();
117 
118  if (output)
119  {
120  TFile* fout = TFile::Open(gSystem->ExpandPathName(output), "RECREATE");
121  hc->Write("bpevo");
122  delete fout;
123  }
124 
125  return hc;
126 }
127 
128 //_________________________________________________________________________________________________
129 void MakeConfigFileForBPEVOda(const char* filename = "mchbpevo.conf") {
132 
133  std::ofstream out(gSystem->ExpandPathName(filename));
134 
135  AssertMapping();
136 
137  out << "# Histogram the hit counts in 10 and 60 seconds bins" << std::endl
138  << std::endl;
139  out << "+timeResolutions: 10s" << std::endl;
140  out << "+timeResolutions: 60s" << std::endl << std::endl;
141 
142  out << "# maxDuration set to 5 hours" << std::endl << std::endl;
143  out << "maxDuration: 18000" << std::endl << std::endl;
144 
145  out << "# number of events needed to take a decision" << std::endl
146  << std::endl;
147  out << "nofEventsRequiredForDecision: 1000" << std::endl << std::endl;
148 
149  out << "# occupancy threshold for alarm" << std::endl << std::endl;
150  out << "occupancyThreshold: 0.5" << std::endl << std::endl;
151 
152  TIter next(AliMpDDLStore::Instance()->CreateBusPatchIterator());
153  AliMpBusPatch* bp;
154 
155  std::map<int, int> buspatches;
156 
157  while ((bp = static_cast<AliMpBusPatch*>(next()))) {
158  int npads(0);
159 
160  Int_t detElemId = AliMpDDLStore::Instance()->GetDEfromBus(bp->GetId());
161 
163  detElemId);
164 
165  for (Int_t i = 0; i < bp->GetNofManus(); ++i) {
166  Int_t manuId = bp->GetManuId(i);
167 
168  npads += de->NofChannelsInManu(manuId);
169  }
170 
171  buspatches[bp->GetId()] = npads;
172  }
173 
174  std::map<int, int>::const_iterator it;
175  int total(0);
176 
177  out << "# +bp: buspatchId:nofPadsInThatBusPatch" << std::endl << std::endl;
178 
179  for (it = buspatches.begin(); it != buspatches.end(); ++it) {
180  total += it->second;
181  out << "+bp: " << it->first << ":" << it->second << std::endl;
182  }
183 
184  assert(total == 1064008);
185 }
186 
187 //_________________________________________________________________________________________________
188 void MakeConfigCodeForBPEVOda(const char* buspatchmapname = "buspatches") {
191 
192  AssertMapping();
193 
194  TIter next(AliMpDDLStore::Instance()->CreateBusPatchIterator());
195  AliMpBusPatch* bp;
196 
197  int total(0);
198 
199  while ((bp = static_cast<AliMpBusPatch*>(next()))) {
200  int npads(0);
201 
202  Int_t detElemId = AliMpDDLStore::Instance()->GetDEfromBus(bp->GetId());
203 
205  detElemId);
206 
207  for (Int_t i = 0; i < bp->GetNofManus(); ++i) {
208  Int_t manuId = bp->GetManuId(i);
209 
210  npads += de->NofChannelsInManu(manuId);
211  }
212 
213  std::cout << Form("%s[%d]=%d;", buspatchmapname, bp->GetId(), npads)
214  << std::endl;
215 
216  total += npads;
217  }
218  assert(total == 1064008);
219 }
220 
221 
222 //_________________________________________________________________________________________________
223 void plot(const std::vector<TH1*>& v, Double_t min, Double_t max, Double_t factor)
224 {
225  if ( !v.size() ) return;
226 
227  gStyle->SetOptTitle(0);
228  gStyle->SetOptStat(0);
229 
230  if (!v[0]) return;
231 
232  TString title = v[0]->GetName();
233  title.ReplaceAll("RIGHT","");
234  title.ReplaceAll("LEFT","");
235 
236  TLegend* legend = new TLegend(0.8,0.1,0.99,0.90,title.Data());
237  legend->SetTextSize(0.08);
238  gPad->SetRightMargin(0.22);
239 
240  for ( std::vector<TH1*>::size_type i = 0; i < v.size(); ++i )
241  {
242  if (!v[i]) continue;
243 
244  TString tf = v[i]->GetXaxis()->GetTimeFormat();
245 
246  tf.ReplaceAll(":%s","");
247  v[i]->GetXaxis()->SetTimeFormat(tf.Data());
248 
249  v[i]->Scale(factor);
250  v[i]->SetMaximum(max);
251  v[i]->SetMinimum(min);
252  v[i]->GetYaxis()->SetLabelSize(0.10);
253  v[i]->GetXaxis()->SetLabelSize(0.10);
254  v[i]->GetYaxis()->SetRange(min,max);
255  Double_t integralError;
256  Double_t integral = v[i]->IntegralAndError(1,v[i]->GetNbinsX(),integralError);
257  integralError /= integral;
258  integral /= v[i]->GetNbinsX();
259  integralError *= integral;
260  TLegendEntry* entry = legend->AddEntry(v[i],Form("Mean %5.2g #pm %5.2g %%",integral,integralError));
261  entry->SetTextSize(0.04);
262  if (i)
263  {
264  v[i]->Draw("same");
265  }
266  else
267  {
268  v[i]->Draw();
269  }
270  }
271  legend->Draw();
272 
273 }
274 
275 //_________________________________________________________________________________________________
276 void PlotStationOccupancies(AliMergeableCollection& hc,
277  int timeResolution = 60) {
278 
279  TCanvas* c = new TCanvas("station", "station");
280  c->Divide(1, 5);
281 
282  for (int i = 1; i <= 5; ++i) {
283 
284  c->cd(i);
285 
286  std::vector<TH1*> v;
287 
288  v.push_back(hc.Histo(Form("/STATION/OCC/%ds/STATION%dRIGHT", timeResolution, i)));
289  v.push_back(hc.Histo(Form("/STATION/OCC/%ds/STATION%dLEFT", timeResolution, i)));
290  v.push_back(hc.Histo(Form("/STATION/OCC/%ds/STATION%d", timeResolution, i)));
291 
292  gPad->SetLogy();
293 
294 
295  plot(v,1E-2,10,100.0);
296 
297  }
298 }
299 
300 //_________________________________________________________________________________________________
301 void PlotChamberOccupancies(AliMergeableCollection& hc,
302  int stationId = 3,
303  int timeResolution = 60) {
304 
305  TCanvas* c = new TCanvas("chamber", "chamber");
306 
307  c->Divide(1, 2);
308 
309  for (int i = 1; i <= 2; ++i) {
310 
311  c->cd(i);
312 
313  std::vector<TH1*> v;
314 
315  int chamberId = (stationId-1)*2 + i;
316 
317  v.push_back(hc.Histo(Form("/CHAMBER/OCC/%ds/CHAMBER%dRIGHT", timeResolution, chamberId)));
318  v.push_back(hc.Histo(Form("/CHAMBER/OCC/%ds/CHAMBER%dLEFT", timeResolution, chamberId)));
319  v.push_back(hc.Histo(Form("/CHAMBER/OCC/%ds/CHAMBER%d", timeResolution, chamberId)));
320 
321  gPad->SetLogy();
322 
323  plot(v,1E-2,10,100.0);
324  }
325 }
326 
327 //_________________________________________________________________________________________________
328 void PlotDDLOccupancies(AliMergeableCollection& hc,
329  int timeResolution = 60) {
330 
331  TCanvas* c = new TCanvas("ddl", "ddl");
332 
333  c->Divide(4, 5);
334 
335  for ( int i = 1; i < 20; ++i ) {
336 
337  c->cd(i);
338 
339  std::vector<TH1*> v;
340 
341  v.push_back(hc.Histo(Form("/DDL/OCC/%ds/DDL%04d", timeResolution, 2560+i)));
342 
343  gPad->SetLogy();
344 
345  plot(v,1E-2,10,100.0);
346  }
347 }
348 
349 //_________________________________________________________________________________________________
350 void PlotDEOccupancies(AliMergeableCollection& hc,
351  int chamberId = 5,
352  int timeResolution = 60) {
353 
354  AssertMapping();
355 
356  TCanvas* c = new TCanvas("de", "de");
357 
358  c->Divide(5, 5);
359 
360  AliMpDEIterator it;
361 
362  it.First(chamberId-1);
363 
364  int i(1);
365 
366  while (!it.IsDone()) {
367 
368  c->cd(i);
369 
370  int detElemId = it.CurrentDEId();
371 
372  std::vector<TH1*> v;
373 
374  v.push_back(hc.Histo(Form("/DE/OCC/%ds/DE%04d", timeResolution, detElemId)));
375 
376  it.Next();
377 
378  gPad->SetLogy();
379 
380  plot(v,1E-2,10,100.0);
381 
382  ++i;
383  }
384 }
385 
386 //_________________________________________________________________________________________________
387 void PlotBusPatchOccupancies(AliMergeableCollection& hc,
388  int ddlId = 2575,
389  int timeResolution = 60) {
390 
391  AssertMapping();
392 
393  TCanvas* c = new TCanvas("bp", "bp");
394 
395  AliMpDDL* ddl = AliMpDDLStore::Instance()->GetDDL(ddlId - AliDAQ::DdlIDOffset("MUONTRK"));
396 
397  int nbuspatches = ddl->GetNofBusPatches();
398 
399  c->DivideSquare(nbuspatches);
400 
401  for ( int i = 1; i < nbuspatches; ++i ) {
402 
403  c->cd(i);
404 
405  std::vector<TH1*> v;
406 
407  v.push_back(hc.Histo(Form("/BUSPATCH/OCC/%ds/BP%04d", timeResolution, ddl->GetBusPatchId(i))));
408 
409  gPad->SetLogy();
410 
411  plot(v,1E-1,100,100.0);
412  }
413 }
414 
Int_t GetNofManus() const
TFile * Open(const char *filename, Long64_t &nevents)
The iterator over detection elements.
TStyle * gStyle
void PlotDEOccupancies(AliMergeableCollection &hc, int chamberId=5, int timeResolution=60)
void MakeConfigFileForBPEVOda(const char *filename="mchbpevo.conf")
Bool_t IsDone() const
TFile f("CalibObjects.root")
Int_t GetDEfromBus(Int_t busPatchId) const
Int_t CurrentDEId() const
AliMpDetElement * GetDetElement(Int_t detElemId, Bool_t warn=true) const
Int_t GetManuId(Int_t index) const
The class defines the electronics properties of detection element.
The class defined electronics properties of DDL.
Definition: AliMpDDL.h:20
void PlotChamberOccupancies(AliMergeableCollection &hc, int stationId=3, int timeResolution=60)
static AliMpDDLStore * Instance(Bool_t warn=true)
Int_t GetNofBusPatches() const
Definition: AliMpDDL.cxx:203
AliMpDDL * GetDDL(Int_t ddlId, Bool_t warn=true) const
AliMergeableCollection * BPEVO(Int_t runNumber, const char *output="mchbepevo.complete.root")
TFile * fout
Definition: PlotSys.C:39
Int_t GetId() const
Return the unique Id.
Definition: AliMpBusPatch.h:89
void PlotDDLOccupancies(AliMergeableCollection &hc, int timeResolution=60)
void PlotBusPatchOccupancies(AliMergeableCollection &hc, int ddlId=2575, int timeResolution=60)
void plot(const std::vector< TH1 * > &v, Double_t min, Double_t max, Double_t factor)
Int_t NofChannelsInManu(Int_t manuId) const
The class defines the properties of BusPatch.
Definition: AliMpBusPatch.h:21
void MakeConfigCodeForBPEVOda(const char *buspatchmapname="buspatches")
Int_t GetBusPatchId(Int_t index) const
Definition: AliMpDDL.cxx:211
Utility class to massage the output of the MCHBPEVO DA.
static Bool_t LoadAll(Bool_t warn=false)
Definition: AliMpCDB.cxx:196
void AssertMapping()
void PlotStationOccupancies(AliMergeableCollection &hc, int timeResolution=60)