ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/UNL/macros/MyHbbAnalyzer.C
Revision: 1.6
Committed: Fri Oct 21 14:39:57 2011 UTC (13 years, 6 months ago) by malbouis
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.5: +2 -12 lines
Log Message:
correct Naj cut

File Contents

# Content
1 #define MyHbbAnalyzer_cxx
2 #include "MyHbbAnalyzer.h"
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+
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
16 // Root > t.Show(16); // Read and show values of entry 16
17 // Root > t.Loop(); // Loop on all entries
18 //
19
20 // This is the loop skeleton where:
21 // jentry is the global entry number in the chain
22 // ientry is the entry number in the current Tree
23 // Note that the argument to GetEntry must be:
24 // jentry for TChain::GetEntry
25 // ientry for TTree::GetEntry and TBranch::GetEntry
26 //
27 // To read only selected branches, Insert statements like:
28 // METHOD1:
29 // fChain->SetBranchStatus("*",0); // disable all branches
30 // fChain->SetBranchStatus("branchname",1); // activate branchname
31 // METHOD2: replace line
32 // fChain->GetEntry(jentry); //read all branches
33 //by b_branchname->GetEntry(ientry); //read only this branch
34 if (fChain == 0) return;
35
36 Long64_t nentries = fChain->GetEntriesFast();
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 = 10000;
73 for (Long64_t jentry=0; jentry<nentries;jentry++) {
74 Long64_t ientry = LoadTree(jentry);
75 if (ientry < 0) break;
76 nb = fChain->GetEntry(jentry); nbytes += nb;
77 // if (Cut(ientry) < 0) continue;
78
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 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 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 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 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 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 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 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 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 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_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 << 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 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 double pTZ = computePtZ();
260 if (pTZ <= 100) passed_pTZ = false;
261 return passed_pTZ;
262 }// pT_Z cut
263
264 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, 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
276 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, 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;
284 return passedCSV2;
285 }// CSV2
286
287 bool MyHbbAnalyzer::Select_DPhi(){
288 bool passedDPhi = true;
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;
297 return passedNaj;
298 }// Naj
299
300 bool MyHbbAnalyzer::Select_Mjj(){
301 bool passedMjj = true;
302 double Mjj = computeMjj();
303 if (Mjj < 95 || Mjj > 125) passedMjj = false;
304 return passedMjj;
305 }
306 void MyHbbAnalyzer::jetIndices(int &index1, int &index2, double bDisc[9]){
307 int fJetIdx = -1, sJetIdx = -1;
308 double highestBdiscSum = 0;
309 for (int i = 0; i != nJets; ++i){
310 if (bDisc[i] < 0) continue;
311 for (int j = i+1; j != nJets; ++j){
312 if (bDisc[j] < 0) continue;
313 if (bDisc[i] + bDisc[j] > highestBdiscSum){
314 highestBdiscSum = bDisc[i] + bDisc[j];
315 fJetIdx = i;
316 sJetIdx = j;
317 }// if greater sum
318 }// second loop
319 }// 1st loop
320 swap(fJetIdx, index1);
321 swap(sJetIdx, index2);
322 }// jetIndices
323
324 void MyHbbAnalyzer::initializeArr(int evt_count[10000][2], int const max_cut_iter){
325 for (int i = 0; i != max_cut_iter; i++){
326 for (int ijet = 0; ijet != 2; ijet++){
327 evt_count[i][ijet] = 0;
328 }// jets
329 }// cuts
330 }// initializeArr: end function
331
332 void MyHbbAnalyzer::initializeArr2D(int evt_count[100][100], int const max_cut_iter){
333 for (int i = 0; i != max_cut_iter; i++){
334 for (int j = 0; j != max_cut_iter; j++){
335 evt_count[i][j] = 0;
336 }// cuts
337 }// cuts
338
339 }// initializeArr2D: end function
340
341 void MyHbbAnalyzer::basicMuonDistributions(string cut, double ph_weight){
342 // var definition
343 string ph_process;
344
345 if (isZH()) ph_process = "ZH";
346 if (isZZ()) ph_process = "ZZ";
347 if (isDY()) ph_process = "DY";
348
349 for (int i=0; i != nMuons; ++i) {
350 string muonpT = Form("muonPt_%i", i) + ph_process + cut;
351 string muoneta = Form("muonEta_%i", i) + ph_process + cut;
352 string muonphi = Form("muonPhi_%i", i) + ph_process + cut;
353
354 //cout << muonpT << endl;
355
356 fillhisto(muonpT, muonPt[i], ph_weight, "#mu pT", 400, 0, 400);
357 fillhisto(muoneta, muonEta[i], ph_weight, "#mu #eta", 100, -5, 5);
358 fillhisto(muonphi, muonPhi[i], ph_weight, "#mu #phi", 80, -4, 4);
359 }
360
361 double mz = computeMmumu();
362 string Mz = "MZ+" + ph_process + cut;
363 fillhisto(Mz, mz, ph_weight, "", 400, 0, 400);
364
365 }// basicMuonDistributions
366
367 void MyHbbAnalyzer::basicJetDistributions(string cut, double ph_weight){
368 // var definition
369 string ph_process;
370
371 if (isZH()){ ph_process = "ZH";}
372 if (isZZ()){ ph_process = "ZZ";}
373 if (isDY()){ ph_process = "DY";}
374
375 for (int i=0; i != nJets; ++i) {
376 string jetpT = Form("jetPt_%i", i) + ph_process + cut;
377 string jeteta = Form("jetEta_%i", i) + ph_process + cut;
378 string jetphi = Form("jetPhi_%i", i) + ph_process + cut;
379
380 fillhisto(jetpT, jetPt[i], ph_weight, "jet pT", 400, 0, 400);
381 fillhisto(jeteta, jetEta[i], ph_weight, "jet #eta", 100, -5, 5);
382 fillhisto(jetphi, jetPhi[i], ph_weight, "jet #phi", 80, -4, 4);
383 }
384
385 // Mjj
386 double Mjj = computeMjj();
387 string mjj = "Mjj_" + ph_process + cut;
388 fillhisto(mjj, Mjj, ph_weight, "", 1000, 0, 1000);
389
390 }// basicJetDistributions
391
392 double MyHbbAnalyzer::computeMjj(){
393 int firstJetIdx, secondJetIdx;
394 jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
395 // reject event with less than 2 CSV tags
396 if (firstJetIdx == -1 || secondJetIdx == -1) return false;
397 TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
398 TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
399
400 TLorentzVector jj = firstJet_p4 + secondJet_p4;
401
402 return jj.M();
403 }// Mjj
404
405 double MyHbbAnalyzer::computeMmumu(){
406
407 double Mz = 0.;
408 TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
409 TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
410 TLorentzVector Z = firstMu_p4 + secondMu_p4;
411 Mz = Z.M();
412
413 return Mz;
414
415 }// Mmumu
416
417 double MyHbbAnalyzer::computePtjj(){
418 double ptjj = 0.;
419 // get jet indices that maximize the b-tag discriminant
420 int firstJetIdx, secondJetIdx;
421 jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
422 // reject event with less than 2 CSV tags
423 if (firstJetIdx == -1 || secondJetIdx == -1) return false;
424 // get the dijet system
425 TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
426 TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
427 TLorentzVector jj = firstJet_p4 + secondJet_p4;
428 ptjj = jj.Pt();
429
430 return ptjj;
431 }// Ptjj
432
433 double MyHbbAnalyzer::computePtZ(){
434 double ptz = 0.;
435 TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
436 TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
437 TLorentzVector Z = firstMu_p4 + secondMu_p4;
438 ptz = Z.Pt();
439
440 return ptz;
441 }// PtZ
442
443 double MyHbbAnalyzer::computeDPhi(){
444
445 double dphi = 99.;
446 TLorentzVector firstMu_p4 = TLorentzVector(muonPx[0], muonPy[0], muonPz[0], muonP[0]);
447 TLorentzVector secondMu_p4 = TLorentzVector(muonPx[1], muonPy[1], muonPz[1], muonP[1]);
448 TLorentzVector Z = firstMu_p4 + secondMu_p4;
449 int firstJetIdx, secondJetIdx;
450 jetIndices(firstJetIdx, secondJetIdx, bDisc_CSV);
451 // reject event with less than 2 CSV tags
452 if (firstJetIdx == -1 || secondJetIdx == -1) return false;
453 TLorentzVector firstJet_p4 = TLorentzVector(jetPx[firstJetIdx], jetPy[firstJetIdx], jetPz[firstJetIdx], jetP[firstJetIdx]);
454 TLorentzVector secondJet_p4 = TLorentzVector(jetPx[secondJetIdx], jetPy[secondJetIdx], jetPz[secondJetIdx], jetP[secondJetIdx]);
455 TLorentzVector jj = firstJet_p4 + secondJet_p4;
456 dphi = Z.DeltaPhi(jj);
457
458 return dphi;
459
460 }// dPhi
461 double MyHbbAnalyzer::DeltaPhi(const double phi1, const double phi2){
462 double delta = phi1-phi2;
463 if (delta > TMath::Pi()) delta -= 2.0*TMath::Pi();
464 if (delta < (-TMath::Pi())) delta += 2.0*TMath::Pi();
465 return fabs(delta);
466 }