test/macros/RelValMacro.C

0001 #include "TFile.h"
0002 #include "TString.h"
0003 #include "TCanvas.h"
0004 #include "TH1.h"
0005 #include "TH2.h"
0006 #include "TProfile.h"
0007 #include "TDirectory.h"
0008 #include "TPaveText.h"
0009 #include "TPaveStats.h"
0010 #include "TROOT.h"
0011 #include "TStyle.h"
0012 #include "TColor.h"
0013 #include "TLegend.h"
0014 #include "TKey.h"
0015 #include "TClass.h"
0016
0017 #include <iostream>
0018 #include <map>
0019 #include <string>
0020 #include <cstring>
0021 #include <cstdio>
0022
0023 #include <fstream>
0024
0025 #include "rootlogon.h"
0026
0027 #include <Python.h>
0028 #include <boost/python.hpp>
0029 #include <vector>
0030
0031 template <class T1, class T2>
0032 void prn(T1 s1, T2 s2) {
0033   std::cout << "\t>> " << s1 << ": " << s2 << std::endl;
0034 }
0035
0036 void RelValMacro(std::string seriesOfTubes);
0037 void ProcessRelVal(TFile* ref_file,
0038                    TFile* val_file,
0039                    std::string ref_vers,
0040                    std::string val_vers,
0041                    std::string histName,
0042                    std::string outLabel,
0043                    int nRebin,
0044                    double xAxisMin,
0045                    double xAxisMax,
0046                    double yAxisMin,
0047                    double yAxisMax,
0048                    std::string dimSwitch,
0049                    std::string statSwitch,
0050                    std::string chi2Switch,
0051                    std::string logSwitch,
0052                    std::string ratioFlag,
0053                    int refCol,
0054                    int valCol,
0055                    std::string xAxisTitle,
0056                    std::string normFlag,
0057                    std::string histName2 = "");
0058 template <class T>
0059 void setObjProps(T obj);
0060
0061 class DirectoryFinder {
0062 private:
0063   std::map<std::string, TDirectory*> ptdMap;
0064   TDirectory* findDirectory(TDirectory* target, std::string& s, int dig = 2);
0065
0066 public:
0067   TDirectory* operator()(TDirectory* target, std::string& s);
0068 } dfRef, dfVal;
0069
0070 void RelValMacro(std::string seriesOfTubes) {
0071   //Split the string passed from the python3 driver
0072   std::stringstream ss(seriesOfTubes);
0073   std::string item;
0074   std::vector<std::string> props;
0075   while (getline(ss, item, '|')) {
0076     props.push_back(item);
0077   }
0078   std::string ref_vers = props[0];
0079   std::string val_vers = props[1];
0080   std::string rfname = props[2];
0081   std::string vfname = props[3];
0082   std::string histName = props[4];
0083   std::string ofileName = props[5];
0084   int nRebin = std::stoi(props[6]);
0085   double xAxisMin = std::stod(props[7]);
0086   double xAxisMax = std::stod(props[8]);
0087   double yAxisMin = std::stod(props[9]);
0088   double yAxisMax = std::stod(props[10]);
0089   std::string dimFlag = props[11];
0090   std::string statFlag = props[12];
0091   std::string chi2Flag = props[13];
0092   std::string logFlag = props[14];
0093   std::string ratioFlag = props[15];
0094   int refCol = std::stoi(props[16]);
0095   int valCol = std::stoi(props[17]);
0096   std::string xAxisTitle = props[18];
0097   std::string histName2 = props[19];
0098   std::string normFlag = props[20];
0099
0100   if (strcmp(histName.c_str(), "HcalDigiTask/HcalDigiTask_signal_amplitude_HE") == 0) {
0101     std::cout << "==================" << std::endl;
0102     std::cout << xAxisMin << "  " << xAxisMax << std::endl;
0103   }
0104   //Warning!!! This rootlogon hacks the root color pallate.  This should probably be rewritten.
0105   setColors();
0106
0107   TFile* Ref_File = new TFile(rfname.c_str());
0108   TFile* Val_File = new TFile(vfname.c_str());
0109
0110   if (Ref_File && Val_File) {
0111     if (histName2 == "none")
0112       histName2 = "";
0113
0114     //Make plot
0115     ProcessRelVal(Ref_File,
0116                   Val_File,
0117                   ref_vers,
0118                   val_vers,
0119                   histName,
0120                   ofileName,
0121                   nRebin,
0122                   xAxisMin,
0123                   xAxisMax,
0124                   yAxisMin,
0125                   yAxisMax,
0126                   dimFlag,
0127                   statFlag,
0128                   chi2Flag,
0129                   logFlag,
0130                   ratioFlag,
0131                   refCol,
0132                   valCol,
0133                   xAxisTitle,
0134                   histName2,
0135                   normFlag);
0136   } else {
0137     if (!Ref_File)
0138       std::cout << "Input root file \"" << rfname << "\" not found!!!" << std::endl;
0139     if (!Val_File)
0140       std::cout << "Input root file \"" << vfname << "\" not found!!!" << std::endl;
0141   }
0142
0143   //    ProcessSubDetCT(Ref_File, Val_File, RelValStream, CT_nHist1, CT_nHist2, CT_nProf, CT_nHistTot, ref_vers, val_vers, harvest);
0144
0145   return;
0146 }
0147
0148 void ProcessRelVal(TFile* ref_file,
0149                    TFile* val_file,
0150                    std::string ref_vers,
0151                    std::string val_vers,
0152                    std::string histName,
0153                    std::string outLabel,
0154                    int nRebin,
0155                    double xAxisMin,
0156                    double xAxisMax,
0157                    double yAxisMin,
0158                    double yAxisMax,
0159                    std::string dimSwitch,
0160                    std::string statSwitch,
0161                    std::string chi2Switch,
0162                    std::string logSwitch,
0163                    std::string ratioFlag,
0164                    int refCol,
0165                    int valCol,
0166                    std::string xAxisTitle,
0167                    std::string histName2,
0168                    std::string normFlag) {
0169   std::string NormHist = "HcalRecHitTask/N_HB";
0170
0171   //split directory off histName
0172   int slashLoc = histName.rfind("/");
0173   std::string histDir = histName.substr(0, slashLoc);
0174   if (slashLoc < histName.size() - 1)
0175     histName = histName.substr(slashLoc + 1, histName.size());
0176
0177   int slashLocN = NormHist.rfind("/");
0178   std::string histDirN = NormHist.substr(0, slashLocN);
0179   if (slashLocN < NormHist.size() - 1)
0180     NormHist = NormHist.substr(slashLocN + 1, NormHist.size());
0181
0182   std::cout << "Processing \"" << histDir << "/" << histName << "\"" << std::endl;
0183
0184   //Get objects from TFiles
0185   TDirectory* refTD = dfRef(ref_file, histDir);
0186   TObject* refObj = 0;
0187
0188   TDirectory* refTDN = dfRef(ref_file, histDirN);
0189   TObject* refObjN = 0;
0190
0191   if (refTD) {
0192     refObj = refTD->Get(histName.c_str());
0193     if (refObj)
0194       refObj = refObj->Clone();
0195   } else {
0196     std::cout << "Cannot find directory \"" << histDir << "\" in file \"" << ref_file->GetName() << "\"" << std::endl;
0197     return;
0198   }
0199   if (!refObj) {
0200     std::cout << "Cannot find histogram \"" << histDir << "/" << histName << "\" in file \"" << ref_file->GetName()
0201               << "\"" << std::endl;
0202     return;
0203   }
0204
0205   if (refTDN) {
0206     refObjN = refTDN->Get(NormHist.c_str());
0207     if (refObjN)
0208       refObjN = refObjN->Clone();
0209   } else {
0210     std::cout << "Cannot find directory \"" << histDirN << "\" in file \"" << ref_file->GetName() << "\"" << std::endl;
0211   }
0212   if (!refObjN) {
0213     std::cout << "Cannot find histogram \"" << histDirN << "/" << NormHist << "\" in file \"" << ref_file->GetName()
0214               << "\"" << std::endl;
0215   }
0216
0217   TDirectory* valTD = dfVal(val_file, histDir);
0218   TObject* valObj = 0;
0219   TDirectory* valTDN = dfVal(val_file, histDirN);
0220   TObject* valObjN = 0;
0221   if (valTD) {
0222     valObj = valTD->Get(histName.c_str());
0223     if (valObj)
0224       valObj = valObj->Clone();
0225   } else {
0226     std::cout << "Cannot find directory \"" << histDir << "\" in file \"" << val_file->GetName() << "\"" << std::endl;
0227     return;
0228   }
0229   if (!valObj) {
0230     std::cout << "Cannot find histogram \"" << histDir << "/" << histName << "\" in file \"" << val_file->GetName()
0231               << "\"" << std::endl;
0232     return;
0233   }
0234
0235   if (valTDN) {
0236     valObjN = valTDN->Get(NormHist.c_str());
0237     if (valObjN)
0238       valObjN = valObjN->Clone();
0239   } else {
0240     std::cout << "Cannot find directory \"" << histDirN << "\" in file \"" << val_file->GetName() << "\"" << std::endl;
0241   }
0242   if (!valObjN) {
0243     std::cout << "Cannot find histogram \"" << histDirN << "/" << NormHist << "\" in file \"" << val_file->GetName()
0244               << "\"" << std::endl;
0245   }
0246
0247   //Try to continue processing even if N_HB is missing
0248   //We only care if the ratio flag is set
0249   //If we can't find any way to normalize the plots, unset the ratioflag
0250   if (std::stoi(ratioFlag) == 1) {
0251     if (!refTDN && !valTDN) {
0252       std::cout << "Cannot find directory \"" << histDirN << "\" in either file \"" << std::endl;
0253       ratioFlag = "0";
0254     }
0255
0256     if (!refObjN && !valObjN) {
0257       std::cout << "Cannot find histogram \"" << histDirN << "/" << NormHist << "\" in either file \"" << std::endl;
0258       ratioFlag = "0";
0259     } else if (!valObjN) {
0260       valObjN = refObjN->Clone();
0261       std::cout << "Using histogram \"" << NormHist << "from file \"" << ref_file->GetName() << std::endl;
0262     } else if (!refObjN) {
0263       refObjN = valObjN->Clone();
0264       std::cout << "Using histogram \"" << NormHist << "from file \"" << val_file->GetName() << std::endl;
0265     }
0266   }  // Make sure we can normalize ratio plots
0267
0268   std::cout << "Loaded \"" << histDir << "/" << histName << "\"" << std::endl;
0269
0270   //Format canvas
0271   TCanvas* myc = 0;
0272   if (dimSwitch.compare("PRwide") == 0) {
0273     gStyle->SetPadLeftMargin(0.06);
0274     gStyle->SetPadRightMargin(0.03);
0275     myc = new TCanvas("myc", "", 1200, 600);
0276   } else
0277     myc = new TCanvas("myc", "", 800, 600);
0278   //    gStyle->SetOptStat(0);
0279   myc->SetGrid();
0280
0281   TPad *pad1, *pad2;
0282
0283   // Ratio Flag
0284
0285   float nRef = 1, nVal = 1;
0286
0287   std::cout << "Ratio Flag: " << std::stoi(ratioFlag) << std::endl;
0288
0289   if (std::stoi(ratioFlag) == 1) {
0290     std::cout << "Histogram will include ratio" << std::endl;
0291
0292     TH1* refN_HB = (TH1*)refObjN;
0293     TH1* valN_HB = (TH1*)valObjN;
0294
0295     nRef = refN_HB->Integral();
0296     nVal = valN_HB->Integral();
0297
0298     // Divide canvas into two pads
0299     //    myc->Divide(1,2,0,0);
0300     pad1 = new TPad("pad1", "pad1", 0.0, 0.3, 1.0, 1.0, 0);
0301     pad1->SetBottomMargin(1);                                 // Upper and lower plots are joined (0) or separate (1)
0302     pad1->SetGridx();                                         // Vertical grid
0303     pad1->SetFillColor(kCyan - 10);                           //spandey
0304     pad2 = new TPad("pad2", "pad2", 0.0, 0.03, 1.0, 0.3, 0);  //spandey updated pad size
0305     pad2->SetTopMargin(0);
0306     pad2->SetBottomMargin(0.2);
0307     pad2->SetGridx();                // vertical grid
0308     pad2->SetGridy();                // horizontal grid
0309     pad2->SetFillColor(kCyan - 10);  //spandey
0310
0311     pad1->Draw();
0312     pad2->Draw();
0313
0314     //    float pad2width = pad2->GetWw();
0315     //    float pad2height = pad2->GetWh() * pad2->GetAbsHNCD();
0316     //    float x2pixels = 10;
0317     //    float y2pixels = 10;
0318     //    float x2size = x2pixels / pad2width;
0319     //    float y2size = y2pixels / pad2height;
0320
0321     //Format pads
0322     //    myc->cd(1);
0323     //    pad1->cd();
0324     if (logSwitch.compare("Log") == 0 && dimSwitch.compare("2D") == 0) {
0325       pad1->SetLogy(0);
0326       pad1->SetLogz(1);
0327     } else if (logSwitch.compare("Log") == 0) {
0328       pad1->SetLogy(1);
0329     }
0330     //        pad2->cd();
0331     pad2->SetGridy();
0332
0333     //        pad1->cd();
0334   }
0335
0336   std::string xTitleCheck = xAxisTitle;
0337   xTitleCheck = xTitleCheck.substr(1, 7);
0338
0339   if (dimSwitch.compare("1D") == 0) {
0340     //Get histograms from objects
0341     TH1* ref_hist1 = (TH1*)refObj;
0342     TH1* val_hist1 = (TH1*)valObj;
0343
0344     // change what is embedded in DQM histos
0345     setObjProps(ref_hist1);
0346     setObjProps(val_hist1);
0347
0348     //Rebin histograms -- has to be done first
0349     if (nRebin != 1) {
0350       ref_hist1->Rebin(nRebin);
0351       val_hist1->Rebin(nRebin);
0352     }
0353
0354     TH1* ratio_hist1;
0355
0356     // Ratio Flag
0357     if (std::stoi(ratioFlag) == 1) {
0358       //Let's normalize the val plot to have the same number of events as the ref plot
0359       //But only if normFlag isn't tripped
0360       if (normFlag.compare("Norm") == 0)
0361         val_hist1->Scale(nRef / nVal);
0362
0363       //Create Copies (Clones) to use in Ratio Plot
0364       TH1* ref_hist1_clone = (TH1*)ref_hist1->Clone("ref_hist1_clone");
0365       TH1* val_hist1_clone = (TH1*)val_hist1->Clone("val_hist1_clone");
0366
0367       //Prepare clones for correct uncertainties
0368       ref_hist1_clone->Sumw2();
0369       val_hist1_clone->Sumw2();
0370
0371       // Normalize (scale = n_ref/n_val)
0372       //float n_ref = ref_hist1_clone->Integral();
0373       //float n_val = val_hist1_clone->Integral();
0374       //float scale = n_ref/n_val;
0375       //val_hist1_clone->Scale(scale);
0376
0377       //Create ratio histogram (val - ref)/ref
0378       ratio_hist1 = (TH1*)val_hist1_clone;
0379       ratio_hist1->Sumw2();
0380       ratio_hist1->Add(ref_hist1_clone, -1.);
0381       ratio_hist1->Divide(ref_hist1_clone);
0382
0383       //            //Format Ratio Plot
0384       //            float pad2width = pad2->GetWw();
0385       //            float pad2height = pad2->GetWh() * pad2->GetAbsHNDC();
0386       //            float x2pixels = 100;
0387       //            float y2pixels = 15;
0388       //            float x2size = x2pixels / pad2width;
0389       //            float y2size = y2pixels / pad2height;
0390       //
0391       //            TAxis* x2axis = ratio_hist1->GetXaxis();
0392       //            TAxis* y2axis = ratio_hist1->GetYaxis();
0393       //
0394       //            x2axis->SetTitleOffset(2);
0395       //            x2axis->SetTitleSize(0.15);
0396       //            x2axis->SetLabelSize(x2size);
0397       //
0398       //            y2axis->SetTitleOffset(0.3);
0399       //            y2axis->SetTitleSize(0.12);
0400       //            y2axis->SetRangeUser(0,2.5);
0401       //            y2axis->SetLabelSize(y2size);
0402
0403       // Sanitizing axis inputs
0404       //Min/Max Convetion: Default AxisMin = 0. Default AxisMax = -1.
0405       //xAxis
0406       if (xAxisMin == 0)
0407         xAxisMin = ref_hist1->GetXaxis()->GetXmin();
0408       if (xAxisMax < 0)
0409         xAxisMax = ref_hist1->GetXaxis()->GetXmax();
0410
0411       //Sanitize xAxis inputs
0412       if (xAxisMin < ref_hist1->GetXaxis()->GetXmin())
0413         xAxisMin = ref_hist1->GetXaxis()->GetXmin();
0414       if (xAxisMax > ref_hist1->GetXaxis()->GetXmax())
0415         xAxisMax = ref_hist1->GetXaxis()->GetXmax();
0416
0417       ratio_hist1->SetTitle("");
0418       ratio_hist1->SetLineStyle(1);
0419       ratio_hist1->SetMarkerStyle(1);
0420       ratio_hist1->SetMarkerSize(0.02);
0421
0422       //Format Ratio Plot
0423       //lets get schwifty
0424       float pad2width = pad2->GetWw();
0425       float pad2height = pad2->GetWh() * pad2->GetAbsHNDC();
0426       float x2pixels = 100;
0427       float y2pixels = 15;
0428       float x2size = x2pixels / pad2width;
0429       float y2size = y2pixels / pad2height;
0430
0431       TAxis* x2axis = ratio_hist1->GetXaxis();
0432       TAxis* y2axis = ratio_hist1->GetYaxis();
0433
0434       x2axis->SetTitleOffset(1.0);  // Important for seeing x-axis title!
0435       x2axis->SetTitleSize(0.1);    //spandey
0436       x2axis->SetLabelSize(x2size * 0.64);
0437       x2axis->SetRangeUser(xAxisMin, xAxisMax);
0438
0439       y2axis->SetTitle("(val - ref)/ref");
0440       y2axis->SetTitleOffset(0.3);
0441       y2axis->SetTitleSize(0.12);
0442       //        y2axis->SetRangeUser(0,2.5);
0443       y2axis->SetLabelSize(y2size);
0444       y2axis->SetNdivisions(4);
0445
0446       ratio_hist1->SetStats(kFALSE);
0447     }
0448
0449     //Set the colors, styles, titles, stat boxes and format axes for the histograms
0450     ref_hist1->SetStats(kTRUE);
0451     val_hist1->SetStats(kTRUE);
0452
0453     if (statSwitch.compare("Stat") != 0 && statSwitch.compare("Statrv") != 0) {
0454       ref_hist1->SetStats(kFALSE);
0455       val_hist1->SetStats(kFALSE);
0456     }
0457
0458     //Min/Max Convetion: Default AxisMin = 0. Default AxisMax = -1.
0459     //xAxis
0460     if (xAxisMin == 0)
0461       xAxisMin = ref_hist1->GetXaxis()->GetXmin();
0462     if (xAxisMax < 0)
0463       xAxisMax = ref_hist1->GetXaxis()->GetXmax();
0464
0465     //Sanitize xAxis inputs
0466     if (xAxisMin < ref_hist1->GetXaxis()->GetXmin())
0467       xAxisMin = ref_hist1->GetXaxis()->GetXmin();
0468     if (xAxisMax > ref_hist1->GetXaxis()->GetXmax())
0469       xAxisMax = ref_hist1->GetXaxis()->GetXmax();
0470
0471     if (xAxisMax > 0 || xAxisMin != 0) {
0472       ref_hist1->GetXaxis()->SetRangeUser(xAxisMin, xAxisMax);
0473       val_hist1->GetXaxis()->SetRangeUser(xAxisMin, xAxisMax);
0474     }
0475     //yAxis
0476     if (yAxisMin != 0)
0477       ref_hist1->SetMinimum(yAxisMin);
0478     if (yAxisMax > 0)
0479       ref_hist1->SetMaximum(yAxisMax);
0480     else if (ref_hist1->GetMaximum() < val_hist1->GetMaximum() && val_hist1->GetMaximum() > 0) {
0481       if (logSwitch.compare("Log") == 0)
0482         ref_hist1->SetMaximum(2 * val_hist1->GetMaximum());
0483       else
0484         ref_hist1->SetMaximum(1.05 * val_hist1->GetMaximum());
0485     }
0486
0487     //Title
0488     //        if (xTitleCheck != "NoTitle") ref_hist1->GetXaxis()->SetTitle(xAxisTitle.c_str());
0489     ref_hist1->GetXaxis()->SetTitle("");
0490     if (xTitleCheck != "NoTitle" && std::stoi(ratioFlag) == 1)
0491       ratio_hist1->GetXaxis()->SetTitle(xAxisTitle.c_str());
0492     if (xTitleCheck != "NoTitle" && std::stoi(ratioFlag) != 1)
0493       ref_hist1->GetXaxis()->SetTitle(xAxisTitle.c_str());
0494
0495     //Different histo colors and styles
0496     ref_hist1->SetTitle("");
0497     ref_hist1->SetLineColor(refCol);
0498     ref_hist1->SetLineStyle(1);
0499     ref_hist1->SetMarkerSize(0.02);
0500
0501     val_hist1->SetTitle("");
0502     val_hist1->SetLineColor(valCol);
0503     val_hist1->SetLineStyle(2);
0504     val_hist1->SetMarkerSize(0.02);
0505     if (statSwitch.compare("Stat") != 0 && statSwitch.compare("Statrv") != 0) {
0506       ref_hist1->SetLineWidth(2);
0507       val_hist1->SetLineWidth(2);
0508     }
0509
0510     //Legend
0511     TLegend* leg = new TLegend(0.50, 0.91, 0.84, 0.99, "", "brNDC");
0512     leg->SetBorderSize(2);
0513     leg->SetFillStyle(1001);
0514     leg->AddEntry(ref_hist1, ("CMSSW_" + ref_vers).c_str(), "l");
0515     leg->AddEntry(val_hist1, ("CMSSW_" + val_vers).c_str(), "l");
0516
0517     //It's time to draw (#yolo)!
0518     if (chi2Switch.compare("Chi2") == 0) {
0519       // Title Time
0520
0521       //Draw and save histograms
0522       if (std::stoi(ratioFlag) == 1) {
0523         pad1->cd();
0524       }
0525       ref_hist1->SetFillColor(40);  //42 Originally, now 40 which is light brown
0526       ref_hist1->Draw("hist");
0527       val_hist1->SetLineStyle(1);
0528       if (statSwitch.compare("Statrv") == 0)
0529         val_hist1->Draw("sames e0");
0530       else
0531         val_hist1->Draw("same e0");
0532
0533       // Ratio Flag
0534       if (std::stoi(ratioFlag) == 1) {
0535         //Draw ratio
0536         pad2->cd();
0537         //pad1->cd();
0538         ratio_hist1->Draw();
0539         pad1->cd();
0540         //pad2->cd();
0541       }
0542
0543       //Get p-value from chi2 test
0544       const float NCHI2MIN = 0.01;
0545
0546       float pval;
0547       char tempbuff[30];
0548
0549       pval = ref_hist1->Chi2Test(val_hist1);
0550
0551       sprintf(tempbuff, "Chi2 p-value: %6.3E", pval);
0552
0553       TPaveText* ptchi2 = new TPaveText(0.05, 0.92, 0.35, 0.99, "NDC");
0554
0555       if (pval > NCHI2MIN)
0556         ptchi2->SetFillColor(kGreen);
0557       else
0558         ptchi2->SetFillColor(kRed);
0559
0560       ptchi2->SetTextSize(0.03);
0561       ptchi2->AddText(tempbuff);
0562       ptchi2->Draw();
0563     } else {
0564       // Title Time
0565
0566       //Draw and save histograms
0567       if (std::stoi(ratioFlag) == 1) {
0568         pad1->cd();
0569       }
0570       ref_hist1->Draw("hist");
0571       if (statSwitch.compare("Statrv") == 0)
0572         val_hist1->Draw("hist sames");
0573       else
0574         val_hist1->Draw("hist same");
0575
0576       // Ratio Flag
0577       if (std::stoi(ratioFlag) == 1) {
0578         //Draw ratio
0579         pad2->cd();
0580         ratio_hist1->Draw();
0581         pad1->cd();
0582       }
0583     }
0584
0585     //Stat Box where required
0586     if (statSwitch.compare("Stat") == 0 || statSwitch.compare("Statrv") == 0) {
0587       TPaveStats* ptstats_r = new TPaveStats(0.85, 0.86, 0.98, 0.98, "brNDC");
0588       ptstats_r->SetTextColor(refCol);
0589       ref_hist1->GetListOfFunctions()->Add(ptstats_r);
0590       ptstats_r->SetParent(ref_hist1->GetListOfFunctions());
0591       TPaveStats* ptstats_v = new TPaveStats(0.85, 0.74, 0.98, 0.86, "brNDC");
0592       ptstats_v->SetTextColor(valCol);
0593       val_hist1->GetListOfFunctions()->Add(ptstats_v);
0594       ptstats_v->SetParent(val_hist1->GetListOfFunctions());
0595
0596       ptstats_r->Draw();
0597       ptstats_v->Draw();
0598     }
0599
0600     leg->Draw();
0601
0602     myc->SaveAs(outLabel.c_str());
0603   }
0604   //Profiles not associated with histograms
0605   else if (dimSwitch.compare("PR") == 0 || dimSwitch.compare("PRwide") == 0) {
0606     //Get profiles from objects
0607     TProfile* ref_prof = (TProfile*)refObj;
0608     TProfile* val_prof = (TProfile*)valObj;
0609
0610     // HACK to change what is embedded in DQM histos
0611     setObjProps(ref_prof);
0612     setObjProps(val_prof);
0613
0614     //Legend
0615     TLegend* leg = new TLegend(0.50, 0.91, 0.84, 0.99, "", "brNDC");
0616     leg->SetBorderSize(2);
0617     leg->SetFillStyle(1001);
0618
0619     //Ordinary profiles
0620     if (dimSwitch.compare("PR") == 0) {
0621       ref_prof->SetTitle("");
0622       ref_prof->SetErrorOption("");
0623
0624       val_prof->SetTitle("");
0625       val_prof->SetErrorOption("");
0626
0627       ref_prof->GetXaxis()->SetTitle(xAxisTitle.c_str());
0628
0629       if (statSwitch.compare("Stat") != 0 && statSwitch.compare("Statrv") != 0) {
0630         ref_prof->SetStats(kFALSE);
0631         val_prof->SetStats(kFALSE);
0632       }
0633
0634       ref_prof->SetLineColor(41);
0635       ref_prof->SetLineStyle(1);
0636       ref_prof->SetLineWidth(1);
0637       ref_prof->SetMarkerColor(41);
0638       ref_prof->SetMarkerStyle(21);
0639       ref_prof->SetMarkerSize(0.8);
0640
0641       val_prof->SetLineColor(43);
0642       val_prof->SetLineStyle(1);
0643       val_prof->SetLineWidth(1);
0644       val_prof->SetMarkerColor(43);
0645       val_prof->SetMarkerStyle(22);
0646       val_prof->SetMarkerSize(1.0);
0647
0648       if (ref_prof->GetMaximum() < val_prof->GetMaximum() && val_prof->GetMaximum() > 0) {
0649         if (logSwitch.compare("Log") == 0)
0650           ref_prof->SetMaximum(2 * val_prof->GetMaximum());
0651         else
0652           ref_prof->SetMaximum(1.05 * val_prof->GetMaximum());
0653       }
0654
0655       ref_prof->Draw("hist pl");
0656       val_prof->Draw("hist pl same");
0657
0658       leg->AddEntry(ref_prof, ("CMSSW_" + ref_vers).c_str(), "pl");
0659       leg->AddEntry(val_prof, ("CMSSW_" + val_vers).c_str(), "pl");
0660     }  //Wide profiles
0661     else if (dimSwitch.compare("PRwide") == 0) {
0662       char temp[128];
0663       sprintf(temp, "%s_px_v", ref_prof->GetName());
0664       TH1* ref_fp = ref_prof->ProjectionX();
0665       TH1* val_fp = val_prof->ProjectionX(temp);
0666
0667       ref_fp->SetTitle("");
0668       val_fp->SetTitle("");
0669
0670       ref_fp->GetXaxis()->SetTitle(xAxisTitle.c_str());
0671
0672       if (statSwitch.compare("Stat") != 0 && statSwitch.compare("Statrv") != 0) {
0673         ref_fp->SetStats(kFALSE);
0674         val_fp->SetStats(kFALSE);
0675       }
0676
0677       int nbins = ref_fp->GetNbinsX();
0678       for (int j = 1; j < nbins; j++) {
0679         ref_fp->SetBinError(j, 0.);
0680         val_fp->SetBinError(j, 0.);
0681       }
0682       ref_fp->SetLineWidth(0);
0683       ref_fp->SetLineColor(0);  // 5 yellow
0684       ref_fp->SetLineStyle(1);
0685       ref_fp->SetMarkerColor(2);
0686       ref_fp->SetMarkerStyle(20);
0687       ref_fp->SetMarkerSize(0.5);
0688
0689       val_fp->SetLineWidth(0);
0690       val_fp->SetLineColor(0);  // 45 blue
0691       val_fp->SetLineStyle(2);
0692       val_fp->SetMarkerColor(4);
0693       val_fp->SetMarkerStyle(22);
0694       val_fp->SetMarkerSize(0.5);
0695
0696       if (ref_fp->GetMaximum() < val_fp->GetMaximum() && val_fp->GetMaximum() > 0) {
0697         if (logSwitch.compare("Log") == 0)
0698           ref_fp->SetMaximum(2 * val_fp->GetMaximum());
0699         else
0700           ref_fp->SetMaximum(1.05 * val_fp->GetMaximum());
0701       }
0702
0703       ref_fp->Draw("p9");
0704       val_fp->Draw("p9same");
0705
0706       leg->AddEntry(ref_fp, ("CMSSW_" + ref_vers).c_str(), "lp");
0707       leg->AddEntry(val_fp, ("CMSSW_" + val_vers).c_str(), "lp");
0708     }
0709
0710     leg->Draw("");
0711
0712     myc->SaveAs(outLabel.c_str());
0713   }  //Timing Histograms (special: read two lines at once)
0714   else if (dimSwitch.compare("TM") == 0) {
0715     //split directory off histName
0716     int slashLoc2 = histName2.rfind("/");
0717     std::string histDir2 = histName2.substr(0, slashLoc2);
0718     if (slashLoc2 < histName2.size() - 1)
0719       histName2 = histName2.substr(slashLoc2 + 1, histName2.size());
0720
0721     //Get objects from TFiles
0722     TDirectory* refTD2 = dfRef(ref_file, histDir2);
0723     TObject* refObj2 = refTD->Get(histName2.c_str())->Clone();
0724     TDirectory* valTD2 = dfVal(val_file, histDir2);
0725     TObject* valObj2 = valTD->Get(histName2.c_str())->Clone();
0726
0727     TH2* ref_hist2 = (TH2*)refObj;
0728     TProfile* ref_prof = (TProfile*)refObj2;
0729
0730     ref_hist2->SetMarkerStyle(21);
0731     ref_prof->SetMarkerStyle(21);
0732     ref_hist2->SetMarkerSize(0.02);
0733     ref_prof->SetMarkerSize(0.02);
0734
0735     TH2* val_hist2 = (TH2F*)valObj;
0736     TProfile* val_prof = (TProfile*)valObj2;
0737
0738     val_hist2->SetMarkerStyle(21);
0739     val_prof->SetMarkerStyle(21);
0740     val_hist2->SetMarkerSize(0.02);
0741     val_prof->SetMarkerSize(0.02);
0742
0743     // HACK to change what is embedded in DQM histos
0744     setObjProps(ref_hist2);
0745     setObjProps(val_hist2);
0746
0747     //Min/Max Convention: Default AxisMin = 0. Default AxisMax = -1.
0748     //xAxis
0749     if (xAxisMin == 0)
0750       xAxisMin = ref_hist2->GetXaxis()->GetXmin();
0751     if (xAxisMax < 0)
0752       xAxisMax = ref_hist2->GetXaxis()->GetXmax();
0753
0754     if (xAxisMax > 0 || xAxisMin != 0) {
0755       ref_hist2->GetXaxis()->SetRangeUser(xAxisMin, xAxisMax);
0756       val_hist2->GetXaxis()->SetRangeUser(xAxisMin, xAxisMax);
0757     }
0758     //yAxis
0759     if (yAxisMin != 0)
0760       ref_hist2->SetMinimum(yAxisMin);
0761     if (yAxisMax > 0)
0762       ref_hist2->SetMaximum(yAxisMax);
0763     else if (ref_hist2->GetMaximum() < val_hist2->GetMaximum() && val_hist2->GetMaximum() > 0) {
0764       if (logSwitch == "Log")
0765         ref_hist2->SetMaximum(2 * val_hist2->GetMaximum());
0766       else
0767         ref_hist2->SetMaximum(1.05 * val_hist2->GetMaximum());
0768     }
0769
0770     //AF
0771     if (yAxisMax > 0 || yAxisMin != 0) {
0772       ref_hist2->GetYaxis()->SetRangeUser(yAxisMin, yAxisMax);
0773       val_hist2->GetYaxis()->SetRangeUser(yAxisMin, yAxisMax);
0774     }
0775
0776     //Legend
0777     TLegend* leg = new TLegend(0.50, 0.91, 0.84, 0.99, "", "brNDC");
0778     leg->SetBorderSize(2);
0779     leg->SetFillStyle(1001);
0780
0781     ref_hist2->GetXaxis()->SetTitle(xAxisTitle.c_str());
0782     ref_hist2->SetStats(kFALSE);
0783
0784     ref_hist2->SetTitle("");
0785     val_hist2->SetTitle("");
0786
0787     ref_hist2->SetMarkerColor(refCol);  // rose
0788     ref_hist2->Draw();
0789     ref_prof->SetLineColor(41);
0790     ref_prof->Draw("same");
0791
0792     val_hist2->SetMarkerColor(valCol);
0793     val_hist2->Draw("same");
0794     val_prof->SetLineColor(45);
0795     val_prof->Draw("same");
0796
0797     leg->AddEntry(ref_prof, ("CMSSW_" + ref_vers).c_str(), "pl");
0798     leg->AddEntry(val_prof, ("CMSSW_" + val_vers).c_str(), "pl");
0799
0800     leg->Draw("");
0801
0802     myc->SaveAs(outLabel.c_str());
0803
0804     if (refObj2)
0805       delete refObj2;
0806     if (valObj2)
0807       delete valObj2;
0808   } else if (dimSwitch.compare("2D") == 0) {
0809     myc->SetGrid(0, 0);
0810
0811     //Get histograms from objects
0812     TH2* ref_hist2D = (TH2*)refObj;
0813     TH2* val_hist2D = (TH2*)valObj;
0814
0815     ref_hist2D->SetStats(kFALSE);
0816     val_hist2D->SetStats(kFALSE);
0817
0818     // HACK to change what is embedded in DQM histos
0819     setObjProps(ref_hist2D);
0820     setObjProps(val_hist2D);
0821
0822     ref_hist2D->SetTitle("");
0823     val_hist2D->SetTitle("");
0824
0825     // special zoom on HB/HE depth1
0826     //if (n2D == 1) {
0827     //    ref_hist2D->GetXaxis()->SetRangeUser(-29., 28.);
0828     //    val_hist2D->GetXaxis()->SetRangeUser(-29., 28.);
0829     //}
0830
0831     //Min/Max Convetion: Default AxisMin = 0. Default AxisMax = -1.
0832     //xAxis
0833     if (xAxisMax > 0 || xAxisMin != 0) {
0834       ref_hist2D->GetXaxis()->SetRangeUser(xAxisMin, xAxisMax);
0835       val_hist2D->GetXaxis()->SetRangeUser(xAxisMin, xAxisMax);
0836     }
0837     //yAxis
0838     if (yAxisMax > 0 || yAxisMin != 0) {
0839       ref_hist2D->GetYaxis()->SetRangeUser(yAxisMin, yAxisMax);
0840       val_hist2D->GetYaxis()->SetRangeUser(yAxisMin, yAxisMax);
0841     }
0842     //Set bin minimum to 0
0843     ref_hist2D->SetMinimum(0.0);
0844     val_hist2D->SetMinimum(0.0);
0845
0846     TLegend* leg1 = new TLegend(0.50, 0.91, 0.84, 0.99, "", "brNDC");
0847     leg1->SetBorderSize(2);
0848     leg1->SetFillStyle(1001);
0849     leg1->AddEntry(ref_hist2D, ("CMSSW_" + ref_vers).c_str(), "l");
0850
0851     if (xTitleCheck != "NoTitle")
0852       ref_hist2D->GetXaxis()->SetTitle(xAxisTitle.c_str());
0853     ref_hist2D->Draw("colz");
0854     leg1->Draw();
0855     myc->SaveAs(("ref_" + outLabel).c_str());
0856
0857     TLegend* leg2 = new TLegend(0.50, 0.91, 0.84, 0.99, "", "brNDC");
0858     leg2->SetBorderSize(2);
0859     leg2->SetFillStyle(1001);
0860     leg2->AddEntry(val_hist2D, ("CMSSW_" + val_vers).c_str(), "l");
0861
0862     if (xTitleCheck != "NoTitle")
0863       val_hist2D->GetXaxis()->SetTitle(xAxisTitle.c_str());
0864     val_hist2D->Draw("colz");
0865     leg2->Draw();
0866     myc->SaveAs(("val_" + outLabel).c_str());
0867   }
0868
0869   if (myc)
0870     delete myc;
0871   if (refObj)
0872     delete refObj;
0873   if (valObj)
0874     delete valObj;
0875
0876   return;
0877 }
0878
0879 TDirectory* DirectoryFinder::operator()(TDirectory* target, std::string& s) {
0880   if (ptdMap.find(s) == ptdMap.end())
0881     return (ptdMap[s] = findDirectory(target, s));
0882   else
0883     return ptdMap[s];
0884 }
0885
0886 TDirectory* DirectoryFinder::findDirectory(TDirectory* target, std::string& s, int dig) {
0887   TDirectory* retval = 0;
0888
0889   // loop over all keys in this directory
0890   TIter nextkey(target->GetListOfKeys());
0891   TKey *key, *oldkey = 0;
0892   while ((key = (TKey*)nextkey())) {
0893     //std::cout << "Found " << key->ReadObj()->GetName() << std::endl;
0894
0895     //keep only the highest cycle number for each key
0896     if (oldkey && !strcmp(oldkey->GetName(), key->GetName()))
0897       continue;
0898
0899     // read object from file
0900     //  target->cd();
0901     TObject* obj = key->ReadObj();
0902     //  obj->Print();
0903
0904     if (obj->IsA()->InheritsFrom(TDirectory::Class())) {
0905       // it's a subdirectory
0906       //std::cout << "Found subdirectory " << obj->GetName() << std::endl;
0907
0908       if (strcmp(s.c_str(), obj->GetName()) == 0)
0909         return (TDirectory*)obj;
0910
0911       if ((retval = findDirectory((TDirectory*)obj, s, dig - 1)))
0912         break;
0913
0914     } else if (dig < 1) {
0915       break;
0916     }
0917   }
0918
0919   return retval;
0920 }
0921
0922 template <class T>
0923 void setObjProps(T obj) {
0924   obj->GetXaxis()->SetLabelSize(0.04);
0925   obj->GetYaxis()->SetLabelSize(0.04);
0926   obj->GetXaxis()->SetTitleSize(0.045);
0927
0928   obj->GetXaxis()->SetTickLength(-0.015);
0929   obj->GetYaxis()->SetTickLength(-0.015);
0930
0931   obj->GetXaxis()->SetLabelOffset(0.02);
0932   obj->GetYaxis()->SetLabelOffset(0.02);
0933
0934   obj->GetXaxis()->SetTitleOffset(1.3);
0935 }
0936
0937 BOOST_PYTHON_MODULE(RelValMacro) {
0938   using namespace boost::python;
0939   def("RelValMacro", RelValMacro);
0940 }