| 1 |
#include <memory>
|
| 2 |
#include <TROOT.h>
|
| 3 |
#include <TString.h>
|
| 4 |
#include <TStyle.h>
|
| 5 |
#include <TH1.h>
|
| 6 |
#include <TCanvas.h>
|
| 7 |
#include <TLegend.h>
|
| 8 |
#include <TFile.h>
|
| 9 |
#include <TTree.h>
|
| 10 |
#include <Riostream.h>
|
| 11 |
// #include "TLorentzVector.h"
|
| 12 |
// #include <map>
|
| 13 |
// #include <vector>
|
| 14 |
// #include <iostream>
|
| 15 |
using namespace std;
|
| 16 |
|
| 17 |
void comparePDF(){
|
| 18 |
gROOT->SetStyle("CMS");
|
| 19 |
// gStyle->SetOptStat("n");
|
| 20 |
// gStyle->SetOptTitle(1);
|
| 21 |
// gStyle->SetOptFit(11111);
|
| 22 |
gStyle->SetPadTopMargin(0.11);
|
| 23 |
gStyle->SetPadRightMargin(0.05);
|
| 24 |
gStyle->SetPadLeftMargin(0.11);
|
| 25 |
gStyle->SetPadBottomMargin(0.11);
|
| 26 |
|
| 27 |
const double wttbar = 0.0081;
|
| 28 |
const int nmembers = 45;
|
| 29 |
TH1D *h1 = new TH1D[nmembers];
|
| 30 |
TH1D *h2 = new TH1D[nmembers];
|
| 31 |
TH1D *hm3 = new TH1D[nmembers-1];
|
| 32 |
TH1D *hm3prime = new TH1D[nmembers-1];
|
| 33 |
double ratio[2] = {0.,0.};
|
| 34 |
double systMax1 = 0.,systMax2 = 0.;
|
| 35 |
|
| 36 |
TFile *f0=TFile::Open("TopAnalysis_TTJets-madgraph_Fall08_all_all_pdf.root");
|
| 37 |
|
| 38 |
TCanvas *cv = new TCanvas("cv","M3 & M3' systematics",1000,800);
|
| 39 |
cv->Divide(2,2,0.001,0.001);
|
| 40 |
cv->cd(1);
|
| 41 |
for(int i=0; i < nmembers; i++){
|
| 42 |
TString suffix="PDF_";
|
| 43 |
suffix+=i;
|
| 44 |
|
| 45 |
h1[i] = (TH1D*)f0->Get(TString("Mass/HadronicTop_mass_M3_")+suffix);
|
| 46 |
h1[i].SetLineColor(i+1);
|
| 47 |
h1[i].Scale(wttbar);
|
| 48 |
|
| 49 |
if(i==0){
|
| 50 |
h1[i].GetXaxis()->SetTitle("M3 [GeV/c^{2}]");
|
| 51 |
h1[i].GetXaxis()->SetLabelFont(42);
|
| 52 |
h1[i].GetXaxis()->SetLabelSize(0.04);
|
| 53 |
h1[i].GetXaxis()->SetTitleFont(42);
|
| 54 |
h1[i].GetXaxis()->SetTitleSize(0.04);
|
| 55 |
h1[i].GetXaxis()->SetTitleOffset(1.2);
|
| 56 |
|
| 57 |
h1[i].GetYaxis()->SetTitle("Events ( L_{Int}= 20 pb^{-1} )");
|
| 58 |
h1[i].SetMaximum(50.);
|
| 59 |
h1[i].GetYaxis()->SetLabelFont(42);
|
| 60 |
h1[i].GetYaxis()->SetLabelSize(0.04);
|
| 61 |
h1[i].GetYaxis()->SetTitleFont(42);
|
| 62 |
h1[i].GetYaxis()->SetTitleSize(0.04);
|
| 63 |
h1[i].GetYaxis()->SetTitleOffset(1.2);
|
| 64 |
h1[i].Draw();
|
| 65 |
}
|
| 66 |
else
|
| 67 |
h1[i].Draw("same");
|
| 68 |
}
|
| 69 |
|
| 70 |
cv->cd(2);
|
| 71 |
for(int i=0; i < nmembers; i++){
|
| 72 |
TString suffix="PDF_";
|
| 73 |
suffix+=i;
|
| 74 |
|
| 75 |
h2[i] = (TH1D*)f0->Get(TString("Mass/HadronicTop_mass_M3prime_")+suffix);
|
| 76 |
h2[i].SetLineColor(i+1);
|
| 77 |
h2[i].Scale(wttbar);
|
| 78 |
|
| 79 |
if(i==0){
|
| 80 |
h2[i].GetXaxis()->SetTitle("M3' (#chi^{2}) [GeV/c^{2}]");
|
| 81 |
h2[i].GetXaxis()->SetLabelFont(42);
|
| 82 |
h2[i].GetXaxis()->SetLabelSize(0.04);
|
| 83 |
h2[i].GetXaxis()->SetTitleFont(42);
|
| 84 |
h2[i].GetXaxis()->SetTitleSize(0.04);
|
| 85 |
h2[i].GetXaxis()->SetTitleOffset(1.2);
|
| 86 |
|
| 87 |
h2[i].GetYaxis()->SetTitle("Events ( L_{Int}= 20 pb^{-1} )");
|
| 88 |
h2[i].SetMaximum(100.);
|
| 89 |
h2[i].GetYaxis()->SetLabelFont(42);
|
| 90 |
h2[i].GetYaxis()->SetLabelSize(0.04);
|
| 91 |
h2[i].GetYaxis()->SetTitleFont(42);
|
| 92 |
h2[i].GetYaxis()->SetTitleSize(0.04);
|
| 93 |
h2[i].GetYaxis()->SetTitleOffset(1.2);
|
| 94 |
h2[i].Draw();
|
| 95 |
}
|
| 96 |
else
|
| 97 |
h2[i].Draw("same");
|
| 98 |
}
|
| 99 |
|
| 100 |
int nbin = h1[0].GetXaxis()->FindBin(500.)-1;
|
| 101 |
// cout<<"nbin = "<<nbin<<endl;
|
| 102 |
const int nbin1 = nbin+1;
|
| 103 |
double xMax=500.;
|
| 104 |
double xMin=100.;
|
| 105 |
double delta1[nbin1];
|
| 106 |
double delta2[nbin1];
|
| 107 |
double delta3[nbin1];
|
| 108 |
double delta4[nbin1];
|
| 109 |
|
| 110 |
for (int i=1;i<nmembers;i++){
|
| 111 |
TString name1 = "hm3_";
|
| 112 |
TString name2 = "hm3prime_";
|
| 113 |
name1+=i;
|
| 114 |
name2+=i;
|
| 115 |
hm3[i] = new TH1D(name1,"M3 [GeV/c^{2}]",nbin,xMin,xMax);
|
| 116 |
hm3prime[i] = new TH1D(name2,"M3' (#chi^{2}) [GeV/c^{2}]",nbin,xMin,xMax);
|
| 117 |
|
| 118 |
for (int j=1;j<=nbin;j++){ // exclude over/underflow bins
|
| 119 |
ratio[0] = h1[i].GetBinContent(j)/h1[0].GetBinContent(j);
|
| 120 |
ratio[1] = h2[i].GetBinContent(j)/h2[0].GetBinContent(j);
|
| 121 |
|
| 122 |
hm3[i].SetBinContent(j,ratio[0]);
|
| 123 |
hm3prime[i].SetBinContent(j,ratio[1]);
|
| 124 |
}
|
| 125 |
|
| 126 |
cv->cd(3);
|
| 127 |
if (i==1){
|
| 128 |
hm3[i].GetXaxis()->SetTitle("M3 [GeV/c^{2}]");
|
| 129 |
hm3[i].GetXaxis()->SetTitleSize(0.04);
|
| 130 |
hm3[i].GetXaxis()->SetTitleOffset(1.2);
|
| 131 |
hm3[i].GetYaxis()->SetTitle("Ratio");
|
| 132 |
hm3[i].GetYaxis()->SetTitleSize(0.04);
|
| 133 |
hm3[i].GetYaxis()->SetTitleOffset(1.2);
|
| 134 |
hm3[i].SetMaximum(1.2);
|
| 135 |
hm3[i].SetMinimum(0.8);
|
| 136 |
hm3[i].Draw();
|
| 137 |
}
|
| 138 |
else{
|
| 139 |
hm3[i].SetLineColor(i+1);
|
| 140 |
hm3[i].Draw("same");
|
| 141 |
}
|
| 142 |
|
| 143 |
cv->cd(4);
|
| 144 |
if (i==1){
|
| 145 |
hm3prime[i].GetXaxis()->SetTitle("M3' (#chi^{2}) [GeV/c^{2}]");
|
| 146 |
hm3prime[i].GetXaxis()->SetTitleSize(0.04);
|
| 147 |
hm3prime[i].GetXaxis()->SetTitleOffset(1.2);
|
| 148 |
hm3prime[i].GetYaxis()->SetTitle("Ratio");
|
| 149 |
hm3prime[i].GetYaxis()->SetTitleSize(0.04);
|
| 150 |
hm3prime[i].GetYaxis()->SetTitleOffset(1.2);
|
| 151 |
hm3prime[i].SetMaximum(1.2);
|
| 152 |
hm3prime[i].SetMinimum(0.8);
|
| 153 |
hm3prime[i].Draw();
|
| 154 |
}
|
| 155 |
else{
|
| 156 |
hm3prime[i].SetLineColor(i+1);
|
| 157 |
hm3prime[i].Draw("same");
|
| 158 |
}
|
| 159 |
}
|
| 160 |
|
| 161 |
TH1D* hm3SystP = new TH1D("M3SystP","M3 [GeV/c^{2}]",nbin,xMin,xMax);
|
| 162 |
TH1D* hm3SystM = new TH1D("M3SystM","M3 [GeV/c^{2}]",nbin,xMin,xMax);
|
| 163 |
TH1D* hm3primeSystP = new TH1D("M3primeSystP","M3' (#chi^{2}) [GeV/c^{2}]",nbin,xMin,xMax);
|
| 164 |
TH1D* hm3primeSystM = new TH1D("M3primeSystM","M3' (#chi^{2}) [GeV/c^{2}]",nbin,xMin,xMax);
|
| 165 |
|
| 166 |
// Calculate syst. errors
|
| 167 |
for (int i=1;i <= nbin;i++){
|
| 168 |
for (int j=1;j<nmembers;j++){
|
| 169 |
delta1[i-1] += (hm3[j].GetBinContent(i)>1)?pow(hm3[j].GetBinContent(i)-1,2):0;
|
| 170 |
delta2[i-1] += (hm3[j].GetBinContent(i)<1)?pow(hm3[j].GetBinContent(i)-1,2):0;
|
| 171 |
}
|
| 172 |
delta1[i-1] = sqrt(delta1[i-1]);
|
| 173 |
hm3SystP->SetBinContent(i,1+delta1[i-1]);
|
| 174 |
delta2[i-1] = sqrt(delta2[i-1]);
|
| 175 |
hm3SystM->SetBinContent(i,1-delta2[i-1]);
|
| 176 |
}
|
| 177 |
cv->cd(3);
|
| 178 |
hm3SystP->SetLineWidth(2);
|
| 179 |
hm3SystP->Draw("same");
|
| 180 |
hm3SystM->SetLineWidth(2);
|
| 181 |
hm3SystM->Draw("same");
|
| 182 |
|
| 183 |
for (int i=1;i <= nbin;i++){
|
| 184 |
for (int j=1;j<nmembers;j++){
|
| 185 |
delta3[i-1] += (hm3prime[j].GetBinContent(i)>1)?pow(hm3prime[j].GetBinContent(i)-1,2):0;
|
| 186 |
delta4[i-1] += (hm3prime[j].GetBinContent(i)<1)?pow(hm3prime[j].GetBinContent(i)-1,2):0;
|
| 187 |
}
|
| 188 |
delta3[i-1] = sqrt(delta3[i-1]);
|
| 189 |
hm3primeSystP->SetBinContent(i,1+delta3[i-1]);
|
| 190 |
delta4[i-1] = sqrt(delta4[i-1]);
|
| 191 |
hm3primeSystM->SetBinContent(i,1-delta4[i-1]);
|
| 192 |
}
|
| 193 |
cv->cd(4);
|
| 194 |
hm3primeSystP->SetLineWidth(2);
|
| 195 |
hm3primeSystP->Draw("same");
|
| 196 |
hm3primeSystM->SetLineWidth(2);
|
| 197 |
hm3primeSystM->Draw("same");
|
| 198 |
|
| 199 |
}
|
