Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-10-08 05:11:29

0001 #include <iostream>
0002 #include <vector>
0003 #include <string>
0004 #include <map>
0005 #include <sstream>
0006 
0007 #include "TFile.h"
0008 #include "TH1F.h"
0009 #include "THStack.h"
0010 #include "TCanvas.h"
0011 #include "TLegend.h"
0012 #include "TPaveText.h"
0013 #include "TSystem.h"
0014 
0015 void setCanvasStyle(TCanvas* c, const bool logScale);
0016 void setHistoStyle(TH1* h);
0017 void setHistoStackStyle(TH1* h, const unsigned int lineColor);
0018 void setLegendStyle(TLegend* l, const unsigned int nColumns);
0019 void setPaveTextStyle(TPaveText* t, const bool isHorizontal = true);
0020 void fillNormFactorMaps();
0021 double findNormFactor(const std::string currentPlotType, const std::string currentPart, const bool stackOption);
0022 void makePlots(std::string inputFileName, std::string outputFileName);
0023 
0024 // Map with <name of tracker part, count of channels in the part>
0025 // It is static so it can be read by the functions fillNormFactorMaps() and findNormFactor(...)
0026 static std::map<std::string, unsigned int> modulesStackNormFactors;
0027 static std::map<std::string, unsigned int> modulesNoStackNormFactors;
0028 static std::map<std::string, unsigned int> fibersStackNormFactors;
0029 static std::map<std::string, unsigned int> fibersNoStackNormFactors;
0030 static std::map<std::string, unsigned int> APVsStackNormFactors;
0031 static std::map<std::string, unsigned int> APVsNoStackNormFactors;
0032 static std::map<std::string, unsigned int> stripsStackNormFactors;
0033 static std::map<std::string, unsigned int> stripsNoStackNormFactors;
0034 
0035 int main(int argc, char* argv[]) {
0036   if (argc == 3) {
0037     char* inputFileName = argv[1];
0038     char* outputFileName = argv[2];
0039 
0040     std::cout << "ready to make plots from " << inputFileName << " to " << outputFileName << std::endl;
0041 
0042     int returncode = 0;
0043     makePlots(inputFileName, outputFileName);
0044 
0045     return returncode;
0046 
0047   } else {
0048     std::cout << "Too few arguments: " << argc << std::endl;
0049     return -1;
0050   }
0051 
0052   return -9;
0053 }
0054 
0055 void makePlots(std::string inputFileName, std::string outputFileName) {
0056   //
0057 
0058   // Open input and output file
0059   TFile* inputFile = new TFile(inputFileName.c_str(), "READ");
0060   TFile* outputFile = new TFile(outputFileName.c_str(), "RECREATE");
0061 
0062   std::ostringstream oss;
0063   std::string histoName;
0064   std::vector<std::string> plotType;
0065   plotType.push_back("BadModules");
0066   plotType.push_back("BadFibers");
0067   plotType.push_back("BadAPVs");
0068   plotType.push_back("BadStrips");
0069   plotType.push_back("BadStripsFromAPVs");
0070   plotType.push_back("AllBadStrips");
0071   std::vector<std::string> subDetName;
0072   subDetName.push_back("TIB");
0073   subDetName.push_back("TID+");
0074   subDetName.push_back("TID-");
0075   subDetName.push_back("TOB");
0076   subDetName.push_back("TEC+");
0077   subDetName.push_back("TEC-");
0078 
0079   // Standard plot options for THStack and for standalone histograms
0080   const bool stackHistograms = true;
0081 
0082   std::string plotStackOptions;
0083   if (stackHistograms)
0084     plotStackOptions = "";
0085   else
0086     plotStackOptions = "nostack p";
0087   const std::string plotHistoOptions("p");
0088 
0089   // Defer the filling of the normFactor maps to this function
0090   // Conceptually trivial but lengthy
0091   fillNormFactorMaps();
0092 
0093   //   // Finds number of channels from above map
0094   //   std::string completePartName = subDetName;
0095   //   if(partName.compare("") != 0)
0096   //     completePartName += " " + partName;
0097   //   if(partNumber != 0)
0098   //   {
0099   //     oss.str("");
0100   //     oss << partNumber;
0101   //     completePartName += " " + oss.str();
0102   //   }
0103   //
0104   //   // Total number of channels in currently processed map
0105   //   const unsigned int nModulesInPart = allModulesTK[completePartName.c_str()];
0106   //   const unsigned int nFibersInPart = allFibersTK[completePartName.c_str()];
0107   //   const unsigned int nAPVsInPart = allAPVsTK[completePartName.c_str()];
0108   //   const unsigned int nStripsInPart = allStripsTK[completePartName.c_str()];
0109 
0110   TH1F* hTracker;
0111   TH1F* hTIB;
0112   TH1F* hTID;
0113   TH1F* hTOB;
0114   TH1F* hTEC;
0115   TH1F* histo;
0116   TH1F* histo2;
0117   THStack* histoStack;
0118   TCanvas* c1;
0119   TCanvas* c2;
0120   TLegend* legend;
0121   TLegend* legend2;
0122   TPaveText* textX;
0123   TPaveText* textY;
0124   std::string entryLabel;
0125   std::vector<TH1F*> hLayers;
0126   //   unsigned int normFactor;
0127 
0128   for (std::vector<std::string>::iterator itPlot = plotType.begin(); itPlot != plotType.end(); itPlot++) {
0129     // Put together the Tracker histograms with the TIB, TID, TOB and TEC ones
0130 
0131     histoName = "h" + *itPlot + "Tracker";
0132     hTracker = (TH1F*)inputFile->Get(histoName.c_str());
0133     if (hTracker) {
0134       hTracker->Scale(1 / findNormFactor(*itPlot, "Tracker", stackHistograms));
0135     } else {
0136       std::cout << histoName << " not found" << std::endl;
0137     }
0138 
0139     histoStack = new THStack(histoName.c_str(), histoName.c_str());
0140 
0141     histoName = "h" + *itPlot + "TIB";
0142     hTIB = (TH1F*)inputFile->Get(histoName.c_str());
0143     if (hTIB) {
0144       hTIB->Scale(1 / findNormFactor(*itPlot, "TIB", stackHistograms));
0145     } else {
0146       std::cout << histoName << " not found" << std::endl;
0147     }
0148 
0149     histoName = "h" + *itPlot + "TID";
0150     hTID = (TH1F*)inputFile->Get(histoName.c_str());
0151     if (hTID) {
0152       hTID->Scale(1 / findNormFactor(*itPlot, "TID", stackHistograms));
0153     } else {
0154       std::cout << histoName << " not found" << std::endl;
0155     }
0156 
0157     histoName = "h" + *itPlot + "TOB";
0158     hTOB = (TH1F*)inputFile->Get(histoName.c_str());
0159     if (hTOB) {
0160       hTOB->Scale(1 / findNormFactor(*itPlot, "TOB", stackHistograms));
0161     } else {
0162       std::cout << histoName << " not found" << std::endl;
0163     }
0164 
0165     histoName = "h" + *itPlot + "TEC";
0166     hTEC = (TH1F*)inputFile->Get(histoName.c_str());
0167     if (hTEC) {
0168       hTEC->Scale(1 / findNormFactor(*itPlot, "TEC", stackHistograms));
0169     } else {
0170       std::cout << histoName << " not found" << std::endl;
0171     }
0172 
0173     c1 = new TCanvas(("c" + *itPlot + "Tracker").c_str(), "", 1200, 600);
0174     setCanvasStyle(c1, false);
0175     //     hTracker->Draw();
0176     if (hTracker)
0177       setHistoStackStyle(hTracker, 1);
0178     //     hTIB->Draw("same");
0179     if (hTIB) {
0180       setHistoStackStyle(hTIB, 2);
0181       histoStack->Add(hTIB);
0182     }
0183     //     hTID->Draw("same");
0184     if (hTID) {
0185       setHistoStackStyle(hTID, 3);
0186       histoStack->Add(hTID);
0187     }
0188     //     hTOB->Draw("same");
0189     if (hTOB) {
0190       setHistoStackStyle(hTOB, 4);
0191       histoStack->Add(hTOB);
0192     }
0193     //     hTEC->Draw("same");
0194     if (hTEC) {
0195       setHistoStackStyle(hTEC, 6);
0196       histoStack->Add(hTEC);
0197     }
0198     // Bug in ROOT? If we plot a stack with the "stack" option and the Y axis is in log scale,
0199     // but there are no entries in any of the histograms of the stack, then ROOT crashes
0200     // Workaround: at this stage, check that at least one histogram has >0 entries,
0201     // otherwise, switch back to linear Y scale
0202     // Curiously, there is no crash if the "nostack" option is chosen...
0203     double histoStackMaximum = histoStack->GetMaximum();
0204     if (histoStackMaximum == 0) {
0205       c1->SetLogy(0);
0206     }
0207     histoStack->Draw(plotStackOptions.c_str());
0208     if (histoStack->GetYaxis()->GetXmax() > 0.9)
0209       histoStack->GetYaxis()->SetRangeUser(0., 0.1);
0210     legend = new TLegend(0.4, 0.9, 0.9, 1);
0211     legend->AddEntry(hTIB, "TIB");
0212     legend->AddEntry(hTID, "TID");
0213     legend->AddEntry(hTOB, "TOB");
0214     legend->AddEntry(hTEC, "TEC");
0215     setLegendStyle(legend, 2);
0216     legend->Draw();
0217     textX = new TPaveText();
0218     textY = new TPaveText();
0219     setPaveTextStyle(textX);
0220     setPaveTextStyle(textY, false);
0221     textX->Draw();
0222     textY->Draw();
0223     gSystem->ProcessEvents();
0224     c1->Update();
0225     outputFile->cd();
0226     c1->Write();
0227     c1->SaveAs((*itPlot + "Tracker.png").c_str());
0228 
0229     delete histoStack;
0230     delete textX;
0231     delete textY;
0232     delete legend;
0233     legend = nullptr;
0234     delete c1;
0235 
0236     // Put together the histograms for the different layers of the detectors
0237     for (std::vector<std::string>::iterator itSub = subDetName.begin(); itSub != subDetName.end(); itSub++) {
0238       unsigned int nLayers = 0;
0239       std::string layerName;
0240       //       std::cout << "itSub = " << (*itSub).c_str() << std::endl;
0241       if ((*itSub) == "TIB") {
0242         nLayers = 4;
0243         layerName = "Layer";
0244         legend = new TLegend(0.4, 0.9, 0.9, 1);
0245         setLegendStyle(legend, 2);
0246       } else if ((*itSub) == "TID+" || (*itSub) == "TID-") {
0247         nLayers = 3;
0248         layerName = "Disk";
0249         legend = new TLegend(0.35, 0.9, 0.9, 1);
0250         setLegendStyle(legend, 3);
0251       } else if ((*itSub) == "TOB") {
0252         nLayers = 6;
0253         layerName = "Layer";
0254         legend = new TLegend(0.35, 0.9, 0.9, 1);
0255         setLegendStyle(legend, 3);
0256       } else if ((*itSub) == "TEC+" || (*itSub) == "TEC-") {
0257         nLayers = 9;
0258         layerName = "Disk";
0259         legend = new TLegend(0.35, 0.9, 1, 1);
0260         setLegendStyle(legend, 5);
0261       }
0262 
0263       if (legend == nullptr) {
0264         std::cerr << "Unknown itSub type " << *itSub << std::endl;
0265         assert(false);
0266       }
0267 
0268       c1 = new TCanvas(("c" + *itPlot + *itSub).c_str(), "", 1200, 600);
0269       setCanvasStyle(c1, false);
0270       //       if((*itSub).compare("TEC+")==0 || (*itSub).compare("TEC-")==0)
0271       //       {
0272       //         histoName = "h" + *itPlot + "TEC";
0273       //       }
0274       //       else
0275       //       {
0276       histoName = "h" + *itPlot + *itSub;
0277       //       }
0278       //       hSubDet = (TH1F*)inputFile->Get(histoName.c_str());
0279       //       setHistoStackStyle(hSubDet,1);
0280       //       hSubDet->Draw();
0281       histoStack = new THStack(histoName.c_str(), histoName.c_str());
0282 
0283       for (unsigned int iLayer = 1; iLayer <= nLayers; iLayer++) {
0284         oss.str("");
0285         oss << iLayer;
0286         //         // TIB and TOB have no plus/minus side division
0287         //         // While TEC has it but I plot them separately
0288         //         // TID has plus/minus side division but I plot everything in a single plot
0289         //         if((*itSub).compare("TID")==0)
0290         //         {
0291         //           histoName = "h" + *itPlot + *itSub + "+" + layerName + oss.str();
0292         // //           std::cout << "histoName = " << histoName.c_str() << std::endl;
0293         //           histo = (TH1F*)inputFile->Get(histoName.c_str());
0294         //
0295         //           // First: plot histogram separately
0296         //           setHistoStyle(histo);
0297         //           c2 = new TCanvas(("c" + *itPlot + *itSub + "+" + layerName + oss.str()).c_str(), "", 1200, 600);
0298         //           setCanvasStyle(c2, true);
0299         //           histo->Draw(plotHistoOptions.c_str());
0300         //           legend2 = new TLegend(0.6,0.92,0.9,0.97);
0301         //           legend2->AddEntry(histo,(*itSub + "+" + layerName + oss.str()).c_str());
0302         //           setLegendStyle(legend2, 1);
0303         //           legend2->Draw();
0304         //           gSystem->ProcessEvents();
0305         //           c2->Update();
0306         //           outputFile->cd();
0307         //           c2->Write();
0308         //           c2->SaveAs((*itPlot + *itSub + "+" + layerName + oss.str() + ".png").c_str());
0309         //           delete legend2;
0310         //           delete c2;
0311         //
0312         //           // Second: add histogram to THStack
0313         //           hLayers.push_back(histo);
0314         //           setHistoStackStyle(hLayers.back(), iLayer);
0315         //           histoStack->Add(hLayers.back());
0316         //           entryLabel = *itSub + "+ " + layerName + " " + oss.str();
0317         //           legend->AddEntry(hLayers.back(), entryLabel.c_str());
0318         //
0319         //
0320         //           histoName = "h" + *itPlot + *itSub + "-" + layerName + oss.str();
0321         // //           std::cout << "histoName = " << histoName.c_str() << std::endl;
0322         //           histo = (TH1F*)inputFile->Get(histoName.c_str());
0323         //
0324         //           // First: plot histogram separately
0325         //           setHistoStyle(histo);
0326         //           c2 = new TCanvas(("c" + *itPlot + *itSub + "-" + layerName + oss.str()).c_str(), "", 1200, 600);
0327         //           setCanvasStyle(c2, true);
0328         //           histo->Draw(plotHistoOptions.c_str());
0329         //           legend2 = new TLegend(0.6,0.92,0.9,0.97);
0330         //           legend2->AddEntry(histo,(*itSub + "-" + layerName + oss.str()).c_str());
0331         //           setLegendStyle(legend2, 1);
0332         //           legend2->Draw();
0333         //           gSystem->ProcessEvents();
0334         //           c2->Update();
0335         //           outputFile->cd();
0336         //           c2->Write();
0337         //           c2->SaveAs((*itPlot + *itSub + "-" + layerName + oss.str() + ".png").c_str());
0338         //           delete legend2;
0339         //           delete c2;
0340         //
0341         //           // Second: add histogram to THStack
0342         //           hLayers.push_back(histo);
0343         //           setHistoStackStyle(hLayers.back(), iLayer+nLayers);
0344         //           histoStack->Add(hLayers.back());
0345         //           entryLabel = *itSub + "- " + layerName + " " + oss.str();
0346         //           legend->AddEntry(hLayers.back(), entryLabel.c_str());
0347         // //           hLayers.back()->Draw("same");
0348         //         }
0349         //         else
0350         //         {
0351         histoName = "h" + *itPlot + *itSub + layerName + oss.str();
0352         //           std::cout << "histoName = " << histoName.c_str() << std::endl;
0353         histo = (TH1F*)inputFile->Get(histoName.c_str());
0354         if (histo) {
0355           histo2 = new TH1F(*histo);
0356           histo->Scale(1 / findNormFactor(*itPlot, *itSub + " " + layerName + " " + oss.str(), false));
0357           histo2->Scale(1 / findNormFactor(*itPlot, *itSub + " " + layerName + " " + oss.str(), stackHistograms));
0358           // First: plot histogram separately
0359           setHistoStyle(histo);
0360           c2 = new TCanvas(("c" + *itPlot + *itSub + layerName + oss.str()).c_str(), "", 1200, 600);
0361           setCanvasStyle(c2, true);
0362           histo->Draw(plotHistoOptions.c_str());
0363           double histoMaximum = histo->GetMaximum();
0364           // Otherwise it does not draw the pad
0365           if (histoMaximum == 0) {
0366             c2->SetLogy(0);
0367           }
0368           legend2 = new TLegend(0.6, 0.92, 0.9, 0.97);
0369           legend2->AddEntry(histo, (*itSub + layerName + oss.str()).c_str());
0370           setLegendStyle(legend2, 1);
0371           legend2->Draw();
0372           textX = new TPaveText();
0373           textY = new TPaveText();
0374           setPaveTextStyle(textX);
0375           setPaveTextStyle(textY, false);
0376           textX->Draw();
0377           textY->Draw();
0378           gSystem->ProcessEvents();
0379           c2->Update();
0380           outputFile->cd();
0381           c2->Write();
0382           c2->SaveAs((*itPlot + *itSub + layerName + oss.str() + ".png").c_str());
0383           delete textX;
0384           delete textY;
0385           delete legend2;
0386           delete c2;
0387 
0388           // Second: add histogram to THStack
0389           setHistoStyle(histo2);
0390           hLayers.push_back(histo2);
0391           setHistoStackStyle(hLayers.back(), iLayer);
0392           histoStack->Add(hLayers.back());
0393           entryLabel = *itSub + " " + layerName + " " + oss.str();
0394           legend->AddEntry(hLayers.back(), entryLabel.c_str());
0395           //           hLayers.back()->Draw("same");
0396           //         }
0397           //           delete histo2;
0398         } else {
0399           std::cout << histoName << " not found" << std::endl;
0400         }
0401       }
0402       histoStack->Draw(plotStackOptions.c_str());
0403 
0404       // Bug in ROOT? If we plot a stack with the "stack" option and the Y axis is in log scale,
0405       // but there are no entries in any of the histograms of the stack, then ROOT crashes
0406       // Workaround: at this stage, check that at least one histogram has >0 entries,
0407       // otherwise, switch back to linear Y scale
0408       // Curiously, there is no crash if the "nostack" option is chosen...
0409       double histoStackMaximum = histoStack->GetMaximum();
0410       if (histoStackMaximum == 0) {
0411         c1->SetLogy(0);
0412       }
0413       if (histoStackMaximum > 0.01)
0414         histoStack->SetMaximum(0.01);
0415       textX = new TPaveText();
0416       textY = new TPaveText();
0417       setPaveTextStyle(textX);
0418       setPaveTextStyle(textY, false);
0419       textX->Draw();
0420       textY->Draw();
0421       legend->Draw();
0422       gSystem->ProcessEvents();
0423       c1->Update();
0424       outputFile->cd();
0425       c1->Write();
0426       c1->SaveAs((*itPlot + *itSub + ".png").c_str());
0427       delete histoStack;
0428       delete textX;
0429       delete textY;
0430       delete legend;
0431       legend = nullptr;
0432       delete c1;
0433     }
0434   }
0435 
0436   inputFile->Close();
0437   outputFile->Close();
0438 }
0439 
0440 void setCanvasStyle(TCanvas* c, const bool logScale) {
0441   c->SetFillColor(0);
0442   c->SetFrameBorderMode(0);
0443   if (logScale)
0444     c->SetLogy(1);
0445   else
0446     c->SetLogy(0);
0447   c->SetCanvasSize(1200, 600);
0448   c->SetWindowSize(1200, 600);
0449 }
0450 
0451 void setHistoStyle(TH1* h) {
0452   h->SetLineStyle(0);
0453   h->SetMarkerStyle(3);
0454   h->SetMarkerSize(1);
0455   h->SetMarkerColor(1);
0456   h->SetStats(kFALSE);
0457   //   h->GetYaxis()->SetTitle("Fraction of total");
0458   //   h->GetXaxis()->SetTitle("IOV");
0459   //   h->GetXaxis()->SetTitleOffset(-0.3);
0460   // Avoid having too many bins with labels
0461   if (h->GetNbinsX() > 25)
0462     for (int i = 1; i < h->GetNbinsX() - 1; i++)
0463       if ((i % (h->GetNbinsX() / 25 + 1)))
0464         h->GetXaxis()->SetBinLabel(i + 1, "");
0465 }
0466 
0467 void setHistoStackStyle(TH1* h, const unsigned int lineColor) {
0468   h->SetLineStyle(0);
0469   //   h->SetDrawOption("e1p");
0470   // Best marker types are 20-23 - use them with different colors
0471   h->SetMarkerStyle(lineColor % 4 + 20);
0472   h->SetMarkerSize(1);
0473   h->SetMarkerColor(lineColor);
0474   h->SetLineColor(lineColor);
0475   h->SetFillColor(lineColor);
0476   h->SetLineWidth(2);
0477   h->SetStats(kFALSE);
0478   h->GetYaxis()->SetTitle("Fraction of total");
0479   //   h->GetXaxis()->SetTitle("IOV");
0480   //   h->GetXaxis()->SetTitleOffset(1.2);
0481   // Avoid having too many bins with labels
0482   if (h->GetNbinsX() > 25)
0483     for (int i = 1; i < h->GetNbinsX() - 1; i++)
0484       if (i % (h->GetNbinsX() / 25 + 1))
0485         h->GetXaxis()->SetBinLabel(i + 1, "");
0486 }
0487 
0488 void setLegendStyle(TLegend* l, const unsigned int nColumns) {
0489   l->SetNColumns(nColumns);
0490   l->SetFillColor(0);
0491 }
0492 
0493 void setPaveTextStyle(TPaveText* t, const bool isHorizontal) {
0494   t->SetLineStyle(0);
0495   t->SetFillColor(0);
0496   t->SetFillStyle(0);
0497   t->SetBorderSize(0);
0498   if (isHorizontal) {
0499     t->SetX1NDC(0.905);
0500     t->SetX2NDC(0.975);
0501     t->SetY1NDC(0.062);
0502     t->SetY2NDC(0.095);
0503     t->AddText("IOV");
0504   } else {
0505     t->SetX1NDC(0.03);
0506     t->SetX2NDC(0.05);
0507     t->SetY1NDC(0.33);
0508     t->SetY2NDC(0.68);
0509     TText* t1 = t->AddText("Fraction of total");
0510     t1->SetTextAngle(90.);
0511   }
0512 }
0513 
0514 void fillNormFactorMaps() {
0515   // Number of modules, fibers, APVs, strips for each tracker part
0516   std::vector<std::string> tkParts;
0517   tkParts.push_back("Tracker");
0518   tkParts.push_back("TIB");
0519   tkParts.push_back("TID");
0520   tkParts.push_back("TOB");
0521   tkParts.push_back("TEC");
0522   tkParts.push_back("TIB Layer 1");
0523   tkParts.push_back("TIB Layer 2");
0524   tkParts.push_back("TIB Layer 3");
0525   tkParts.push_back("TIB Layer 4");
0526   tkParts.push_back("TID+ Disk 1");
0527   tkParts.push_back("TID+ Disk 2");
0528   tkParts.push_back("TID+ Disk 3");
0529   tkParts.push_back("TID- Disk 1");
0530   tkParts.push_back("TID- Disk 2");
0531   tkParts.push_back("TID- Disk 3");
0532   tkParts.push_back("TOB Layer 1");
0533   tkParts.push_back("TOB Layer 2");
0534   tkParts.push_back("TOB Layer 3");
0535   tkParts.push_back("TOB Layer 4");
0536   tkParts.push_back("TOB Layer 5");
0537   tkParts.push_back("TOB Layer 6");
0538   tkParts.push_back("TEC+ Disk 1");
0539   tkParts.push_back("TEC+ Disk 2");
0540   tkParts.push_back("TEC+ Disk 3");
0541   tkParts.push_back("TEC+ Disk 4");
0542   tkParts.push_back("TEC+ Disk 5");
0543   tkParts.push_back("TEC+ Disk 6");
0544   tkParts.push_back("TEC+ Disk 7");
0545   tkParts.push_back("TEC+ Disk 8");
0546   tkParts.push_back("TEC+ Disk 9");
0547   tkParts.push_back("TEC- Disk 1");
0548   tkParts.push_back("TEC- Disk 2");
0549   tkParts.push_back("TEC- Disk 3");
0550   tkParts.push_back("TEC- Disk 4");
0551   tkParts.push_back("TEC- Disk 5");
0552   tkParts.push_back("TEC- Disk 6");
0553   tkParts.push_back("TEC- Disk 7");
0554   tkParts.push_back("TEC- Disk 8");
0555   tkParts.push_back("TEC- Disk 9");
0556 
0557   std::vector<unsigned int> nModulesStack;
0558   nModulesStack.push_back(15148);
0559   nModulesStack.push_back(15148);
0560   nModulesStack.push_back(15148);
0561   nModulesStack.push_back(15148);
0562   nModulesStack.push_back(15148);
0563   nModulesStack.push_back(2724);
0564   nModulesStack.push_back(2724);
0565   nModulesStack.push_back(2724);
0566   nModulesStack.push_back(2724);
0567   nModulesStack.push_back(408);
0568   nModulesStack.push_back(408);
0569   nModulesStack.push_back(408);
0570   nModulesStack.push_back(408);
0571   nModulesStack.push_back(408);
0572   nModulesStack.push_back(408);
0573   nModulesStack.push_back(5208);
0574   nModulesStack.push_back(5208);
0575   nModulesStack.push_back(5208);
0576   nModulesStack.push_back(5208);
0577   nModulesStack.push_back(5208);
0578   nModulesStack.push_back(5208);
0579   nModulesStack.push_back(3200);
0580   nModulesStack.push_back(3200);
0581   nModulesStack.push_back(3200);
0582   nModulesStack.push_back(3200);
0583   nModulesStack.push_back(3200);
0584   nModulesStack.push_back(3200);
0585   nModulesStack.push_back(3200);
0586   nModulesStack.push_back(3200);
0587   nModulesStack.push_back(3200);
0588   nModulesStack.push_back(3200);
0589   nModulesStack.push_back(3200);
0590   nModulesStack.push_back(3200);
0591   nModulesStack.push_back(3200);
0592   nModulesStack.push_back(3200);
0593   nModulesStack.push_back(3200);
0594   nModulesStack.push_back(3200);
0595   nModulesStack.push_back(3200);
0596   nModulesStack.push_back(3200);
0597 
0598   std::vector<unsigned int> nModulesNoStack;
0599   nModulesNoStack.push_back(15148);
0600   nModulesNoStack.push_back(2724);
0601   nModulesNoStack.push_back(816);
0602   nModulesNoStack.push_back(5208);
0603   nModulesNoStack.push_back(6400);
0604   nModulesNoStack.push_back(672);
0605   nModulesNoStack.push_back(864);
0606   nModulesNoStack.push_back(540);
0607   nModulesNoStack.push_back(648);
0608   nModulesNoStack.push_back(136);
0609   nModulesNoStack.push_back(136);
0610   nModulesNoStack.push_back(136);
0611   nModulesNoStack.push_back(136);
0612   nModulesNoStack.push_back(136);
0613   nModulesNoStack.push_back(136);
0614   nModulesNoStack.push_back(1008);
0615   nModulesNoStack.push_back(1152);
0616   nModulesNoStack.push_back(648);
0617   nModulesNoStack.push_back(720);
0618   nModulesNoStack.push_back(792);
0619   nModulesNoStack.push_back(888);
0620   nModulesNoStack.push_back(408);
0621   nModulesNoStack.push_back(408);
0622   nModulesNoStack.push_back(408);
0623   nModulesNoStack.push_back(360);
0624   nModulesNoStack.push_back(360);
0625   nModulesNoStack.push_back(360);
0626   nModulesNoStack.push_back(312);
0627   nModulesNoStack.push_back(312);
0628   nModulesNoStack.push_back(272);
0629   nModulesNoStack.push_back(408);
0630   nModulesNoStack.push_back(408);
0631   nModulesNoStack.push_back(408);
0632   nModulesNoStack.push_back(360);
0633   nModulesNoStack.push_back(360);
0634   nModulesNoStack.push_back(360);
0635   nModulesNoStack.push_back(312);
0636   nModulesNoStack.push_back(312);
0637   nModulesNoStack.push_back(272);
0638 
0639   //
0640   std::vector<unsigned int> nFibersStack;
0641   nFibersStack.push_back(36392);
0642   nFibersStack.push_back(36392);
0643   nFibersStack.push_back(36392);
0644   nFibersStack.push_back(36392);
0645   nFibersStack.push_back(36392);
0646   nFibersStack.push_back(6984);
0647   nFibersStack.push_back(6984);
0648   nFibersStack.push_back(6984);
0649   nFibersStack.push_back(6984);
0650   nFibersStack.push_back(1104);
0651   nFibersStack.push_back(1104);
0652   nFibersStack.push_back(1104);
0653   nFibersStack.push_back(1104);
0654   nFibersStack.push_back(1104);
0655   nFibersStack.push_back(1104);
0656   nFibersStack.push_back(12096);
0657   nFibersStack.push_back(12096);
0658   nFibersStack.push_back(12096);
0659   nFibersStack.push_back(12096);
0660   nFibersStack.push_back(12096);
0661   nFibersStack.push_back(12096);
0662   nFibersStack.push_back(7552);
0663   nFibersStack.push_back(7552);
0664   nFibersStack.push_back(7552);
0665   nFibersStack.push_back(7552);
0666   nFibersStack.push_back(7552);
0667   nFibersStack.push_back(7552);
0668   nFibersStack.push_back(7552);
0669   nFibersStack.push_back(7552);
0670   nFibersStack.push_back(7552);
0671   nFibersStack.push_back(7552);
0672   nFibersStack.push_back(7552);
0673   nFibersStack.push_back(7552);
0674   nFibersStack.push_back(7552);
0675   nFibersStack.push_back(7552);
0676   nFibersStack.push_back(7552);
0677   nFibersStack.push_back(7552);
0678   nFibersStack.push_back(7552);
0679   nFibersStack.push_back(7552);
0680 
0681   std::vector<unsigned int> nFibersNoStack;
0682   nFibersNoStack.push_back(36392);
0683   nFibersNoStack.push_back(6984);
0684   nFibersNoStack.push_back(2208);
0685   nFibersNoStack.push_back(12096);
0686   nFibersNoStack.push_back(15104);
0687   nFibersNoStack.push_back(2016);
0688   nFibersNoStack.push_back(2592);
0689   nFibersNoStack.push_back(1080);
0690   nFibersNoStack.push_back(1296);
0691   nFibersNoStack.push_back(368);
0692   nFibersNoStack.push_back(368);
0693   nFibersNoStack.push_back(368);
0694   nFibersNoStack.push_back(368);
0695   nFibersNoStack.push_back(368);
0696   nFibersNoStack.push_back(368);
0697   nFibersNoStack.push_back(2016);
0698   nFibersNoStack.push_back(2304);
0699   nFibersNoStack.push_back(1296);
0700   nFibersNoStack.push_back(1440);
0701   nFibersNoStack.push_back(2376);
0702   nFibersNoStack.push_back(2664);
0703   nFibersNoStack.push_back(992);
0704   nFibersNoStack.push_back(992);
0705   nFibersNoStack.push_back(992);
0706   nFibersNoStack.push_back(848);
0707   nFibersNoStack.push_back(848);
0708   nFibersNoStack.push_back(848);
0709   nFibersNoStack.push_back(704);
0710   nFibersNoStack.push_back(704);
0711   nFibersNoStack.push_back(624);
0712   nFibersNoStack.push_back(992);
0713   nFibersNoStack.push_back(992);
0714   nFibersNoStack.push_back(992);
0715   nFibersNoStack.push_back(848);
0716   nFibersNoStack.push_back(848);
0717   nFibersNoStack.push_back(848);
0718   nFibersNoStack.push_back(704);
0719   nFibersNoStack.push_back(704);
0720   nFibersNoStack.push_back(624);
0721 
0722   //
0723   std::vector<unsigned int> nAPVsStack;
0724   nAPVsStack.push_back(72784);
0725   nAPVsStack.push_back(72784);
0726   nAPVsStack.push_back(72784);
0727   nAPVsStack.push_back(72784);
0728   nAPVsStack.push_back(72784);
0729   nAPVsStack.push_back(13968);
0730   nAPVsStack.push_back(13968);
0731   nAPVsStack.push_back(13968);
0732   nAPVsStack.push_back(13968);
0733   nAPVsStack.push_back(2208);
0734   nAPVsStack.push_back(2208);
0735   nAPVsStack.push_back(2208);
0736   nAPVsStack.push_back(2208);
0737   nAPVsStack.push_back(2208);
0738   nAPVsStack.push_back(2208);
0739   nAPVsStack.push_back(24192);
0740   nAPVsStack.push_back(24192);
0741   nAPVsStack.push_back(24192);
0742   nAPVsStack.push_back(24192);
0743   nAPVsStack.push_back(24192);
0744   nAPVsStack.push_back(24192);
0745   nAPVsStack.push_back(15104);
0746   nAPVsStack.push_back(15104);
0747   nAPVsStack.push_back(15104);
0748   nAPVsStack.push_back(15104);
0749   nAPVsStack.push_back(15104);
0750   nAPVsStack.push_back(15104);
0751   nAPVsStack.push_back(15104);
0752   nAPVsStack.push_back(15104);
0753   nAPVsStack.push_back(15104);
0754   nAPVsStack.push_back(15104);
0755   nAPVsStack.push_back(15104);
0756   nAPVsStack.push_back(15104);
0757   nAPVsStack.push_back(15104);
0758   nAPVsStack.push_back(15104);
0759   nAPVsStack.push_back(15104);
0760   nAPVsStack.push_back(15104);
0761   nAPVsStack.push_back(15104);
0762   nAPVsStack.push_back(15104);
0763 
0764   std::vector<unsigned int> nAPVsNoStack;
0765   nAPVsNoStack.push_back(72784);
0766   nAPVsNoStack.push_back(13968);
0767   nAPVsNoStack.push_back(4416);
0768   nAPVsNoStack.push_back(24192);
0769   nAPVsNoStack.push_back(30208);
0770   nAPVsNoStack.push_back(4032);
0771   nAPVsNoStack.push_back(5184);
0772   nAPVsNoStack.push_back(2160);
0773   nAPVsNoStack.push_back(2592);
0774   nAPVsNoStack.push_back(736);
0775   nAPVsNoStack.push_back(736);
0776   nAPVsNoStack.push_back(736);
0777   nAPVsNoStack.push_back(736);
0778   nAPVsNoStack.push_back(736);
0779   nAPVsNoStack.push_back(736);
0780   nAPVsNoStack.push_back(4032);
0781   nAPVsNoStack.push_back(4608);
0782   nAPVsNoStack.push_back(2592);
0783   nAPVsNoStack.push_back(2880);
0784   nAPVsNoStack.push_back(4752);
0785   nAPVsNoStack.push_back(5328);
0786   nAPVsNoStack.push_back(1984);
0787   nAPVsNoStack.push_back(1984);
0788   nAPVsNoStack.push_back(1984);
0789   nAPVsNoStack.push_back(1696);
0790   nAPVsNoStack.push_back(1696);
0791   nAPVsNoStack.push_back(1696);
0792   nAPVsNoStack.push_back(1408);
0793   nAPVsNoStack.push_back(1408);
0794   nAPVsNoStack.push_back(1248);
0795   nAPVsNoStack.push_back(1984);
0796   nAPVsNoStack.push_back(1984);
0797   nAPVsNoStack.push_back(1984);
0798   nAPVsNoStack.push_back(1696);
0799   nAPVsNoStack.push_back(1696);
0800   nAPVsNoStack.push_back(1696);
0801   nAPVsNoStack.push_back(1408);
0802   nAPVsNoStack.push_back(1408);
0803   nAPVsNoStack.push_back(1248);
0804 
0805   //
0806   std::vector<unsigned int> nStripsStack;
0807   nStripsStack.push_back(9316352);
0808   nStripsStack.push_back(9316352);
0809   nStripsStack.push_back(9316352);
0810   nStripsStack.push_back(9316352);
0811   nStripsStack.push_back(9316352);
0812   nStripsStack.push_back(1787904);
0813   nStripsStack.push_back(1787904);
0814   nStripsStack.push_back(1787904);
0815   nStripsStack.push_back(1787904);
0816   nStripsStack.push_back(282624);
0817   nStripsStack.push_back(282624);
0818   nStripsStack.push_back(282624);
0819   nStripsStack.push_back(282624);
0820   nStripsStack.push_back(282624);
0821   nStripsStack.push_back(282624);
0822   nStripsStack.push_back(3096576);
0823   nStripsStack.push_back(3096576);
0824   nStripsStack.push_back(3096576);
0825   nStripsStack.push_back(3096576);
0826   nStripsStack.push_back(3096576);
0827   nStripsStack.push_back(3096576);
0828   nStripsStack.push_back(1933312);
0829   nStripsStack.push_back(1933312);
0830   nStripsStack.push_back(1933312);
0831   nStripsStack.push_back(1933312);
0832   nStripsStack.push_back(1933312);
0833   nStripsStack.push_back(1933312);
0834   nStripsStack.push_back(1933312);
0835   nStripsStack.push_back(1933312);
0836   nStripsStack.push_back(1933312);
0837   nStripsStack.push_back(1933312);
0838   nStripsStack.push_back(1933312);
0839   nStripsStack.push_back(1933312);
0840   nStripsStack.push_back(1933312);
0841   nStripsStack.push_back(1933312);
0842   nStripsStack.push_back(1933312);
0843   nStripsStack.push_back(1933312);
0844   nStripsStack.push_back(1933312);
0845   nStripsStack.push_back(1933312);
0846 
0847   std::vector<unsigned int> nStripsNoStack;
0848   nStripsNoStack.push_back(9316352);
0849   nStripsNoStack.push_back(1787904);
0850   nStripsNoStack.push_back(565248);
0851   nStripsNoStack.push_back(3096576);
0852   nStripsNoStack.push_back(3866624);
0853   nStripsNoStack.push_back(516096);
0854   nStripsNoStack.push_back(663552);
0855   nStripsNoStack.push_back(276480);
0856   nStripsNoStack.push_back(331776);
0857   nStripsNoStack.push_back(94208);
0858   nStripsNoStack.push_back(94208);
0859   nStripsNoStack.push_back(94208);
0860   nStripsNoStack.push_back(94208);
0861   nStripsNoStack.push_back(94208);
0862   nStripsNoStack.push_back(94208);
0863   nStripsNoStack.push_back(516096);
0864   nStripsNoStack.push_back(589824);
0865   nStripsNoStack.push_back(331776);
0866   nStripsNoStack.push_back(368640);
0867   nStripsNoStack.push_back(608256);
0868   nStripsNoStack.push_back(681984);
0869   nStripsNoStack.push_back(253952);
0870   nStripsNoStack.push_back(253952);
0871   nStripsNoStack.push_back(253952);
0872   nStripsNoStack.push_back(217088);
0873   nStripsNoStack.push_back(217088);
0874   nStripsNoStack.push_back(217088);
0875   nStripsNoStack.push_back(180224);
0876   nStripsNoStack.push_back(180224);
0877   nStripsNoStack.push_back(159744);
0878   nStripsNoStack.push_back(253952);
0879   nStripsNoStack.push_back(253952);
0880   nStripsNoStack.push_back(253952);
0881   nStripsNoStack.push_back(217088);
0882   nStripsNoStack.push_back(217088);
0883   nStripsNoStack.push_back(217088);
0884   nStripsNoStack.push_back(180224);
0885   nStripsNoStack.push_back(180224);
0886   nStripsNoStack.push_back(159744);
0887 
0888   for (unsigned int i = 0; i < tkParts.size(); i++) {
0889     modulesStackNormFactors[tkParts[i]] = nModulesStack[i];
0890     modulesNoStackNormFactors[tkParts[i]] = nModulesNoStack[i];
0891     fibersStackNormFactors[tkParts[i]] = nFibersStack[i];
0892     fibersNoStackNormFactors[tkParts[i]] = nFibersNoStack[i];
0893     APVsStackNormFactors[tkParts[i]] = nAPVsStack[i];
0894     APVsNoStackNormFactors[tkParts[i]] = nAPVsNoStack[i];
0895     stripsStackNormFactors[tkParts[i]] = nStripsStack[i];
0896     stripsNoStackNormFactors[tkParts[i]] = nStripsNoStack[i];
0897   }
0898 
0899   //   for(std::map< std::string, unsigned int>::iterator it = allStripsTK.begin(); it != allStripsTK.end(); it++)
0900   //   {
0901   //     std::cout << it->first.c_str() << " " << it->second << std::endl;
0902   //   }
0903 }
0904 
0905 double findNormFactor(const std::string currentPlotType, const std::string currentPart, const bool stackOption) {
0906   //   std::cout << "findNormFactor(): Finding normalization factor for this part: \"" << currentPart.c_str() << "\".\n";
0907   //   std::cout << "                  Plot type is: \"" << currentPlotType.c_str() << "\".\n";
0908   //   std::cout << "                  stack option is: " << stackOption << std::endl;
0909   if (stackOption) {
0910     if (currentPlotType == "BadModules") {
0911       return modulesStackNormFactors[currentPart];
0912     } else if (currentPlotType == "BadFibers") {
0913       return fibersStackNormFactors[currentPart];
0914     } else if (currentPlotType == "BadAPVs") {
0915       return APVsStackNormFactors[currentPart];
0916     } else if (currentPlotType == "AllBadStrips" || currentPlotType == "BadStripsFromAPVs" ||
0917                currentPlotType == "BadStrips") {
0918       return stripsStackNormFactors[currentPart];
0919     } else {
0920       std::cout << "findNormFactor(): ERROR! Requested a non supported plot type: " << currentPlotType.c_str()
0921                 << std::endl;
0922       std::cout << "                  Add this to the function body or correct the error\n";
0923       return 0;  // This should trigger a divByZero error...
0924     }
0925   } else {
0926     if (currentPlotType == "BadModules") {
0927       return modulesNoStackNormFactors[currentPart];
0928     } else if (currentPlotType == "BadFibers") {
0929       return fibersNoStackNormFactors[currentPart];
0930     } else if (currentPlotType == "BadAPVs") {
0931       return APVsNoStackNormFactors[currentPart];
0932     } else if (currentPlotType == "BadStrips" || currentPlotType == "BadStripsFromAPVs" ||
0933                currentPlotType == "AllBadStrips") {
0934       return stripsNoStackNormFactors[currentPart];
0935     } else {
0936       std::cout << "findNormFactor(): ERROR! Requested a non supported plot type: \"" << currentPlotType.c_str()
0937                 << "\"\n";
0938       std::cout << "                  Add this to the function body or correct the error otherwise.\n";
0939       return 0;  // This should trigger a divByZero error...
0940     }
0941   }
0942 }