ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/UNL/macros/MyHbbAnalyzer.C
(Generate patch)

Comparing UserCode/UNL/macros/MyHbbAnalyzer.C (file contents):
Revision 1.2 by malbouis, Wed Aug 10 18:04:26 2011 UTC vs.
Revision 1.5 by malbouis, Thu Oct 13 15:10:45 2011 UTC

# Line 3 | Line 3
3   #include <TH2.h>
4   #include <TStyle.h>
5   #include <TCanvas.h>
6 + #include <iostream>
7 + #include <fstream>
8  
9   void MyHbbAnalyzer::Loop()
10   {
11   //   In a ROOT session, you can do:
12 < //      Root > .L MyHbbAnalyzer.C
12 > //      Root > .L MyHbbAnalyzer.C+
13   //      Root > MyHbbAnalyzer t
14   //      Root > t.GetEntry(12); // Fill t data members with entry number 12
15   //      Root > t.Show();       // Show values of entry 12
# Line 35 | Line 37 | void MyHbbAnalyzer::Loop()
37  
38     Long64_t nbytes = 0, nb = 0;
39  
40 +  
41 +   ////////
42 +   // get samples weight:
43 +   ////////
44 +   double weight = 1.0;
45 +   double ZH_weight = SetTTreeWeight("ZH");
46 +   double ZZ_weight = SetTTreeWeight("ZZ");
47 +   double DY_weight = SetTTreeWeight("DY");
48 +
49 +   ////////
50 +   // initialize variables for significance
51 +   ////////
52 +  
53     ////////
54     // cut and count variables
55     ////////
56     int allEvts = 0;
57     int preSelection = 0;
58     int pT_jj = 0, pT_Z = 0, CSV1 = 0, CSV2 = 0, DPhi_ZH = 0, N_aj = 0, M_jj = 0;
59 +   int allEvts_zz = 0;
60 +   int preSelection_zz = 0;
61 +   int pT_jj_zz = 0, pT_Z_zz = 0, CSV1_zz = 0, CSV2_zz = 0, DPhi_ZH_zz = 0, N_aj_zz = 0, M_jj_zz = 0;
62 +   int allEvts_z = 0;
63 +   int preSelection_z = 0;
64 +   int pT_jj_z = 0, pT_Z_z = 0, CSV1_z = 0, CSV2_z = 0, DPhi_ZH_z = 0, N_aj_z = 0, M_jj_z = 0;
65 +
66     ////////
67     // loop on entries
68     ////////
69     pct = 0;
70     total_entries = fChain->GetEntries();
71  
72 <   //nentries = 1000;
72 >   //nentries = 10000;
73     for (Long64_t jentry=0; jentry<nentries;jentry++) {
74        Long64_t ientry = LoadTree(jentry);
75        if (ientry < 0) break;
# Line 57 | Line 79 | void MyHbbAnalyzer::Loop()
79        // jobProgress:
80        jobProgress(jentry);
81        
82 +      // set weights
83 +      if (isZH()) weight = ZH_weight;
84 +      if (isZZ()) weight = ZZ_weight;
85 +      if (isDY()) weight = DY_weight;
86 +
87        // count events:
88 <      allEvts++;
89 <      basicMuonDistributions("allEvts");
90 <      basicJetDistributions("allEvts");
91 <      
88 >      if (isZH()) allEvts++;
89 >      if (isZZ()) allEvts_zz++;
90 >      if (isDY()) allEvts_z++;
91 >      basicMuonDistributions("allEvts", weight);
92 >      basicJetDistributions("allEvts", weight);
93 >
94 >
95        // pre-selection cuts:
96        if (!preSelect()) continue;
97 <      preSelection++;
98 <      basicMuonDistributions("preSel");
99 <      basicJetDistributions("preSel");
97 >      if (isZH()) preSelection++;
98 >      if (isZZ()) preSelection_zz++;
99 >      if (isDY()) preSelection_z++;
100 >      basicMuonDistributions("preSel", weight);
101 >      basicJetDistributions("preSel", weight);
102        
103        // pt_jj cut:
104        if (!Select_pTjj()) continue;
105 <      pT_jj++;
106 <      basicMuonDistributions("pT_jj");
107 <      basicJetDistributions("pT_jj");
105 >      if (isZH()) pT_jj++;
106 >      if (isZZ()) pT_jj_zz++;
107 >      if (isDY()) pT_jj_z++;
108 >      basicMuonDistributions("pT_jj", weight);
109 >      basicJetDistributions("pT_jj", weight);
110  
111        // pT_Z cut:
112        if (!Select_pTZ()) continue;
113 <      pT_Z++;
114 <      basicMuonDistributions("pT_Z");
115 <      basicJetDistributions("pT_Z");
116 <
113 >      if (isZH()) pT_Z++;
114 >      if (isZZ()) pT_Z_zz++;
115 >      if (isDY()) pT_Z_z++;
116 >      basicMuonDistributions("pT_Z", weight);
117 >      basicJetDistributions("pT_Z", weight);
118 >      
119        // CSV1 cut:
120        if (!Select_CSV1()) continue;
121 <      CSV1++;
122 <      basicMuonDistributions("CSV1");
123 <      basicJetDistributions("CSV1");
124 <
121 >      if (isZH()) CSV1++;
122 >      if (isZZ()) CSV1_zz++;
123 >      if (isDY()) CSV1_z++;
124 >      basicMuonDistributions("CSV1", weight);
125 >      basicJetDistributions("CSV1", weight);
126 >      
127        // CSV2 cut:
128        if (!Select_CSV2()) continue;
129 <      CSV2++;
130 <      basicMuonDistributions("CSV2");
131 <      basicJetDistributions("CSV2");
132 <
129 >      if (isZH()) CSV2++;
130 >      if (isZZ()) CSV2_zz++;
131 >      if (isDY()) CSV2_z++;
132 >      basicMuonDistributions("CSV2", weight);
133 >      basicJetDistributions("CSV2", weight);
134 >      
135        // DPhi cut:
136        if (!Select_DPhi()) continue;
137 <      DPhi_ZH++;
138 <      basicMuonDistributions("DPhi");
139 <      basicJetDistributions("DPhi");
140 <
137 >      if (isZH()) DPhi_ZH++;
138 >      if (isZZ()) DPhi_ZH_zz++;
139 >      if (isDY()) DPhi_ZH_z++;
140 >      basicMuonDistributions("DPhi", weight);
141 >      basicJetDistributions("DPhi", weight);
142 >            
143        // Naj cut:
144        if (!Select_Naj()) continue;
145 <      N_aj++;
146 <      basicMuonDistributions("Naj");
147 <      basicJetDistributions("Naj");
148 <
145 >      if (isZH()) N_aj++;
146 >      if (isZZ()) N_aj_zz++;
147 >      if (isDY()) N_aj_z++;
148 >      basicMuonDistributions("Naj", weight);
149 >      basicJetDistributions("Naj", weight);
150 >    
151        // Mjj cut:
152        if (!Select_Mjj()) continue;
153 <      M_jj++;
154 <      basicMuonDistributions("Mjj");
155 <      basicJetDistributions("Mjj");
156 <
153 >      if (isZH()) M_jj++;
154 >      if (isZZ()) M_jj_zz++;
155 >      if (isDY()) M_jj_z++;
156 >      basicMuonDistributions("Mjj", weight);
157 >      basicJetDistributions("Mjj", weight);
158 >      
159     }// loop on entries
160  
161 +
162     ////////
163     // write histos to file:
164     ////////
165 <   TFile *f = new TFile("hbb.root","RECREATE");
165 >   TFile *f = new TFile("hbb_CSV.root","RECREATE");
166     for (map<std::string,TH1*>::iterator it=histmap.begin(); it!=histmap.end();it++) {
167       (*it).second->Write();
168     }
169 +   for (map<string,TH2*>::iterator it2=bidimhistmap.begin(); it2!=bidimhistmap.end();it2++) {
170 +     (*it2).second->Write();
171 +     delete (*it2).second;
172 +   }
173 +
174     f->Write();
175     f->Close();
176 +
177    
178     ////////
179     // print stats
180     ////////
181 <   cout << "allEvts: " << allEvts << endl;
182 <   cout << "preSelection: " << preSelection << endl;
183 <   cout << "pT_jj: " << pT_jj << endl;
184 <   cout << "pT_Z: " << pT_Z << endl;
185 <   cout << "CSV1: " << CSV1 << endl;
186 <   cout << "CSV2: " << CSV2 << endl;
187 <   cout << "DPhi: " << DPhi_ZH << endl;
188 <   cout << "Naj: " << N_aj << endl;
189 <   cout << "Mjj: " << M_jj << endl;
190 <   cout << "signal eff. PRE-SELECT: " << (double)preSelection/allEvts << endl;
191 <   cout << "signal eff. pT_jj: " << (double)pT_jj/preSelection << endl;
192 <   cout << "signal eff. pT_Z: " << (double)pT_Z/pT_jj << endl;
193 <   cout << "signal eff. CSV1: " << (double)CSV1/pT_Z << endl;
194 <   cout << "signal eff. CSV2: " << (double)CSV2/CSV1 << endl;
195 <   cout << "signal eff. DPhi: " << (double)DPhi_ZH/CSV2 << endl;
196 <   cout << "signal eff. Naj: " << (double)N_aj/DPhi_ZH << endl;
197 <   cout << "signal eff. Mjj: " << (double)M_jj/N_aj << endl;
198 <  
181 >   cout << endl;
182 >   cout << endl;
183 >   cout << "**** stats (no normalization) ****" << endl;
184 >   cout << "allEvts, ZH: " << allEvts << " , ZZ: " << allEvts_zz << ", DY: " << allEvts_z << endl;
185 >   cout << "preSelection, ZH: " << preSelection << " , ZZ: " << preSelection_zz << ", DY: " << preSelection_z <<  endl;
186 >   cout << "pT_jj, ZH: " << pT_jj << " , ZZ: " << pT_jj_zz << ", DY: " << pT_jj_z <<  endl;
187 >   cout << "pT_Z, ZH: " << pT_Z << " , ZZ: " << pT_Z_zz << ", DY: " << pT_Z_z <<  endl;
188 >   cout << "CSV1, ZH: " << CSV1 << " , ZZ: " << CSV1_zz << ", DY: " << CSV1_z <<  endl;
189 >   cout << "CSV2, ZH: " << CSV2 << " , ZZ: " << CSV2_zz << ", DY: " << CSV2_z <<  endl;
190 >   cout << "DPhi, ZH: " << DPhi_ZH << " , ZZ: " << DPhi_ZH_zz << ", DY: " << DPhi_ZH_z <<  endl;
191 >   cout << "Naj, ZH: " << N_aj << " , ZZ: " << N_aj_zz << ", DY: " << N_aj_z <<  endl;
192 >   cout << "Mjj, ZH: " << M_jj << " , ZZ: " << M_jj_zz << ", DY: " << M_jj_z <<  endl;
193 >   cout << endl;
194 >   cout << endl;
195 >   cout << "**** stats normalized to 10 fb-1 ****" << endl;
196 >   cout << "allEvts: " << allEvts*ZH_weight <<  ", ZZ allEvts: " << allEvts_zz*ZZ_weight <<  ", DY allEvts: " << allEvts_z*DY_weight << endl;
197 >   cout << "preSelection: " << preSelection*ZH_weight << ", ZZ preSelection: " << preSelection_zz*ZZ_weight << ", DY preSelection: " << preSelection_z*DY_weight << endl;
198 >   cout << "pT_jj: " << pT_jj*ZH_weight << ", ZZ pT_jj: " << pT_jj_zz*ZZ_weight << ", DY pT_jj: " << pT_jj_z*DY_weight << endl;
199 >   cout << "pT_Z: " << pT_Z*ZH_weight << ", ZZ pT_Z: " << pT_Z_zz*ZZ_weight << ", DY pT_Z: " << pT_Z_z*DY_weight << endl;
200 >   cout << "CSV1: " << CSV1*ZH_weight << ", ZZ CSV1: " << CSV1_zz*ZZ_weight << ", DY CSV1: " << CSV1_z*DY_weight << endl;
201 >   cout << "CSV2: " << CSV2*ZH_weight << ", ZZ CSV2: " << CSV2_zz*ZZ_weight << ", DY CSV2: " << CSV2_z*DY_weight << endl;
202 >   cout << "DPhi: " << DPhi_ZH*ZH_weight << ", ZZ DPhi: " << DPhi_ZH_zz*ZZ_weight << ", DY DPhi: " << DPhi_ZH_z*DY_weight << endl;
203 >   cout << "Naj: " << N_aj*ZH_weight << ", ZZ Naj: " << N_aj_zz*ZZ_weight << ", DY Naj: " << N_aj_z*DY_weight << endl;
204 >   cout << "Mjj: " << M_jj*ZH_weight << ", ZZ Mjj: " << M_jj_zz*ZZ_weight << ", DY Mjj: " << M_jj_z*DY_weight << endl;
205 >   cout << endl;
206 >   cout << "signal eff. PRE-SELECT: " << (double)preSelection/allEvts << ", ZZ eff. PRE-SELECT: " << (double)preSelection_zz/allEvts_zz << ", DY eff. PRE-SELECT: " << (double)preSelection_z/allEvts_z << endl;
207 >   cout << "signal eff. pT_jj: " << (double)pT_jj/preSelection << ", ZZ eff. pT_jj: " << (double)pT_jj_zz/preSelection_zz << ", DY eff. pT_jj: " << (double)pT_jj_z/preSelection_z << endl;
208 >   cout << "signal eff. pT_Z: " << (double)pT_Z/pT_jj <<  ", ZZ signal eff. pT_Z: " << (double)pT_Z_zz/pT_jj_zz <<  ", DY signal eff. pT_Z: " << (double)pT_Z_z/pT_jj_z << endl;
209 >   cout << "signal eff. CSV1: " << (double)CSV1/pT_Z << ", ZZ eff. CSV1: " << (double)CSV1_zz/pT_Z_zz << ", DY eff. CSV1: " << (double)CSV1_z/pT_Z_z << endl;
210 >   cout << "signal eff. CSV2: " << (double)CSV2/CSV1 << ", ZZ eff. CSV2: " << (double)CSV2_zz/CSV1_zz << ", DY eff. CSV2: " << (double)CSV2_z/CSV1_z << endl;
211 >   cout << "signal eff. DPhi: " << (double)DPhi_ZH/CSV2 << ", ZZ eff. DPhi: " << (double)DPhi_ZH_zz/CSV2_zz << ", DY eff. DPhi: " << (double)DPhi_ZH_z/CSV2_z << endl;
212 >   cout << "signal eff. Naj: " << (double)N_aj/DPhi_ZH << ", ZZ eff. Naj: " << (double)N_aj_zz/DPhi_ZH_zz << ", DY eff. Naj: " << (double)N_aj_z/DPhi_ZH_z << endl;
213 >   cout << "signal eff. Mjj: " << (double)M_jj/N_aj << ", ZZ eff. Mjj: " << (double)M_jj_zz/N_aj_zz << ", DY eff. Mjj: " << (double)M_jj_z/N_aj_z << endl;
214 >   cout << endl;
215 >   cout << "**** signal only ****" << endl;
216 >   cout << "signal eff. PRE-SELECT: " << (double)preSelection*100/allEvts << " +- " << relativeError((double)preSelection, (double)allEvts)*100 << endl;
217 >   cout << "signal eff. pT_jj: " << (double)pT_jj*100/preSelection << " +- " << relativeError((double)pT_jj, (double)preSelection)*100 << endl;
218 >   cout << "signal eff. pT_Z: " << (double)pT_Z*100/pT_jj << " +- " << relativeError((double)pT_Z, (double)pT_jj)*100 << endl;
219 >   cout << "signal eff. CSV1: " << (double)CSV1*100/pT_Z << " +- " << relativeError((double)CSV1, (double)pT_Z)*100 << endl;
220 >   cout << "signal eff. CSV2: " << (double)CSV2*100/CSV1 << " +- " << relativeError((double)CSV2, (double)CSV1)*100 << endl;
221 >   cout << "signal eff. DPhi: " << (double)DPhi_ZH*100/CSV2 << " +- " << relativeError((double)DPhi_ZH, (double)CSV2)*100 << endl;
222 >   cout << "signal eff. Naj: " << (double)N_aj*100/DPhi_ZH << " +- " << relativeError((double)N_aj, (double)DPhi_ZH)*100 << endl;
223 >   cout << "signal eff. Mjj: " << (double)M_jj*100/N_aj << " +- " << relativeError((double)M_jj, (double)N_aj)*100 << endl;
224 >   cout << "Total signal eff.: " << (double)M_jj*100/allEvts << " +- " << relativeError((double)M_jj, (double)allEvts)*100 << endl;
225 >   cout << endl;
226 >   cout << "**** TOTAL SIGNIFICANCE, 95 < M_jj < 125 ****" << endl;
227 >   cout << "S = S/sqrt{B} = " << ((double)M_jj*ZH_weight)/sqrt((M_jj_zz*ZZ_weight) + (M_jj_z*DY_weight)) << endl;
228 >   cout << endl;
229 >   cout << "**** TOTAL SIGNIFICANCE before Mjj cut ****" << endl;
230 >   cout << "S = S/sqrt{B} = " << ((double)N_aj*ZH_weight)/sqrt((N_aj_zz*ZZ_weight) + (N_aj_z*DY_weight)) << endl;
231 >   cout << endl;
232 >   cout << endl;
233  
234   }// Loop
235   bool MyHbbAnalyzer::preSelect(){
236    bool passedPreSel = true;
237 +  // trigger requirement:
238 +  if (!HLT_IsoMu17_v5_trig) return false;
239    if (nMuons < 2) passedPreSel = false;
240    if (nJets < 2) passedPreSel = false;
241 +  int firstJetIdx, secondJetIdx;
242 +  jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
243 +  // reject events with less than 2 CSV tags
244 +  if (firstJetIdx == -1 || secondJetIdx == -1) return false;
245 +  double mz = computeMmumu();
246 +  if (mz <= 75. || mz >= 105.) passedPreSel = false;
247    return passedPreSel;
248   }// preSelect
249  
250   bool MyHbbAnalyzer::Select_pTjj(){
251    bool passedSel_pTjj = true;
252 <  // get jet indices that maximize the CSV discriminant
253 <  int firstJetIdx, secondJetIdx;
159 <  jetIndices(firstJetIdx, secondJetIdx);
160 <  // reject event with less than 2 CSV tags
161 <  if (firstJetIdx == -1 || secondJetIdx == -1) return false;
162 <  // get the dijet system
163 <  TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
164 <  TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
165 <  TLorentzVector jj = firstJet_p4 + secondJet_p4;
166 <  if (jj.Pt() <= 100) passedSel_pTjj = false;
252 >  double pTjj = computePtjj();
253 >  if (pTjj <= 100) passedSel_pTjj = false;
254    return passedSel_pTjj;
255   }// pT_jj cut
256  
257   bool MyHbbAnalyzer::Select_pTZ(){
258    bool passed_pTZ = true;
259 <  TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
260 <  TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
174 <  TLorentzVector Z = firstMu_p4 + secondMu_p4;
175 <  if (Z.Pt() <= 100) passed_pTZ = false;
259 >  double pTZ = computePtZ();
260 >  if (pTZ <= 100) passed_pTZ = false;
261    return passed_pTZ;
262   }// pT_Z cut
263  
# Line 180 | Line 265 | bool MyHbbAnalyzer::Select_CSV1(){
265    // max CSV tagged jet should satisfy the tight Operating Point (CSVT > 0.898)
266    bool passedCSV1 = true;
267    int firstJetIdx, secondJetIdx;
268 <  jetIndices(firstJetIdx, secondJetIdx);
268 >  jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
269    double bDiscMax;
270    bDiscMax = bDisc_CSV[firstJetIdx] > bDisc_CSV[secondJetIdx] ? bDisc_CSV[firstJetIdx]:bDisc_CSV[secondJetIdx];
271    if (bDiscMax <= 0.898) passedCSV1 = false;
272 +  //if (bDiscMax <= 0.866) passedCSV1 = false;
273    return passedCSV1;
274   }// CSV1
275  
# Line 191 | Line 277 | bool MyHbbAnalyzer::Select_CSV2(){
277    // min CSV tagged jet should satisfy CSV > 0.5
278    bool passedCSV2 = true;
279    int firstJetIdx, secondJetIdx;
280 <  jetIndices(firstJetIdx, secondJetIdx);
280 >  jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
281    double bDiscMin;
282    bDiscMin = bDisc_CSV[firstJetIdx] < bDisc_CSV[secondJetIdx] ? bDisc_CSV[firstJetIdx]:bDisc_CSV[secondJetIdx];
283    if (bDiscMin <= 0.5) passedCSV2 = false;
# Line 200 | Line 286 | bool MyHbbAnalyzer::Select_CSV2(){
286  
287   bool MyHbbAnalyzer::Select_DPhi(){
288    bool passedDPhi = true;
289 <  TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
290 <  TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
205 <  TLorentzVector Z = firstMu_p4 + secondMu_p4;
206 <  int firstJetIdx, secondJetIdx;
207 <  jetIndices(firstJetIdx, secondJetIdx);
208 <  TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
209 <  TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
210 <  TLorentzVector jj = firstJet_p4 + secondJet_p4;
211 <  if (fabs(Z.DeltaPhi(jj)) <= 2.90) passedDPhi = false;
289 >  double DPhi = computeDPhi();
290 >  if (fabs(DPhi) <= 2.90) passedDPhi = false;
291    return passedDPhi;
292   }// DPhi
293  
294   bool MyHbbAnalyzer::Select_Naj(){
295    bool passedNaj = true;
296 <  if (nJets > 3) passedNaj = false;
296 >  ////if (nJets > 3) passedNaj = false;
297 >  int bjets = 0;
298 >  for(int j = 0; j != nJets; ++j){
299 >    if (bDisc_CSV[j] < 0.244) continue;
300 >    bjets++;
301 >  }// jet loop
302 >  if (bjets > 3) passedNaj = false;
303    return passedNaj;
304   }// Naj
305  
306   bool  MyHbbAnalyzer::Select_Mjj(){
307    bool passedMjj = true;
308 <  int firstJetIdx, secondJetIdx;
309 <  jetIndices(firstJetIdx, secondJetIdx);
225 <  TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
226 <  TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
227 <  TLorentzVector jj = firstJet_p4 + secondJet_p4;
228 <  if (jj.M() < 95 || jj.M() > 125) passedMjj = false;
308 >  double Mjj = computeMjj();
309 >  if (Mjj < 95 || Mjj > 125) passedMjj = false;
310    return passedMjj;
311   }
312 < void MyHbbAnalyzer::jetIndices(int &index1, int &index2){
312 > void MyHbbAnalyzer::jetIndices(int &index1, int &index2, double bDisc[9]){
313    int fJetIdx = -1, sJetIdx = -1;
314    double highestBdiscSum = 0;
315    for (int i = 0; i != nJets; ++i){
316 <    if (bDisc_CSV[i] < 0) continue;
316 >    if (bDisc[i] < 0) continue;
317      for (int j = i+1; j != nJets; ++j){
318 <      if (bDisc_CSV[j] < 0) continue;
319 <      if (bDisc_CSV[i] + bDisc_CSV[j] > highestBdiscSum){
320 <        highestBdiscSum = bDisc_CSV[i] + bDisc_CSV[j];
318 >      if (bDisc[j] < 0) continue;
319 >      if (bDisc[i] + bDisc[j] > highestBdiscSum){
320 >        highestBdiscSum = bDisc[i] + bDisc[j];
321          fJetIdx = i;
322          sJetIdx = j;
323        }// if greater sum
# Line 246 | Line 327 | void MyHbbAnalyzer::jetIndices(int &inde
327    swap(sJetIdx, index2);
328   }// jetIndices
329  
330 < void MyHbbAnalyzer::basicMuonDistributions(string cut){
330 > void MyHbbAnalyzer::initializeArr(int evt_count[10000][2], int const max_cut_iter){
331 >  for (int i = 0; i != max_cut_iter; i++){
332 >    for (int ijet = 0; ijet != 2; ijet++){
333 >     evt_count[i][ijet] = 0;
334 >    }// jets
335 >  }// cuts
336 > }// initializeArr: end function
337 >
338 > void MyHbbAnalyzer::initializeArr2D(int evt_count[100][100], int const max_cut_iter){
339 >  for (int i = 0; i != max_cut_iter; i++){
340 >    for (int j = 0; j != max_cut_iter; j++){
341 >      evt_count[i][j] = 0;
342 >    }// cuts
343 >  }// cuts
344 >
345 > }// initializeArr2D: end function
346 >
347 > void MyHbbAnalyzer::basicMuonDistributions(string cut, double ph_weight){
348    // var definition
349    string ph_process;
350    
# Line 255 | Line 353 | void MyHbbAnalyzer::basicMuonDistributio
353    if (isDY()) ph_process = "DY";
354        
355    for (int i=0; i != nMuons; ++i) {
356 <    string muonpT = "muonPt_" + ph_process + cut;
357 <    string muoneta = "muonEta_" + ph_process + cut;
358 <    string muonphi = "muonPhi_" + ph_process + cut;
356 >    string muonpT = Form("muonPt_%i", i) + ph_process + cut;
357 >    string muoneta = Form("muonEta_%i", i) + ph_process + cut;
358 >    string muonphi = Form("muonPhi_%i", i) + ph_process + cut;
359  
360      //cout << muonpT << endl;
361            
362 <    fillhisto(muonpT, muonPt[i], 1., "#mu pT", 400, 0, 400);
363 <    fillhisto(muoneta, muonEta[i], 1., "#mu #eta", 100, -5, 5);
364 <    fillhisto(muonphi, muonPhi[i], 1., "#mu #phi", 80, -4, 4);
362 >    fillhisto(muonpT, muonPt[i], ph_weight, "#mu pT", 400, 0, 400);
363 >    fillhisto(muoneta, muonEta[i], ph_weight, "#mu #eta", 100, -5, 5);
364 >    fillhisto(muonphi, muonPhi[i], ph_weight, "#mu #phi", 80, -4, 4);
365    }
366 +
367 +  double mz = computeMmumu();
368 +  string Mz = "MZ+" + ph_process + cut;
369 +  fillhisto(Mz, mz, ph_weight, "", 400, 0, 400);
370 +  
371   }// basicMuonDistributions
372  
373 < void MyHbbAnalyzer::basicJetDistributions(string cut){
373 > void MyHbbAnalyzer::basicJetDistributions(string cut, double ph_weight){
374    // var definition
375    string ph_process;
376    
377 <  if (isZH()) ph_process = "ZH";
378 <  if (isZZ()) ph_process = "ZZ";
379 <  if (isDY()) ph_process = "DY";
377 >  if (isZH()){ ph_process = "ZH";}
378 >  if (isZZ()){ ph_process = "ZZ";}
379 >  if (isDY()){ ph_process = "DY";}
380        
381    for (int i=0; i != nJets; ++i) {
382 <    string jetpT = "jetPt_" + ph_process + cut;
383 <    string jeteta = "jetEta_" + ph_process + cut;
384 <    string jetphi = "jetPhi_" + ph_process + cut;
282 <
283 <    //cout << muonpT << endl;
382 >    string jetpT = Form("jetPt_%i", i) + ph_process + cut;
383 >    string jeteta = Form("jetEta_%i", i) + ph_process + cut;
384 >    string jetphi = Form("jetPhi_%i", i) + ph_process + cut;
385            
386 <    fillhisto(jetpT, jetPt[i], 1., "jet pT", 400, 0, 400);
387 <    fillhisto(jeteta, jetEta[i], 1., "jet #eta", 100, -5, 5);
388 <    fillhisto(jetphi, jetPhi[i], 1., "jet #phi", 80, -4, 4);
386 >    fillhisto(jetpT, jetPt[i], ph_weight, "jet pT", 400, 0, 400);
387 >    fillhisto(jeteta, jetEta[i], ph_weight, "jet #eta", 100, -5, 5);
388 >    fillhisto(jetphi, jetPhi[i], ph_weight, "jet #phi", 80, -4, 4);
389    }
390 +
391 +  // Mjj
392 +  double Mjj = computeMjj();
393 +  string mjj = "Mjj_" + ph_process + cut;
394 +  fillhisto(mjj, Mjj, ph_weight, "", 1000, 0, 1000);
395 +  
396   }// basicJetDistributions
397 +
398 + double MyHbbAnalyzer::computeMjj(){
399 +  int firstJetIdx, secondJetIdx;
400 +  jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
401 +  // reject event with less than 2 CSV tags
402 +  if (firstJetIdx == -1 || secondJetIdx == -1) return false;
403 +  TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
404 +  TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
405 +
406 +  TLorentzVector jj = firstJet_p4 + secondJet_p4;
407 +
408 +  return jj.M();
409 + }// Mjj
410 +
411 + double MyHbbAnalyzer::computeMmumu(){
412 +
413 +  double Mz = 0.;
414 +  TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
415 +  TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
416 +  TLorentzVector Z = firstMu_p4 + secondMu_p4;
417 +  Mz = Z.M();
418 +
419 +  return Mz;
420 +
421 + }// Mmumu
422 +
423 + double MyHbbAnalyzer::computePtjj(){
424 +  double ptjj = 0.;
425 +  // get jet indices that maximize the b-tag discriminant
426 +  int firstJetIdx, secondJetIdx;
427 +  jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
428 +  // reject event with less than 2 CSV tags
429 +  if (firstJetIdx == -1 || secondJetIdx == -1) return false;
430 +  // get the dijet system
431 +  TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
432 +  TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
433 +  TLorentzVector jj = firstJet_p4 + secondJet_p4;
434 +  ptjj = jj.Pt();
435 +
436 +  return ptjj;
437 + }// Ptjj
438 +
439 + double MyHbbAnalyzer::computePtZ(){
440 +  double ptz = 0.;
441 +  TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
442 +  TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
443 +  TLorentzVector Z = firstMu_p4 + secondMu_p4;
444 +  ptz = Z.Pt();
445 +
446 +  return ptz;
447 + }// PtZ
448 +
449 + double MyHbbAnalyzer::computeDPhi(){
450 +
451 +  double dphi = 99.;
452 +  TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
453 +  TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
454 +  TLorentzVector Z = firstMu_p4 + secondMu_p4;
455 +  int firstJetIdx, secondJetIdx;
456 +  jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
457 +  // reject event with less than 2 CSV tags
458 +  if (firstJetIdx == -1 || secondJetIdx == -1) return false;
459 +  TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
460 +  TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
461 +  TLorentzVector jj = firstJet_p4 + secondJet_p4;
462 +  dphi = Z.DeltaPhi(jj);
463 +
464 +  return dphi;
465 +
466 + }// dPhi
467 + double MyHbbAnalyzer::DeltaPhi(const double phi1, const double phi2){
468 +  double delta = phi1-phi2;
469 +  if (delta > TMath::Pi()) delta -= 2.0*TMath::Pi();
470 +  if (delta < (-TMath::Pi())) delta += 2.0*TMath::Pi();
471 +  return fabs(delta);
472 + }
473 +
474 +
475 +
476 +

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines