Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:29:59 

0001 //////////////////////////////////////////////////////////////////////////////
0002 //
0003 // Usage:
0004 // .L MakeZDCPlots.C+g
0005 //
0006 //   To make plot of various quantities which are created by ZDCSimHitStudy
0007 //   from upto three different settings
0008 //
0009 //     makeHitStudyPlots(file1, tag1, file2, tag2, file3, tag3, todomin,
0010 //                       todomax, gtitle, ratio, save, dirnm)
0011 //
0012 //   where (for makeHitStudyPlots)
0013 //     fileN    std::string   Name of the Nth ROOT file
0014 //                            The fist must exist; blank means beng ignored
0015 //     tag1     std::string   Tag for the Nth file
0016 //     todomin  int           Minimum type # of histograms to be plotted [0]
0017 //     todomax  int           Maximum type # of histograms to be plotted [5]
0018 //     gtitle   std::string   Overall Titile [""]
0019 //     ratio    bool          if the ratio to be plotted [true]
0020 //                            (with respect to the first file)
0021 //     save     bool          If the canvas is to be saved [false]
0022 //     dirnm    std::string   Name of the directory [zdcSimHitStudy]
0023 //
0024 //////////////////////////////////////////////////////////////////////////////
0025 
0026 #include <TCanvas.h>
0027 #include <TChain.h>
0028 #include <TFile.h>
0029 #include <TFitResult.h>
0030 #include <TFitResultPtr.h>
0031 #include <TH1D.h>
0032 #include <TH2D.h>
0033 #include <TLegend.h>
0034 #include <TPaveStats.h>
0035 #include <TPaveText.h>
0036 #include <TProfile.h>
0037 #include <TROOT.h>
0038 #include <TStyle.h>
0039 #include <vector>
0040 #include <string>
0041 #include <iomanip>
0042 #include <iostream>
0043 #include <fstream>
0044 
0045 void makeHitStudyPlots(std::string file1 = "Forward.root",
0046                        std::string tag1 = "New",
0047                        std::string file2 = "Standard.root",
0048                        std::string tag2 = "Old",
0049                        std::string file3 = "",
0050                        std::string tag3 = "",
0051                        std::string gtitle = "",
0052                        int todomin = 0,
0053                        int todomax = 5,
0054                        bool ratio = true,
0055                        bool save = false,
0056                        std::string dirnm = "zdcSimHitStudy") {
0057   const int plots = 6;
0058   std::string names[plots] = {"ETot", "ETotT", "Edep", "Hits", "Indx", "Time"};
0059   int logy[plots] = {1, 1, 1, 0, 1, 1};
0060   int rebin[plots] = {10, 10, 10, 2, 1, 10};
0061   int xmax[plots] = {5, 5, 3, 5, 2, 2};
0062   int colors[5] = {1, 2, 4, 6, 46};
0063   int marker[5] = {20, 21, 22, 23, 24};
0064   int styles[5] = {1, 2, 3, 4, 5};
0065   bool debug(false);
0066 
0067   gStyle->SetCanvasBorderMode(0);
0068   gStyle->SetCanvasColor(kWhite);
0069   gStyle->SetPadColor(kWhite);
0070   gStyle->SetFillColor(kWhite);
0071   if (ratio) {
0072     gStyle->SetOptStat(10);
0073     gStyle->SetOptFit(1);
0074   } else {
0075     gStyle->SetOptStat(1110);
0076   }
0077   TFile* file[3];
0078   int nfile(0);
0079   std::string tag(""), tags[3];
0080   if (file1 != "") {
0081     file[nfile] = new TFile(file1.c_str());
0082     if (file[nfile]) {
0083       tags[nfile] = tag1;
0084       ++nfile;
0085       tag += tag1;
0086     }
0087   }
0088   if (file2 != "") {
0089     file[nfile] = new TFile(file2.c_str());
0090     if (file[nfile]) {
0091       tags[nfile] = tag2;
0092       ++nfile;
0093       tag += tag2;
0094     }
0095   }
0096   if (file3 != "") {
0097     file[nfile] = new TFile(file3.c_str());
0098     if (file[nfile]) {
0099       tags[nfile] = tag3;
0100       ++nfile;
0101       tag += tag3;
0102     }
0103   }
0104   if ((todomin < 0) || (todomin >= plots))
0105     todomin = 0;
0106   if (todomax < todomin) {
0107     todomax = todomin;
0108   } else if (todomax >= plots) {
0109     todomax = plots - 1;
0110   }
0111   std::cout << "Use " << nfile << " files from " << file1 << "," << file2 << " and " << file3 << " and look for "
0112             << todomin << ":" << todomax << std::endl;
0113   for (int todo = todomin; todo <= todomax; ++todo) {
0114     double y1(0.90), dy(0.12);
0115     double y2 = y1 - dy * nfile - 0.01;
0116     double y3 = y1 - 0.04 * (nfile - 1);
0117     TLegend* leg = (ratio) ? (new TLegend(0.10, y3, 0.65, y1)) : (new TLegend(0.65, y2 - nfile * 0.04, 0.90, y2));
0118     leg->SetBorderSize(1);
0119     leg->SetFillColor(10);
0120     TH1D* hist0[nfile];
0121     for (int i1 = 0; i1 < nfile; ++i1) {
0122       TDirectory* dir = static_cast<TDirectory*>(file[i1]->FindObjectAny(dirnm.c_str()));
0123       hist0[i1] = static_cast<TH1D*>(dir->FindObjectAny(names[todo].c_str()));
0124       if (debug)
0125         std::cout << "Reads histogram " << hist0[i1]->GetName() << " for " << hist0[i1]->GetTitle() << " from "
0126                   << file[i1] << std::endl;
0127     }
0128     std::string title = hist0[0]->GetTitle();
0129     if (debug)
0130       std::cout << "Histogram title " << title << std::endl;
0131     int istart = (ratio) ? 1 : 0;
0132     char name[100], namec[100];
0133     std::vector<TH1D*> hist1;
0134     for (int i1 = istart; i1 < nfile; ++i1) {
0135       if (ratio) {
0136         int nbin = hist0[0]->GetNbinsX();
0137         int nbinR = nbin / rebin[todo];
0138         double xlow = hist0[0]->GetXaxis()->GetBinLowEdge(1);
0139         double xhigh = hist0[0]->GetXaxis()->GetBinLowEdge(nbin) + hist0[0]->GetXaxis()->GetBinWidth(nbin);
0140         sprintf(name, "%s%sVs%sRatio", names[todo].c_str(), tags[0].c_str(), tags[i1].c_str());
0141         TH1D* histr = new TH1D(name, hist0[0]->GetTitle(), nbinR, xlow, xhigh);
0142         histr->SetTitle(gtitle.c_str());
0143         histr->GetXaxis()->SetTitle(title.c_str());
0144         histr->GetYaxis()->SetTitle(name);
0145         histr->GetXaxis()->SetLabelOffset(0.005);
0146         histr->GetXaxis()->SetTitleOffset(1.00);
0147         histr->GetYaxis()->SetLabelOffset(0.005);
0148         histr->GetYaxis()->SetTitleOffset(1.20);
0149         histr->GetYaxis()->SetRangeUser(0.0, xmax[todo]);
0150         histr->SetLineColor(colors[i1 - 1]);
0151         histr->SetLineStyle(styles[i1 - 1]);
0152         histr->SetMarkerStyle(marker[i1 - 1]);
0153         histr->SetMarkerColor(colors[i1 - 1]);
0154         double sumNum(0), sumDen(0), maxDev(0), xlow1(-1), xhigh1(xlow);
0155         for (int j = 0; j < nbinR; ++j) {
0156           double tnum(0), tden(0), rnum(0), rden(0);
0157           for (int i = 0; i < rebin[todo]; ++i) {
0158             int ib = j * rebin[todo] + i + 1;
0159             tnum += hist0[0]->GetBinContent(ib);
0160             tden += hist0[1]->GetBinContent(ib);
0161             rnum += ((hist0[0]->GetBinError(ib)) * (hist0[0]->GetBinError(ib)));
0162             rden += ((hist0[1]->GetBinError(ib)) * (hist0[1]->GetBinError(ib)));
0163           }
0164           if (tden > 0 && tnum > 0) {
0165             if (xlow1 < 0)
0166               xlow1 = hist0[0]->GetXaxis()->GetBinLowEdge(j * rebin[todo] + 1);
0167             xhigh1 = hist0[0]->GetXaxis()->GetBinLowEdge((j + 1) * rebin[todo]) +
0168                      hist0[0]->GetXaxis()->GetBinWidth((j + 1) * rebin[todo]);
0169             double rat = tnum / tden;
0170             double err = rat * sqrt((rnum / (tnum * tnum)) + (rden / (tden * tden)));
0171             if (debug)
0172               std::cout << "Bin " << j << " Ratio " << rat << " +- " << err << std::endl;
0173             histr->SetBinContent(j + 1, rat);
0174             histr->SetBinError(j + 1, err);
0175             double temp1 = (rat > 1.0) ? 1.0 / rat : rat;
0176             double temp2 = (rat > 1.0) ? err / (rat * rat) : err;
0177             sumNum += (fabs(1 - temp1) / (temp2 * temp2));
0178             sumDen += (1.0 / (temp2 * temp2));
0179             if (fabs(1 - temp1) > maxDev)
0180               maxDev = fabs(1 - temp1);
0181           }
0182         }
0183         histr->Fit("pol0", "+QRWLS", "", xlow1, xhigh1);
0184         sumNum = (sumDen > 0) ? (sumNum / sumDen) : 0;
0185         sumDen = (sumDen > 0) ? 1.0 / sqrt(sumDen) : 0;
0186         sprintf(name, "%sVs%sRatio", tags[0].c_str(), tags[i1].c_str());
0187         if (sumNum == 0)
0188           sumDen = 0;
0189         sprintf(name, "%s vs %s (%6.3f +- %6.3f)", tags[0].c_str(), tags[i1].c_str(), sumNum, sumDen);
0190         hist1.push_back(histr);
0191         leg->AddEntry(histr, name, "lp");
0192       } else {
0193         hist1.push_back(hist0[i1]);
0194         hist1.back()->SetTitle(gtitle.c_str());
0195         hist1.back()->GetXaxis()->SetTitle(title.c_str());
0196         hist1.back()->GetYaxis()->SetTitle("Events");
0197         hist1.back()->GetXaxis()->SetLabelOffset(0.005);
0198         hist1.back()->GetXaxis()->SetTitleOffset(1.00);
0199         hist1.back()->GetYaxis()->SetLabelOffset(0.005);
0200         hist1.back()->GetYaxis()->SetTitleOffset(1.40);
0201         hist1.back()->SetLineColor(colors[i1]);
0202         hist1.back()->SetLineStyle(styles[i1]);
0203         leg->AddEntry(hist1.back(), tags[i1].c_str(), "lp");
0204       }
0205     }
0206     sprintf(namec, "c_%s", hist1[0]->GetName());
0207     if (debug)
0208       std::cout << namec << " Canvas made for " << hist1[0]->GetName() << " with " << hist1.size() << " plots"
0209                 << std::endl;
0210     TCanvas* pad = new TCanvas(namec, namec, 500, 500);
0211     pad->SetRightMargin(0.10);
0212     pad->SetTopMargin(0.10);
0213     if ((!ratio) && (logy[todo] > 0))
0214       pad->SetLogy(1);
0215     else
0216       pad->SetLogy(0);
0217     for (unsigned int i1 = 0; i1 < hist1.size(); ++i1) {
0218       if ((!ratio) && (rebin[todo] > 1))
0219         hist1[i1]->Rebin(rebin[todo]);
0220       if (i1 == 0) {
0221         hist1[i1]->Draw();
0222       } else {
0223         hist1[i1]->Draw("sames");
0224       }
0225       if (debug)
0226         std::cout << "Drawing histograms for " << hist1[i1]->GetName() << ":" << i1 << " in canvas " << pad->GetName()
0227                   << std::endl;
0228       pad->Update();
0229       double dy0 = (ratio) ? 0.75 * dy : dy;
0230       TPaveStats* st = ((TPaveStats*)hist1[i1]->GetListOfFunctions()->FindObject("stats"));
0231       if (st != NULL) {
0232         st->SetFillColor(kWhite);
0233         st->SetLineColor(colors[i1]);
0234         st->SetTextColor(colors[i1]);
0235         st->SetY1NDC(y1 - dy0);
0236         st->SetY2NDC(y1);
0237         st->SetX1NDC(0.65);
0238         st->SetX2NDC(0.90);
0239         y1 -= dy0;
0240       }
0241       pad->Modified();
0242       pad->Update();
0243       leg->Draw("same");
0244       pad->Update();
0245       if (save) {
0246         sprintf(name, "%s.pdf", pad->GetName());
0247         pad->Print(name);
0248       }
0249     }
0250   }
0251 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!