Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:24:29

0001 #include "TFile.h"
0002 #include "TTree.h"
0003 #include "TH1F.h"
0004 #include "TCanvas.h"
0005 #include "TApplication.h"
0006 #include "TRefArray.h"
0007 #include "TStyle.h"
0008 #include "TGraph.h"
0009 #include <string>
0010 
0011 int         colrs[5]  = {1, 2, 6, 3, 7};
0012 int         symbol[5] = {24, 29, 25, 27, 26};
0013 double      cutv[5]   = {10., 1., 0.1, 0.01, 0.001};
0014 
0015 void PlotEcal(int ien=10, int birk=0, int cuts=1, int logy=0, int save=0) {
0016 
0017   setStyle();
0018   char fileName[50];
0019   TH1F *hi[5];
0020   int   nh=5;
0021 
0022   std::string flag="";
0023   TFile *file;
0024   if (birk != 2) {
0025     cuts = 0;
0026     if (birk == 0) flag = "No";
0027     for (int i=0; i<5; i++) {
0028       int cut=i+1;
0029       sprintf (fileName,"el%dGeVCut%d%sBirk_2.root",ien,cut,flag.c_str());
0030       std::cout << fileName << "\n";
0031       file = new TFile (fileName);
0032       hi[i] = (TH1F*) file->Get("DQMData/HcalTB06Histo/edecS");
0033     }
0034   } else {
0035     sprintf (fileName,"el%dGeVCut%dBirk_2.root",ien,cuts);
0036     std::cout << fileName << "\n";
0037     file = new TFile (fileName);
0038     hi[1] = (TH1F*) file->Get("DQMData/HcalTB06Histo/edecS");
0039     sprintf (fileName,"el%dGeVCut%dNoBirk_2.root",ien,cuts);
0040     std::cout << fileName << "\n";
0041     file = new TFile (fileName);
0042     hi[0] = (TH1F*) file->Get("DQMData/HcalTB06Histo/edecS");
0043     nh = 2;
0044   }
0045 
0046   double xmin = 0.8*ien;
0047   double xmax = 1.05*ien;
0048   double ymax = 0, ymin = 100;
0049   for (int i=0; i<nh; i++) {
0050     hi[i]->Rebin(2); hi[i]->SetTitle("");
0051     int nx = hi[i]->GetNbinsX();
0052     for (int k=1; k <= nx; k++) {
0053       double xx = hi[i]->GetBinCenter(k);
0054       double yy = hi[i]->GetBinContent(k);
0055       if (xx > xmin && xx < xmax) {
0056         if (yy > ymax) ymax = yy;
0057         if (yy < ymin && yy > 0) ymin = yy;
0058       }
0059     }
0060   }
0061   if (logy == 0) {ymax *= 1.2; ymin *= 0.8;}
0062   else           {ymax *=10.0; ymin *= 0.2; }
0063 
0064   TCanvas *myc = new TCanvas("myc","",800,600); myc->SetLeftMargin(0.15);
0065   if (logy != 0) gPad->SetLogy(1);
0066 
0067   for (int i=0; i<nh; i++) {
0068     hi[i]->SetLineStyle(1);  hi[i]->SetLineWidth(1); 
0069     hi[i]->SetLineColor(colrs[i]); hi[i]->GetXaxis()->SetRangeUser(xmin,xmax);
0070     hi[i]->GetYaxis()->SetRangeUser(ymin,ymax);
0071     hi[i]->GetXaxis()->SetTitle("Energy Depost (GeV)");
0072     hi[i]->GetYaxis()->SetTitle("Events");
0073     if (i == 0) hi[i]->Draw();
0074     else        hi[i]->Draw("same");
0075   }
0076 
0077   char title[60], list[60];
0078   TLegend *leg1;
0079   if (birk == 2) {
0080     leg1 = new TLegend(0.15,0.80,0.62,0.90);
0081     sprintf(title, "Effect of Birks Law for electron with cut %6.3f mm", cutv[cuts-1]);
0082   } else {
0083     leg1 = new TLegend(0.15,0.70,0.55,0.90);
0084     if (birk == 0) sprintf(title, "Birks Law switched Off (electron)");
0085     else           sprintf(title, "Birks Law switched On (electron)");
0086   }
0087   leg1->SetHeader(title); leg1->SetFillColor(0);  leg1->SetTextSize(0.03);
0088   if (birk == 2) {
0089     sprintf(list, "Birks Law switched Off"); leg1->AddEntry(hi[0],list,"F");
0090     sprintf(list, "Birks Law switched On");  leg1->AddEntry(hi[1],list,"F");
0091   } else {
0092     for (int i=0; i<nh; i++) {
0093       sprintf (list, "Production cut set to %6.3f mm", cutv[i]);
0094       leg1->AddEntry(hi[i],list,"F");
0095     }
0096   }
0097   leg1->Draw();
0098 
0099   if (save != 0) {
0100     char fname[100];
0101     if (save > 0) sprintf (fname, "plot%dGeVBirk%dCut%d.eps", ien, birk, cuts);
0102     else          sprintf (fname, "plot%dGeVBirk%dCut%d.gif", ien, birk, cuts);
0103     myc->SaveAs(fname);
0104   }
0105 }
0106 
0107 void PlotMean(int ien=10, int birk=0, int save=0) {
0108 
0109   setStyle();
0110   char fileName[50];
0111   TH1F *hi;
0112   int   nh=5, nen=1, iene[5], nset=1, ienb=0;
0113   double mean[5][5];
0114 
0115   std::string flag="";
0116   TFile *file;
0117   if (ien > 0) {
0118     iene[0] = ien;
0119   } else {
0120     nen = 5;
0121     iene[0] = 1;
0122     iene[1] = 3;
0123     iene[2] = 10;
0124     iene[3] = 30;
0125     iene[4] = 100;
0126     if (ien == -2) ienb = 1;
0127   }
0128 
0129   if (birk != 2) {
0130     if (birk == 0) flag = "No";
0131     for (int k=0; k<nen; k++) {
0132       ien = iene[k];
0133       for (int i=0; i<nh; i++) {
0134     int cut=i+1;
0135     sprintf (fileName,"el%dGeVCut%d%sBirk_2.root",ien,cut,flag.c_str());
0136     std::cout << fileName << "\n";
0137     file = new TFile (fileName);
0138     hi = (TH1F*) file->Get("DQMData/HcalTB06Histo/edecS");
0139     mean[k][i]  = hi->GetMean(1);
0140     mean[k][i] /= (double)(ien);
0141       }
0142     }
0143     nset = nen;
0144   } else {
0145     ien = iene[0];
0146     for (int i=0; i<nh; i++) {
0147       int cut=i+1;
0148       sprintf (fileName,"el%dGeVCut%dBirk_2.root",ien,cut);
0149       std::cout << fileName << "\n";
0150       file = new TFile (fileName);
0151       hi = (TH1F*) file->Get("DQMData/HcalTB06Histo/edecS");
0152       mean[0][i]  = hi->GetMean(1);
0153       mean[0][i] /= (double)(ien);
0154       sprintf (fileName,"el%dGeVCut%dNoBirk_2.root",ien,cut);
0155       std::cout << fileName << "\n";
0156       file = new TFile (fileName);
0157       hi = (TH1F*) file->Get("DQMData/HcalTB06Histo/edecS");
0158       mean[1][i]  = hi->GetMean(1);
0159       mean[1][i] /= (double)(ien);
0160     }
0161     nset = 2;
0162   }
0163 
0164   TGraph *gr[5];
0165   double val[5], ymax=0.95, ymin=0.85;
0166   if (birk == 0 || birk == 2) ymax = 1.05;
0167   if (birk == 0)              ymin = 0.90;
0168   std::cout << nset << " set(s) of " << nh << " means\n";
0169   for (int k=0; k<nset; k++) {
0170     std::cout << "Set " << k;
0171     for (int i=0; i<nh; i++) {
0172       std::cout << " Mean[" << i << "] = " << mean[k][i];
0173       val[i] = mean[k][i];
0174     }
0175     std::cout << "\n";
0176     gr[k] = new TGraph(nh, cutv, val); gr[k]->SetMarkerSize(1.5);
0177     gr[k]->SetTitle(""); gr[k]->SetLineColor(colrs[k]);
0178     gr[k]->SetLineStyle(k+1); gr[k]->SetLineWidth(2);
0179     gr[k]->SetMarkerColor(colrs[k]);  gr[k]->SetMarkerStyle(symbol[k]);
0180     gr[k]->GetXaxis()->SetTitle("Production cut value (mm)");
0181     gr[k]->GetYaxis()->SetTitle("Mean Energy/Incident Energy");
0182     gr[k]->GetXaxis()->SetRangeUser(0.0005, 20.0);
0183     gr[k]->GetYaxis()->SetRangeUser(ymin, ymax);
0184   }
0185 
0186   TCanvas *myc = new TCanvas("myc","",800,600); myc->SetLeftMargin(0.15);
0187   gPad->SetLogx(1);
0188   gr[ienb]->GetXaxis()->SetRangeUser(0.0005, 20.0);
0189   gr[ienb]->Draw("alp");
0190   for (int k=ienb+1; k<nset; k++) 
0191     gr[k]->Draw("lp");
0192 
0193   TLegend *leg1;
0194   if (birk == 2) leg1 = new TLegend(0.15,0.80,0.50,0.90);
0195   else           leg1 = new TLegend(0.15,0.75,0.37,0.90);
0196   leg1->SetFillColor(0);  leg1->SetTextSize(0.03); 
0197   char list[60];
0198   if (birk == 2) {
0199     ien = iene[0];
0200     sprintf(list, "Birks Law On  at %d GeV",ien); leg1->AddEntry(gr[0],list,"P");
0201     sprintf(list, "Birks Law Off at %d GeV",ien); leg1->AddEntry(gr[1],list,"P");
0202   } else {
0203     if (birk == 0) flag = "Off";
0204     else           flag = "On";
0205     sprintf(list, "Birks Law %s", flag.c_str()); leg1->SetHeader(list);
0206     for (int k=ienb; k<nset; k++) {
0207       sprintf (list, "%d GeV electron",iene[k]); leg1->AddEntry(gr[k],list,"P");
0208     }
0209   }
0210   leg1->Draw();
0211 
0212   if (save != 0) {
0213     char fname[100];
0214     ien = iene[ienb];
0215     if (save > 0) {
0216       if (nen == 1) sprintf (fname, "plotRatio%dGeVBirk%d.eps", ien, birk);
0217       else          sprintf (fname, "plotRatio%dCut%d.eps",     ienb,birk);
0218     } else {
0219       if (nen == 1) sprintf (fname, "plotRatio%dGeVBirk%d.gif", ien, birk);
0220       else          sprintf (fname, "plotRatio%dCut%d.gif",     ienb,birk);
0221     }
0222     myc->SaveAs(fname);
0223   }
0224 }
0225 
0226 void setStyle() {
0227   
0228   gStyle->SetCanvasBorderMode(0); gStyle->SetCanvasColor(kWhite);
0229   gStyle->SetPadColor(kWhite);    gStyle->SetFrameBorderMode(0);
0230   gStyle->SetFrameBorderSize(1);  gStyle->SetFrameFillColor(0);
0231   gStyle->SetFrameFillStyle(0);   gStyle->SetFrameLineColor(1);
0232   gStyle->SetFrameLineStyle(1);   gStyle->SetFrameLineWidth(1);
0233   gStyle->SetTitleOffset(1.6,"Y");  gStyle->SetOptStat(0);
0234   gStyle->SetLegendBorderSize(1);
0235 
0236 }