15 |
|
#include <TPaveStats.h> |
16 |
|
#include <TTree.h> |
17 |
|
|
18 |
< |
TString fileName = "QCDbin4.root"; |
18 |
> |
// Output file |
19 |
> |
TString fileName = "skimmed/QCDbinall_NewRelIso.root"; |
20 |
> |
bool FastSim = false; |
21 |
> |
TString isoname = "New"; |
22 |
> |
// Jet bin: |
23 |
> |
Int_t jbin = 0; // 1-3, 4 for >=4. Put in number less than 5; 0 for all jets |
24 |
> |
double binwidth = 0.01; // 0.02,0.02,0.05,0.1 for New RelIso |
25 |
|
// Histo range |
26 |
|
double xMin = 0.0; |
27 |
< |
double xMax = 2.0; |
22 |
< |
// Jet bin: |
23 |
< |
Int_t jbin = 4; // 1-3, 4 for >=4. Put in number less than 5 |
24 |
< |
double binwidth = 0.1; // 0.02,0.02,0.05,0.1 |
25 |
< |
int nbin = (xMax-xMin)/binwidth; |
27 |
> |
double xMax = 0.0; |
28 |
|
// Specify fitting function: |
27 |
– |
TString fitfunc = "landau"; // this version is cut for Landau |
29 |
|
// Select New or Old reliso: |
29 |
– |
TString isoname = "New"; |
30 |
|
double d0sigCut = 3.; |
31 |
|
|
32 |
|
Int_t fjbin; |
33 |
|
Double_t fna; // Num of QCD events in signal region |
34 |
|
|
35 |
+ |
int nbin; |
36 |
+ |
|
37 |
|
void skim::Loop() |
38 |
|
{ |
39 |
|
// In a ROOT session, you can do: |
44 |
|
// Root > t.Show(); // Show values of entry 12 |
45 |
|
// Root > t.Show(16); // Read and show values of entry 16 |
46 |
|
// Root > t.Loop(); // Loop on all entries |
47 |
+ |
|
48 |
+ |
if ( isoname.Contains("Old") ) { |
49 |
+ |
xMax = 1+binwidth; |
50 |
+ |
nbin = (1.0-xMin)/binwidth+1; |
51 |
+ |
} |
52 |
+ |
else if ( isoname.Contains("New") ) { |
53 |
+ |
xMax = 2.0; |
54 |
+ |
nbin = (xMax-xMin)/binwidth; |
55 |
+ |
} |
56 |
+ |
else { |
57 |
+ |
cout << "No Rel Iso defined!" << endl; |
58 |
+ |
return; |
59 |
+ |
} |
60 |
+ |
cout << "nbin for histo = " << nbin << endl; |
61 |
|
TTree *ftree = new TTree("myTree","myTree"); |
62 |
|
ftree->AutoSave(); |
63 |
|
ftree->Branch("jbin",&fjbin,"jbin/I"); |
74 |
|
// define fit region: |
75 |
|
double AX0 = 0.; |
76 |
|
double AXf = 0.; |
61 |
– |
double BX0 = 0.; |
62 |
– |
double BXf = 0.; |
77 |
|
|
78 |
|
cout<<"Use "<<isoname<<" RelIso"<<endl; |
79 |
|
|
81 |
|
// old |
82 |
|
AX0 = 0.95; |
83 |
|
AXf = 1.0; |
70 |
– |
BX0 = 0.2; |
71 |
– |
BXf = 0.8; |
84 |
|
} |
85 |
|
else if ( isoname.Contains("New") ) { |
86 |
|
// new |
87 |
|
AX0 = 0.0; |
88 |
< |
AXf = 0.053; |
77 |
< |
BX0 = 0.4; |
78 |
< |
BXf = 2.0; |
88 |
> |
AXf = 0.053; // default:0.053 ; loose:0.11 |
89 |
|
} |
90 |
|
|
91 |
< |
// weights |
92 |
< |
double wttbar = 0.0101; |
93 |
< |
double wqcd = 1.3161; |
94 |
< |
double wwjets = 0.0977; |
95 |
< |
double wzjets = 0.078; |
91 |
> |
// Weights: |
92 |
> |
double wttbar,wqcd,wwjets,wzjets; |
93 |
> |
wttbar = wqcd = wwjets = wzjets = 0.; |
94 |
> |
if ( ! FastSim ) { |
95 |
> |
cout << "Full sim sample weights used " << endl; |
96 |
> |
// Full Sim |
97 |
> |
wttbar = 0.0101; |
98 |
> |
wqcd = 1.3161; |
99 |
> |
wwjets = 0.0977; |
100 |
> |
wzjets = 0.078; |
101 |
> |
} |
102 |
> |
else { |
103 |
> |
// Fast Sim |
104 |
> |
cout << "Fast sim sample weights used " << endl; |
105 |
> |
wttbar = 0.0008; |
106 |
> |
wqcd = 0.676; |
107 |
> |
wwjets = 0.0091; |
108 |
> |
wzjets = 0.0094; |
109 |
> |
} |
110 |
|
|
111 |
|
// Na = num of qcd in signal region |
112 |
|
double Na; |
127 |
|
nb = fChain->GetEntry(jentry); nbytes += nb; |
128 |
|
// if (Cut(ientry) < 0) continue; |
129 |
|
bool signalRegion = false; |
106 |
– |
bool fitRegion = false; |
130 |
|
|
131 |
|
TString filename(fChain->GetCurrentFile()->GetName()); |
132 |
|
if (isoname.Contains("Old")) { |
133 |
|
reliso = top_muon_old_reliso[0]; |
134 |
|
if ( reliso > AX0 ) signalRegion = true; |
112 |
– |
else if ( reliso > BX0 && reliso < BXf ) fitRegion = true; |
135 |
|
} |
136 |
|
else if (isoname.Contains("New")) { |
137 |
|
reliso = top_muon_new_reliso[0]; |
138 |
|
if ( reliso < AXf ) signalRegion = true; |
117 |
– |
else if ( reliso > BX0 && reliso < BXf ) fitRegion = true; |
139 |
|
} |
140 |
|
else { |
141 |
|
cout<<"No RelIso defined!!!"<<endl; |
142 |
|
return; |
143 |
|
} |
144 |
< |
|
144 |
> |
|
145 |
|
double d0sig = TMath::Abs(top_muon_d0[0])/top_muon_d0Error[0]; |
146 |
|
bool correctBin = false; // jet bin condition |
147 |
|
bool goodD0sig = false; // d0sig condition |
148 |
|
|
149 |
|
if ( jbin > 4 ) {cout<<"Don't mess around!!! jbin should be less than 5!"<<endl;return;} |
150 |
< |
if (top_njets == jbin && jbin != 4) correctBin = true; |
151 |
< |
else if (top_njets >= jbin && jbin == 4) correctBin = true; |
150 |
> |
if ( jbin == 0 ){ |
151 |
> |
if ( top_njets > jbin ) { correctBin = true; } |
152 |
> |
} |
153 |
> |
else if ( jbin == 4 ){ |
154 |
> |
if ( top_njets >= jbin ) { correctBin = true; } |
155 |
> |
} |
156 |
> |
else { |
157 |
> |
if ( top_njets == jbin ) { correctBin = true; } |
158 |
> |
} |
159 |
> |
|
160 |
|
if ( d0sig < d0sigCut ) goodD0sig = true; |
161 |
|
|
162 |
< |
// Jet bin |
162 |
> |
// Counting |
163 |
|
if( correctBin && goodD0sig ) { |
164 |
< |
if ( filename.Contains("TTJets") ) { |
164 |
> |
if ( filename.Contains("TT") ) { |
165 |
|
h_RelIso_ttbar->Fill(reliso); |
166 |
|
if ( signalRegion ) { Nttbar += wttbar; } |
167 |
|
} |
168 |
< |
if ( filename.Contains("MuPt15") ) { |
168 |
> |
if ( filename.Contains("Mu") ) { |
169 |
|
h_RelIso_qcd->Fill(reliso); |
170 |
|
if ( signalRegion ) { Na += wqcd; Nqcd+= wqcd; } // QCD in signal region |
171 |
|
} |
190 |
|
|
191 |
|
TFile *file0 = new TFile(fileName,"RECREATE"); |
192 |
|
file0->mkdir("histos"); |
193 |
+ |
file0->mkdir("hstack"); |
194 |
|
file0->cd("histos"); |
195 |
|
TCanvas *cv1 = new TCanvas("cv1","cv1",800,600); |
196 |
|
THStack *hs = new THStack("hs","RelIso (stacked)"); |
254 |
|
return; |
255 |
|
} |
256 |
|
h_RelIso_all->SetMinimum(0); |
257 |
< |
|
257 |
> |
|
258 |
> |
file0->cd("hstack"); |
259 |
> |
hs->Write(); |
260 |
> |
|
261 |
|
file0->cd(); |
262 |
|
ftree->Write(); |
230 |
– |
// file0->Write(); |
263 |
|
file0->Close(); |
264 |
|
} |
265 |
|
|