CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
makePlots.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <vector>
3 #include <string>
4 #include <map>
5 #include <sstream>
6 
7 #include "TFile.h"
8 #include "TH1F.h"
9 #include "THStack.h"
10 #include "TCanvas.h"
11 #include "TLegend.h"
12 #include "TPaveText.h"
13 #include "TSystem.h"
14 
15 void setCanvasStyle(TCanvas* c, const bool logScale);
16 void setHistoStyle(TH1* h);
17 void setHistoStackStyle(TH1* h, const unsigned int lineColor);
18 void setLegendStyle(TLegend* l, const unsigned int nColumns);
19 void setPaveTextStyle(TPaveText* t, const bool isHorizontal = true);
20 void fillNormFactorMaps();
21 double findNormFactor(const std::string currentPlotType, const std::string currentPart, const bool stackOption);
23 
24 // Map with <name of tracker part, count of channels in the part>
25 // It is static so it can be read by the functions fillNormFactorMaps() and findNormFactor(...)
26 static std::map<std::string, unsigned int> modulesStackNormFactors;
27 static std::map<std::string, unsigned int> modulesNoStackNormFactors;
28 static std::map<std::string, unsigned int> fibersStackNormFactors;
29 static std::map<std::string, unsigned int> fibersNoStackNormFactors;
30 static std::map<std::string, unsigned int> APVsStackNormFactors;
31 static std::map<std::string, unsigned int> APVsNoStackNormFactors;
32 static std::map<std::string, unsigned int> stripsStackNormFactors;
33 static std::map<std::string, unsigned int> stripsNoStackNormFactors;
34 
35 int main(int argc, char* argv[]) {
36  if (argc == 3) {
37  char* inputFileName = argv[1];
38  char* outputFileName = argv[2];
39 
40  std::cout << "ready to make plots from " << inputFileName << " to " << outputFileName << std::endl;
41 
42  int returncode = 0;
44 
45  return returncode;
46 
47  } else {
48  std::cout << "Too few arguments: " << argc << std::endl;
49  return -1;
50  }
51 
52  return -9;
53 }
54 
56  //
57 
58  // Open input and output file
59  TFile* inputFile = new TFile(inputFileName.c_str(), "READ");
60  TFile* outputFile = new TFile(outputFileName.c_str(), "RECREATE");
61 
62  std::ostringstream oss;
64  std::vector<std::string> plotType;
65  plotType.push_back("BadModules");
66  plotType.push_back("BadFibers");
67  plotType.push_back("BadAPVs");
68  plotType.push_back("BadStrips");
69  plotType.push_back("BadStripsFromAPVs");
70  plotType.push_back("AllBadStrips");
71  std::vector<std::string> subDetName;
72  subDetName.push_back("TIB");
73  subDetName.push_back("TID+");
74  subDetName.push_back("TID-");
75  subDetName.push_back("TOB");
76  subDetName.push_back("TEC+");
77  subDetName.push_back("TEC-");
78 
79  // Standard plot options for THStack and for standalone histograms
80  const bool stackHistograms = true;
81 
82  std::string plotStackOptions;
83  if (stackHistograms)
84  plotStackOptions = "";
85  else
86  plotStackOptions = "nostack p";
87  const std::string plotHistoOptions("p");
88 
89  // Defer the filling of the normFactor maps to this function
90  // Conceptually trivial but lengthy
92 
93  // // Finds number of channels from above map
94  // std::string completePartName = subDetName;
95  // if(partName.compare("") != 0)
96  // completePartName += " " + partName;
97  // if(partNumber != 0)
98  // {
99  // oss.str("");
100  // oss << partNumber;
101  // completePartName += " " + oss.str();
102  // }
103  //
104  // // Total number of channels in currently processed map
105  // const unsigned int nModulesInPart = allModulesTK[completePartName.c_str()];
106  // const unsigned int nFibersInPart = allFibersTK[completePartName.c_str()];
107  // const unsigned int nAPVsInPart = allAPVsTK[completePartName.c_str()];
108  // const unsigned int nStripsInPart = allStripsTK[completePartName.c_str()];
109 
110  TH1F* hTracker;
111  TH1F* hTIB;
112  TH1F* hTID;
113  TH1F* hTOB;
114  TH1F* hTEC;
115  TH1F* histo;
116  TH1F* histo2;
117  THStack* histoStack;
118  TCanvas* c1;
119  TCanvas* c2;
120  TLegend* legend;
121  TLegend* legend2;
122  TPaveText* textX;
123  TPaveText* textY;
124  std::string entryLabel;
125  std::vector<TH1F*> hLayers;
126  // unsigned int normFactor;
127 
128  for (std::vector<std::string>::iterator itPlot = plotType.begin(); itPlot != plotType.end(); itPlot++) {
129  // Put together the Tracker histograms with the TIB, TID, TOB and TEC ones
130 
131  histoName = "h" + *itPlot + "Tracker";
132  hTracker = (TH1F*)inputFile->Get(histoName.c_str());
133  if (hTracker) {
134  hTracker->Scale(1 / findNormFactor(*itPlot, "Tracker", stackHistograms));
135  } else {
136  std::cout << histoName << " not found" << std::endl;
137  }
138 
139  histoStack = new THStack(histoName.c_str(), histoName.c_str());
140 
141  histoName = "h" + *itPlot + "TIB";
142  hTIB = (TH1F*)inputFile->Get(histoName.c_str());
143  if (hTIB) {
144  hTIB->Scale(1 / findNormFactor(*itPlot, "TIB", stackHistograms));
145  } else {
146  std::cout << histoName << " not found" << std::endl;
147  }
148 
149  histoName = "h" + *itPlot + "TID";
150  hTID = (TH1F*)inputFile->Get(histoName.c_str());
151  if (hTID) {
152  hTID->Scale(1 / findNormFactor(*itPlot, "TID", stackHistograms));
153  } else {
154  std::cout << histoName << " not found" << std::endl;
155  }
156 
157  histoName = "h" + *itPlot + "TOB";
158  hTOB = (TH1F*)inputFile->Get(histoName.c_str());
159  if (hTOB) {
160  hTOB->Scale(1 / findNormFactor(*itPlot, "TOB", stackHistograms));
161  } else {
162  std::cout << histoName << " not found" << std::endl;
163  }
164 
165  histoName = "h" + *itPlot + "TEC";
166  hTEC = (TH1F*)inputFile->Get(histoName.c_str());
167  if (hTEC) {
168  hTEC->Scale(1 / findNormFactor(*itPlot, "TEC", stackHistograms));
169  } else {
170  std::cout << histoName << " not found" << std::endl;
171  }
172 
173  c1 = new TCanvas(("c" + *itPlot + "Tracker").c_str(), "", 1200, 600);
174  setCanvasStyle(c1, false);
175  // hTracker->Draw();
176  if (hTracker)
177  setHistoStackStyle(hTracker, 1);
178  // hTIB->Draw("same");
179  if (hTIB) {
180  setHistoStackStyle(hTIB, 2);
181  histoStack->Add(hTIB);
182  }
183  // hTID->Draw("same");
184  if (hTID) {
185  setHistoStackStyle(hTID, 3);
186  histoStack->Add(hTID);
187  }
188  // hTOB->Draw("same");
189  if (hTOB) {
190  setHistoStackStyle(hTOB, 4);
191  histoStack->Add(hTOB);
192  }
193  // hTEC->Draw("same");
194  if (hTEC) {
195  setHistoStackStyle(hTEC, 6);
196  histoStack->Add(hTEC);
197  }
198  // Bug in ROOT? If we plot a stack with the "stack" option and the Y axis is in log scale,
199  // but there are no entries in any of the histograms of the stack, then ROOT crashes
200  // Workaround: at this stage, check that at least one histogram has >0 entries,
201  // otherwise, switch back to linear Y scale
202  // Curiously, there is no crash if the "nostack" option is chosen...
203  double histoStackMaximum = histoStack->GetMaximum();
204  if (histoStackMaximum == 0) {
205  c1->SetLogy(0);
206  }
207  histoStack->Draw(plotStackOptions.c_str());
208  if (histoStack->GetYaxis()->GetXmax() > 0.9)
209  histoStack->GetYaxis()->SetRangeUser(0., 0.1);
210  legend = new TLegend(0.4, 0.9, 0.9, 1);
211  legend->AddEntry(hTIB, "TIB");
212  legend->AddEntry(hTID, "TID");
213  legend->AddEntry(hTOB, "TOB");
214  legend->AddEntry(hTEC, "TEC");
216  legend->Draw();
217  textX = new TPaveText();
218  textY = new TPaveText();
219  setPaveTextStyle(textX);
220  setPaveTextStyle(textY, false);
221  textX->Draw();
222  textY->Draw();
223  gSystem->ProcessEvents();
224  c1->Update();
225  outputFile->cd();
226  c1->Write();
227  c1->SaveAs((*itPlot + "Tracker.png").c_str());
228 
229  delete histoStack;
230  delete textX;
231  delete textY;
232  delete legend;
233  legend = nullptr;
234  delete c1;
235 
236  // Put together the histograms for the different layers of the detectors
237  for (std::vector<std::string>::iterator itSub = subDetName.begin(); itSub != subDetName.end(); itSub++) {
238  unsigned int nLayers = 0;
239  std::string layerName;
240  // std::cout << "itSub = " << (*itSub).c_str() << std::endl;
241  if ((*itSub) == "TIB") {
242  nLayers = 4;
243  layerName = "Layer";
244  legend = new TLegend(0.4, 0.9, 0.9, 1);
246  } else if ((*itSub) == "TID+" || (*itSub) == "TID-") {
247  nLayers = 3;
248  layerName = "Disk";
249  legend = new TLegend(0.35, 0.9, 0.9, 1);
251  } else if ((*itSub) == "TOB") {
252  nLayers = 6;
253  layerName = "Layer";
254  legend = new TLegend(0.35, 0.9, 0.9, 1);
256  } else if ((*itSub) == "TEC+" || (*itSub) == "TEC-") {
257  nLayers = 9;
258  layerName = "Disk";
259  legend = new TLegend(0.35, 0.9, 1, 1);
261  }
262 
263  if (legend == nullptr) {
264  std::cerr << "Unknown itSub type " << *itSub << std::endl;
265  assert(false);
266  }
267 
268  c1 = new TCanvas(("c" + *itPlot + *itSub).c_str(), "", 1200, 600);
269  setCanvasStyle(c1, false);
270  // if((*itSub).compare("TEC+")==0 || (*itSub).compare("TEC-")==0)
271  // {
272  // histoName = "h" + *itPlot + "TEC";
273  // }
274  // else
275  // {
276  histoName = "h" + *itPlot + *itSub;
277  // }
278  // hSubDet = (TH1F*)inputFile->Get(histoName.c_str());
279  // setHistoStackStyle(hSubDet,1);
280  // hSubDet->Draw();
281  histoStack = new THStack(histoName.c_str(), histoName.c_str());
282 
283  for (unsigned int iLayer = 1; iLayer <= nLayers; iLayer++) {
284  oss.str("");
285  oss << iLayer;
286  // // TIB and TOB have no plus/minus side division
287  // // While TEC has it but I plot them separately
288  // // TID has plus/minus side division but I plot everything in a single plot
289  // if((*itSub).compare("TID")==0)
290  // {
291  // histoName = "h" + *itPlot + *itSub + "+" + layerName + oss.str();
292  // // std::cout << "histoName = " << histoName.c_str() << std::endl;
293  // histo = (TH1F*)inputFile->Get(histoName.c_str());
294  //
295  // // First: plot histogram separately
296  // setHistoStyle(histo);
297  // c2 = new TCanvas(("c" + *itPlot + *itSub + "+" + layerName + oss.str()).c_str(), "", 1200, 600);
298  // setCanvasStyle(c2, true);
299  // histo->Draw(plotHistoOptions.c_str());
300  // legend2 = new TLegend(0.6,0.92,0.9,0.97);
301  // legend2->AddEntry(histo,(*itSub + "+" + layerName + oss.str()).c_str());
302  // setLegendStyle(legend2, 1);
303  // legend2->Draw();
304  // gSystem->ProcessEvents();
305  // c2->Update();
306  // outputFile->cd();
307  // c2->Write();
308  // c2->SaveAs((*itPlot + *itSub + "+" + layerName + oss.str() + ".png").c_str());
309  // delete legend2;
310  // delete c2;
311  //
312  // // Second: add histogram to THStack
313  // hLayers.push_back(histo);
314  // setHistoStackStyle(hLayers.back(), iLayer);
315  // histoStack->Add(hLayers.back());
316  // entryLabel = *itSub + "+ " + layerName + " " + oss.str();
317  // legend->AddEntry(hLayers.back(), entryLabel.c_str());
318  //
319  //
320  // histoName = "h" + *itPlot + *itSub + "-" + layerName + oss.str();
321  // // std::cout << "histoName = " << histoName.c_str() << std::endl;
322  // histo = (TH1F*)inputFile->Get(histoName.c_str());
323  //
324  // // First: plot histogram separately
325  // setHistoStyle(histo);
326  // c2 = new TCanvas(("c" + *itPlot + *itSub + "-" + layerName + oss.str()).c_str(), "", 1200, 600);
327  // setCanvasStyle(c2, true);
328  // histo->Draw(plotHistoOptions.c_str());
329  // legend2 = new TLegend(0.6,0.92,0.9,0.97);
330  // legend2->AddEntry(histo,(*itSub + "-" + layerName + oss.str()).c_str());
331  // setLegendStyle(legend2, 1);
332  // legend2->Draw();
333  // gSystem->ProcessEvents();
334  // c2->Update();
335  // outputFile->cd();
336  // c2->Write();
337  // c2->SaveAs((*itPlot + *itSub + "-" + layerName + oss.str() + ".png").c_str());
338  // delete legend2;
339  // delete c2;
340  //
341  // // Second: add histogram to THStack
342  // hLayers.push_back(histo);
343  // setHistoStackStyle(hLayers.back(), iLayer+nLayers);
344  // histoStack->Add(hLayers.back());
345  // entryLabel = *itSub + "- " + layerName + " " + oss.str();
346  // legend->AddEntry(hLayers.back(), entryLabel.c_str());
347  // // hLayers.back()->Draw("same");
348  // }
349  // else
350  // {
351  histoName = "h" + *itPlot + *itSub + layerName + oss.str();
352  // std::cout << "histoName = " << histoName.c_str() << std::endl;
353  histo = (TH1F*)inputFile->Get(histoName.c_str());
354  if (histo) {
355  histo2 = new TH1F(*histo);
356  histo->Scale(1 / findNormFactor(*itPlot, *itSub + " " + layerName + " " + oss.str(), false));
357  histo2->Scale(1 / findNormFactor(*itPlot, *itSub + " " + layerName + " " + oss.str(), stackHistograms));
358  // First: plot histogram separately
360  c2 = new TCanvas(("c" + *itPlot + *itSub + layerName + oss.str()).c_str(), "", 1200, 600);
361  setCanvasStyle(c2, true);
362  histo->Draw(plotHistoOptions.c_str());
363  double histoMaximum = histo->GetMaximum();
364  // Otherwise it does not draw the pad
365  if (histoMaximum == 0) {
366  c2->SetLogy(0);
367  }
368  legend2 = new TLegend(0.6, 0.92, 0.9, 0.97);
369  legend2->AddEntry(histo, (*itSub + layerName + oss.str()).c_str());
370  setLegendStyle(legend2, 1);
371  legend2->Draw();
372  textX = new TPaveText();
373  textY = new TPaveText();
374  setPaveTextStyle(textX);
375  setPaveTextStyle(textY, false);
376  textX->Draw();
377  textY->Draw();
378  gSystem->ProcessEvents();
379  c2->Update();
380  outputFile->cd();
381  c2->Write();
382  c2->SaveAs((*itPlot + *itSub + layerName + oss.str() + ".png").c_str());
383  delete textX;
384  delete textY;
385  delete legend2;
386  delete c2;
387 
388  // Second: add histogram to THStack
389  setHistoStyle(histo2);
390  hLayers.push_back(histo2);
391  setHistoStackStyle(hLayers.back(), iLayer);
392  histoStack->Add(hLayers.back());
393  entryLabel = *itSub + " " + layerName + " " + oss.str();
394  legend->AddEntry(hLayers.back(), entryLabel.c_str());
395  // hLayers.back()->Draw("same");
396  // }
397  // delete histo2;
398  } else {
399  std::cout << histoName << " not found" << std::endl;
400  }
401  }
402  histoStack->Draw(plotStackOptions.c_str());
403 
404  // Bug in ROOT? If we plot a stack with the "stack" option and the Y axis is in log scale,
405  // but there are no entries in any of the histograms of the stack, then ROOT crashes
406  // Workaround: at this stage, check that at least one histogram has >0 entries,
407  // otherwise, switch back to linear Y scale
408  // Curiously, there is no crash if the "nostack" option is chosen...
409  double histoStackMaximum = histoStack->GetMaximum();
410  if (histoStackMaximum == 0) {
411  c1->SetLogy(0);
412  }
413  if (histoStackMaximum > 0.01)
414  histoStack->SetMaximum(0.01);
415  textX = new TPaveText();
416  textY = new TPaveText();
417  setPaveTextStyle(textX);
418  setPaveTextStyle(textY, false);
419  textX->Draw();
420  textY->Draw();
421  legend->Draw();
422  gSystem->ProcessEvents();
423  c1->Update();
424  outputFile->cd();
425  c1->Write();
426  c1->SaveAs((*itPlot + *itSub + ".png").c_str());
427  delete histoStack;
428  delete textX;
429  delete textY;
430  delete legend;
431  legend = nullptr;
432  delete c1;
433  }
434  }
435 
436  inputFile->Close();
437  outputFile->Close();
438 }
439 
440 void setCanvasStyle(TCanvas* c, const bool logScale) {
441  c->SetFillColor(0);
442  c->SetFrameBorderMode(0);
443  if (logScale)
444  c->SetLogy(1);
445  else
446  c->SetLogy(0);
447  c->SetCanvasSize(1200, 600);
448  c->SetWindowSize(1200, 600);
449 }
450 
451 void setHistoStyle(TH1* h) {
452  h->SetLineStyle(0);
453  h->SetMarkerStyle(3);
454  h->SetMarkerSize(1);
455  h->SetMarkerColor(1);
456  h->SetStats(kFALSE);
457  // h->GetYaxis()->SetTitle("Fraction of total");
458  // h->GetXaxis()->SetTitle("IOV");
459  // h->GetXaxis()->SetTitleOffset(-0.3);
460  // Avoid having too many bins with labels
461  if (h->GetNbinsX() > 25)
462  for (int i = 1; i < h->GetNbinsX() - 1; i++)
463  if ((i % (h->GetNbinsX() / 25 + 1)))
464  h->GetXaxis()->SetBinLabel(i + 1, "");
465 }
466 
467 void setHistoStackStyle(TH1* h, const unsigned int lineColor) {
468  h->SetLineStyle(0);
469  // h->SetDrawOption("e1p");
470  // Best marker types are 20-23 - use them with different colors
471  h->SetMarkerStyle(lineColor % 4 + 20);
472  h->SetMarkerSize(1);
473  h->SetMarkerColor(lineColor);
474  h->SetLineColor(lineColor);
475  h->SetFillColor(lineColor);
476  h->SetLineWidth(2);
477  h->SetStats(kFALSE);
478  h->GetYaxis()->SetTitle("Fraction of total");
479  // h->GetXaxis()->SetTitle("IOV");
480  // h->GetXaxis()->SetTitleOffset(1.2);
481  // Avoid having too many bins with labels
482  if (h->GetNbinsX() > 25)
483  for (int i = 1; i < h->GetNbinsX() - 1; i++)
484  if (i % (h->GetNbinsX() / 25 + 1))
485  h->GetXaxis()->SetBinLabel(i + 1, "");
486 }
487 
488 void setLegendStyle(TLegend* l, const unsigned int nColumns) {
489  l->SetNColumns(nColumns);
490  l->SetFillColor(0);
491 }
492 
493 void setPaveTextStyle(TPaveText* t, const bool isHorizontal) {
494  t->SetLineStyle(0);
495  t->SetFillColor(0);
496  t->SetFillStyle(0);
497  t->SetBorderSize(0);
498  if (isHorizontal) {
499  t->SetX1NDC(0.905);
500  t->SetX2NDC(0.975);
501  t->SetY1NDC(0.062);
502  t->SetY2NDC(0.095);
503  t->AddText("IOV");
504  } else {
505  t->SetX1NDC(0.03);
506  t->SetX2NDC(0.05);
507  t->SetY1NDC(0.33);
508  t->SetY2NDC(0.68);
509  TText* t1 = t->AddText("Fraction of total");
510  t1->SetTextAngle(90.);
511  }
512 }
513 
515  // Number of modules, fibers, APVs, strips for each tracker part
516  std::vector<std::string> tkParts;
517  tkParts.push_back("Tracker");
518  tkParts.push_back("TIB");
519  tkParts.push_back("TID");
520  tkParts.push_back("TOB");
521  tkParts.push_back("TEC");
522  tkParts.push_back("TIB Layer 1");
523  tkParts.push_back("TIB Layer 2");
524  tkParts.push_back("TIB Layer 3");
525  tkParts.push_back("TIB Layer 4");
526  tkParts.push_back("TID+ Disk 1");
527  tkParts.push_back("TID+ Disk 2");
528  tkParts.push_back("TID+ Disk 3");
529  tkParts.push_back("TID- Disk 1");
530  tkParts.push_back("TID- Disk 2");
531  tkParts.push_back("TID- Disk 3");
532  tkParts.push_back("TOB Layer 1");
533  tkParts.push_back("TOB Layer 2");
534  tkParts.push_back("TOB Layer 3");
535  tkParts.push_back("TOB Layer 4");
536  tkParts.push_back("TOB Layer 5");
537  tkParts.push_back("TOB Layer 6");
538  tkParts.push_back("TEC+ Disk 1");
539  tkParts.push_back("TEC+ Disk 2");
540  tkParts.push_back("TEC+ Disk 3");
541  tkParts.push_back("TEC+ Disk 4");
542  tkParts.push_back("TEC+ Disk 5");
543  tkParts.push_back("TEC+ Disk 6");
544  tkParts.push_back("TEC+ Disk 7");
545  tkParts.push_back("TEC+ Disk 8");
546  tkParts.push_back("TEC+ Disk 9");
547  tkParts.push_back("TEC- Disk 1");
548  tkParts.push_back("TEC- Disk 2");
549  tkParts.push_back("TEC- Disk 3");
550  tkParts.push_back("TEC- Disk 4");
551  tkParts.push_back("TEC- Disk 5");
552  tkParts.push_back("TEC- Disk 6");
553  tkParts.push_back("TEC- Disk 7");
554  tkParts.push_back("TEC- Disk 8");
555  tkParts.push_back("TEC- Disk 9");
556 
557  std::vector<unsigned int> nModulesStack;
558  nModulesStack.push_back(15148);
559  nModulesStack.push_back(15148);
560  nModulesStack.push_back(15148);
561  nModulesStack.push_back(15148);
562  nModulesStack.push_back(15148);
563  nModulesStack.push_back(2724);
564  nModulesStack.push_back(2724);
565  nModulesStack.push_back(2724);
566  nModulesStack.push_back(2724);
567  nModulesStack.push_back(408);
568  nModulesStack.push_back(408);
569  nModulesStack.push_back(408);
570  nModulesStack.push_back(408);
571  nModulesStack.push_back(408);
572  nModulesStack.push_back(408);
573  nModulesStack.push_back(5208);
574  nModulesStack.push_back(5208);
575  nModulesStack.push_back(5208);
576  nModulesStack.push_back(5208);
577  nModulesStack.push_back(5208);
578  nModulesStack.push_back(5208);
579  nModulesStack.push_back(3200);
580  nModulesStack.push_back(3200);
581  nModulesStack.push_back(3200);
582  nModulesStack.push_back(3200);
583  nModulesStack.push_back(3200);
584  nModulesStack.push_back(3200);
585  nModulesStack.push_back(3200);
586  nModulesStack.push_back(3200);
587  nModulesStack.push_back(3200);
588  nModulesStack.push_back(3200);
589  nModulesStack.push_back(3200);
590  nModulesStack.push_back(3200);
591  nModulesStack.push_back(3200);
592  nModulesStack.push_back(3200);
593  nModulesStack.push_back(3200);
594  nModulesStack.push_back(3200);
595  nModulesStack.push_back(3200);
596  nModulesStack.push_back(3200);
597 
598  std::vector<unsigned int> nModulesNoStack;
599  nModulesNoStack.push_back(15148);
600  nModulesNoStack.push_back(2724);
601  nModulesNoStack.push_back(816);
602  nModulesNoStack.push_back(5208);
603  nModulesNoStack.push_back(6400);
604  nModulesNoStack.push_back(672);
605  nModulesNoStack.push_back(864);
606  nModulesNoStack.push_back(540);
607  nModulesNoStack.push_back(648);
608  nModulesNoStack.push_back(136);
609  nModulesNoStack.push_back(136);
610  nModulesNoStack.push_back(136);
611  nModulesNoStack.push_back(136);
612  nModulesNoStack.push_back(136);
613  nModulesNoStack.push_back(136);
614  nModulesNoStack.push_back(1008);
615  nModulesNoStack.push_back(1152);
616  nModulesNoStack.push_back(648);
617  nModulesNoStack.push_back(720);
618  nModulesNoStack.push_back(792);
619  nModulesNoStack.push_back(888);
620  nModulesNoStack.push_back(408);
621  nModulesNoStack.push_back(408);
622  nModulesNoStack.push_back(408);
623  nModulesNoStack.push_back(360);
624  nModulesNoStack.push_back(360);
625  nModulesNoStack.push_back(360);
626  nModulesNoStack.push_back(312);
627  nModulesNoStack.push_back(312);
628  nModulesNoStack.push_back(272);
629  nModulesNoStack.push_back(408);
630  nModulesNoStack.push_back(408);
631  nModulesNoStack.push_back(408);
632  nModulesNoStack.push_back(360);
633  nModulesNoStack.push_back(360);
634  nModulesNoStack.push_back(360);
635  nModulesNoStack.push_back(312);
636  nModulesNoStack.push_back(312);
637  nModulesNoStack.push_back(272);
638 
639  //
640  std::vector<unsigned int> nFibersStack;
641  nFibersStack.push_back(36392);
642  nFibersStack.push_back(36392);
643  nFibersStack.push_back(36392);
644  nFibersStack.push_back(36392);
645  nFibersStack.push_back(36392);
646  nFibersStack.push_back(6984);
647  nFibersStack.push_back(6984);
648  nFibersStack.push_back(6984);
649  nFibersStack.push_back(6984);
650  nFibersStack.push_back(1104);
651  nFibersStack.push_back(1104);
652  nFibersStack.push_back(1104);
653  nFibersStack.push_back(1104);
654  nFibersStack.push_back(1104);
655  nFibersStack.push_back(1104);
656  nFibersStack.push_back(12096);
657  nFibersStack.push_back(12096);
658  nFibersStack.push_back(12096);
659  nFibersStack.push_back(12096);
660  nFibersStack.push_back(12096);
661  nFibersStack.push_back(12096);
662  nFibersStack.push_back(7552);
663  nFibersStack.push_back(7552);
664  nFibersStack.push_back(7552);
665  nFibersStack.push_back(7552);
666  nFibersStack.push_back(7552);
667  nFibersStack.push_back(7552);
668  nFibersStack.push_back(7552);
669  nFibersStack.push_back(7552);
670  nFibersStack.push_back(7552);
671  nFibersStack.push_back(7552);
672  nFibersStack.push_back(7552);
673  nFibersStack.push_back(7552);
674  nFibersStack.push_back(7552);
675  nFibersStack.push_back(7552);
676  nFibersStack.push_back(7552);
677  nFibersStack.push_back(7552);
678  nFibersStack.push_back(7552);
679  nFibersStack.push_back(7552);
680 
681  std::vector<unsigned int> nFibersNoStack;
682  nFibersNoStack.push_back(36392);
683  nFibersNoStack.push_back(6984);
684  nFibersNoStack.push_back(2208);
685  nFibersNoStack.push_back(12096);
686  nFibersNoStack.push_back(15104);
687  nFibersNoStack.push_back(2016);
688  nFibersNoStack.push_back(2592);
689  nFibersNoStack.push_back(1080);
690  nFibersNoStack.push_back(1296);
691  nFibersNoStack.push_back(368);
692  nFibersNoStack.push_back(368);
693  nFibersNoStack.push_back(368);
694  nFibersNoStack.push_back(368);
695  nFibersNoStack.push_back(368);
696  nFibersNoStack.push_back(368);
697  nFibersNoStack.push_back(2016);
698  nFibersNoStack.push_back(2304);
699  nFibersNoStack.push_back(1296);
700  nFibersNoStack.push_back(1440);
701  nFibersNoStack.push_back(2376);
702  nFibersNoStack.push_back(2664);
703  nFibersNoStack.push_back(992);
704  nFibersNoStack.push_back(992);
705  nFibersNoStack.push_back(992);
706  nFibersNoStack.push_back(848);
707  nFibersNoStack.push_back(848);
708  nFibersNoStack.push_back(848);
709  nFibersNoStack.push_back(704);
710  nFibersNoStack.push_back(704);
711  nFibersNoStack.push_back(624);
712  nFibersNoStack.push_back(992);
713  nFibersNoStack.push_back(992);
714  nFibersNoStack.push_back(992);
715  nFibersNoStack.push_back(848);
716  nFibersNoStack.push_back(848);
717  nFibersNoStack.push_back(848);
718  nFibersNoStack.push_back(704);
719  nFibersNoStack.push_back(704);
720  nFibersNoStack.push_back(624);
721 
722  //
723  std::vector<unsigned int> nAPVsStack;
724  nAPVsStack.push_back(72784);
725  nAPVsStack.push_back(72784);
726  nAPVsStack.push_back(72784);
727  nAPVsStack.push_back(72784);
728  nAPVsStack.push_back(72784);
729  nAPVsStack.push_back(13968);
730  nAPVsStack.push_back(13968);
731  nAPVsStack.push_back(13968);
732  nAPVsStack.push_back(13968);
733  nAPVsStack.push_back(2208);
734  nAPVsStack.push_back(2208);
735  nAPVsStack.push_back(2208);
736  nAPVsStack.push_back(2208);
737  nAPVsStack.push_back(2208);
738  nAPVsStack.push_back(2208);
739  nAPVsStack.push_back(24192);
740  nAPVsStack.push_back(24192);
741  nAPVsStack.push_back(24192);
742  nAPVsStack.push_back(24192);
743  nAPVsStack.push_back(24192);
744  nAPVsStack.push_back(24192);
745  nAPVsStack.push_back(15104);
746  nAPVsStack.push_back(15104);
747  nAPVsStack.push_back(15104);
748  nAPVsStack.push_back(15104);
749  nAPVsStack.push_back(15104);
750  nAPVsStack.push_back(15104);
751  nAPVsStack.push_back(15104);
752  nAPVsStack.push_back(15104);
753  nAPVsStack.push_back(15104);
754  nAPVsStack.push_back(15104);
755  nAPVsStack.push_back(15104);
756  nAPVsStack.push_back(15104);
757  nAPVsStack.push_back(15104);
758  nAPVsStack.push_back(15104);
759  nAPVsStack.push_back(15104);
760  nAPVsStack.push_back(15104);
761  nAPVsStack.push_back(15104);
762  nAPVsStack.push_back(15104);
763 
764  std::vector<unsigned int> nAPVsNoStack;
765  nAPVsNoStack.push_back(72784);
766  nAPVsNoStack.push_back(13968);
767  nAPVsNoStack.push_back(4416);
768  nAPVsNoStack.push_back(24192);
769  nAPVsNoStack.push_back(30208);
770  nAPVsNoStack.push_back(4032);
771  nAPVsNoStack.push_back(5184);
772  nAPVsNoStack.push_back(2160);
773  nAPVsNoStack.push_back(2592);
774  nAPVsNoStack.push_back(736);
775  nAPVsNoStack.push_back(736);
776  nAPVsNoStack.push_back(736);
777  nAPVsNoStack.push_back(736);
778  nAPVsNoStack.push_back(736);
779  nAPVsNoStack.push_back(736);
780  nAPVsNoStack.push_back(4032);
781  nAPVsNoStack.push_back(4608);
782  nAPVsNoStack.push_back(2592);
783  nAPVsNoStack.push_back(2880);
784  nAPVsNoStack.push_back(4752);
785  nAPVsNoStack.push_back(5328);
786  nAPVsNoStack.push_back(1984);
787  nAPVsNoStack.push_back(1984);
788  nAPVsNoStack.push_back(1984);
789  nAPVsNoStack.push_back(1696);
790  nAPVsNoStack.push_back(1696);
791  nAPVsNoStack.push_back(1696);
792  nAPVsNoStack.push_back(1408);
793  nAPVsNoStack.push_back(1408);
794  nAPVsNoStack.push_back(1248);
795  nAPVsNoStack.push_back(1984);
796  nAPVsNoStack.push_back(1984);
797  nAPVsNoStack.push_back(1984);
798  nAPVsNoStack.push_back(1696);
799  nAPVsNoStack.push_back(1696);
800  nAPVsNoStack.push_back(1696);
801  nAPVsNoStack.push_back(1408);
802  nAPVsNoStack.push_back(1408);
803  nAPVsNoStack.push_back(1248);
804 
805  //
806  std::vector<unsigned int> nStripsStack;
807  nStripsStack.push_back(9316352);
808  nStripsStack.push_back(9316352);
809  nStripsStack.push_back(9316352);
810  nStripsStack.push_back(9316352);
811  nStripsStack.push_back(9316352);
812  nStripsStack.push_back(1787904);
813  nStripsStack.push_back(1787904);
814  nStripsStack.push_back(1787904);
815  nStripsStack.push_back(1787904);
816  nStripsStack.push_back(282624);
817  nStripsStack.push_back(282624);
818  nStripsStack.push_back(282624);
819  nStripsStack.push_back(282624);
820  nStripsStack.push_back(282624);
821  nStripsStack.push_back(282624);
822  nStripsStack.push_back(3096576);
823  nStripsStack.push_back(3096576);
824  nStripsStack.push_back(3096576);
825  nStripsStack.push_back(3096576);
826  nStripsStack.push_back(3096576);
827  nStripsStack.push_back(3096576);
828  nStripsStack.push_back(1933312);
829  nStripsStack.push_back(1933312);
830  nStripsStack.push_back(1933312);
831  nStripsStack.push_back(1933312);
832  nStripsStack.push_back(1933312);
833  nStripsStack.push_back(1933312);
834  nStripsStack.push_back(1933312);
835  nStripsStack.push_back(1933312);
836  nStripsStack.push_back(1933312);
837  nStripsStack.push_back(1933312);
838  nStripsStack.push_back(1933312);
839  nStripsStack.push_back(1933312);
840  nStripsStack.push_back(1933312);
841  nStripsStack.push_back(1933312);
842  nStripsStack.push_back(1933312);
843  nStripsStack.push_back(1933312);
844  nStripsStack.push_back(1933312);
845  nStripsStack.push_back(1933312);
846 
847  std::vector<unsigned int> nStripsNoStack;
848  nStripsNoStack.push_back(9316352);
849  nStripsNoStack.push_back(1787904);
850  nStripsNoStack.push_back(565248);
851  nStripsNoStack.push_back(3096576);
852  nStripsNoStack.push_back(3866624);
853  nStripsNoStack.push_back(516096);
854  nStripsNoStack.push_back(663552);
855  nStripsNoStack.push_back(276480);
856  nStripsNoStack.push_back(331776);
857  nStripsNoStack.push_back(94208);
858  nStripsNoStack.push_back(94208);
859  nStripsNoStack.push_back(94208);
860  nStripsNoStack.push_back(94208);
861  nStripsNoStack.push_back(94208);
862  nStripsNoStack.push_back(94208);
863  nStripsNoStack.push_back(516096);
864  nStripsNoStack.push_back(589824);
865  nStripsNoStack.push_back(331776);
866  nStripsNoStack.push_back(368640);
867  nStripsNoStack.push_back(608256);
868  nStripsNoStack.push_back(681984);
869  nStripsNoStack.push_back(253952);
870  nStripsNoStack.push_back(253952);
871  nStripsNoStack.push_back(253952);
872  nStripsNoStack.push_back(217088);
873  nStripsNoStack.push_back(217088);
874  nStripsNoStack.push_back(217088);
875  nStripsNoStack.push_back(180224);
876  nStripsNoStack.push_back(180224);
877  nStripsNoStack.push_back(159744);
878  nStripsNoStack.push_back(253952);
879  nStripsNoStack.push_back(253952);
880  nStripsNoStack.push_back(253952);
881  nStripsNoStack.push_back(217088);
882  nStripsNoStack.push_back(217088);
883  nStripsNoStack.push_back(217088);
884  nStripsNoStack.push_back(180224);
885  nStripsNoStack.push_back(180224);
886  nStripsNoStack.push_back(159744);
887 
888  for (unsigned int i = 0; i < tkParts.size(); i++) {
889  modulesStackNormFactors[tkParts[i]] = nModulesStack[i];
890  modulesNoStackNormFactors[tkParts[i]] = nModulesNoStack[i];
891  fibersStackNormFactors[tkParts[i]] = nFibersStack[i];
892  fibersNoStackNormFactors[tkParts[i]] = nFibersNoStack[i];
893  APVsStackNormFactors[tkParts[i]] = nAPVsStack[i];
894  APVsNoStackNormFactors[tkParts[i]] = nAPVsNoStack[i];
895  stripsStackNormFactors[tkParts[i]] = nStripsStack[i];
896  stripsNoStackNormFactors[tkParts[i]] = nStripsNoStack[i];
897  }
898 
899  // for(std::map< std::string, unsigned int>::iterator it = allStripsTK.begin(); it != allStripsTK.end(); it++)
900  // {
901  // std::cout << it->first.c_str() << " " << it->second << std::endl;
902  // }
903 }
904 
905 double findNormFactor(const std::string currentPlotType, const std::string currentPart, const bool stackOption) {
906  // std::cout << "findNormFactor(): Finding normalization factor for this part: \"" << currentPart.c_str() << "\".\n";
907  // std::cout << " Plot type is: \"" << currentPlotType.c_str() << "\".\n";
908  // std::cout << " stack option is: " << stackOption << std::endl;
909  if (stackOption) {
910  if (currentPlotType == "BadModules") {
911  return modulesStackNormFactors[currentPart];
912  } else if (currentPlotType == "BadFibers") {
913  return fibersStackNormFactors[currentPart];
914  } else if (currentPlotType == "BadAPVs") {
915  return APVsStackNormFactors[currentPart];
916  } else if (currentPlotType == "AllBadStrips" || currentPlotType == "BadStripsFromAPVs" ||
917  currentPlotType == "BadStrips") {
918  return stripsStackNormFactors[currentPart];
919  } else {
920  std::cout << "findNormFactor(): ERROR! Requested a non supported plot type: " << currentPlotType.c_str()
921  << std::endl;
922  std::cout << " Add this to the function body or correct the error\n";
923  return 0; // This should trigger a divByZero error...
924  }
925  } else {
926  if (currentPlotType == "BadModules") {
927  return modulesNoStackNormFactors[currentPart];
928  } else if (currentPlotType == "BadFibers") {
929  return fibersNoStackNormFactors[currentPart];
930  } else if (currentPlotType == "BadAPVs") {
931  return APVsNoStackNormFactors[currentPart];
932  } else if (currentPlotType == "BadStrips" || currentPlotType == "BadStripsFromAPVs" ||
933  currentPlotType == "AllBadStrips") {
934  return stripsNoStackNormFactors[currentPart];
935  } else {
936  std::cout << "findNormFactor(): ERROR! Requested a non supported plot type: \"" << currentPlotType.c_str()
937  << "\"\n";
938  std::cout << " Add this to the function body or correct the error otherwise.\n";
939  return 0; // This should trigger a divByZero error...
940  }
941  }
942 }
static std::map< std::string, unsigned int > modulesStackNormFactors
Definition: makePlots.cc:26
int main(int argc, char *argv[])
Definition: makePlots.cc:35
double findNormFactor(const std::string currentPlotType, const std::string currentPart, const bool stackOption)
Definition: makePlots.cc:905
static std::map< std::string, unsigned int > APVsNoStackNormFactors
Definition: makePlots.cc:31
void fillNormFactorMaps()
Definition: makePlots.cc:514
assert(be >=bs)
static std::map< std::string, unsigned int > stripsNoStackNormFactors
Definition: makePlots.cc:33
void setCanvasStyle(TCanvas *c, const bool logScale)
Definition: makePlots.cc:440
static std::map< std::string, unsigned int > stripsStackNormFactors
Definition: makePlots.cc:32
void setLegendStyle(TLegend *l, const unsigned int nColumns)
Definition: makePlots.cc:488
static std::map< std::string, unsigned int > fibersStackNormFactors
Definition: makePlots.cc:28
void setPaveTextStyle(TPaveText *t, const bool isHorizontal=true)
Definition: makePlots.cc:493
static std::map< std::string, unsigned int > modulesNoStackNormFactors
Definition: makePlots.cc:27
void makePlots(std::string inputFileName, std::string outputFileName)
Definition: makePlots.cc:55
void setHistoStyle(TH1 *h)
Definition: makePlots.cc:451
static std::map< std::string, unsigned int > fibersNoStackNormFactors
Definition: makePlots.cc:29
void setHistoStackStyle(TH1 *h, const unsigned int lineColor)
Definition: makePlots.cc:467
static std::map< std::string, unsigned int > APVsStackNormFactors
Definition: makePlots.cc:30
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4