ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/FGolf/Tools/EfficiencyModelTools.cc
Revision: 1.2
Committed: Mon May 28 22:11:53 2012 UTC (12 years, 11 months ago) by fgolf
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +23 -25 lines
Error occurred while calculating annotation data.
Log Message:
update for new gen particle struct; change btag eta cut from 2.5 to 2.4; update for 2012 lepton selections

File Contents

# Content
1 #include "EfficiencyModelTools.h"
2
3 #include "CORE/CMS2.h"
4
5 #include "CORE/ssSelections.h"
6 #include "CORE/mcSelections.h"
7
8 #include "Math/VectorUtil.h"
9
10 std::vector<unsigned int> efftools::getRecoJets () {
11
12 const float mupt = 20.;
13 const float elpt = 20.;
14
15 std::vector<unsigned int> goodJets_;
16 for (unsigned int jidx = 0; jidx < cms2.pfjets_p4().size(); jidx++) {
17
18 float jet_eta = cms2.pfjets_p4().at(jidx).eta();
19 float jet_pt = cms2.pfjets_p4().at(jidx).pt() * cms2.pfjets_corL1FastL2L3().at(jidx);
20
21 if (fabs(jet_eta) > 2.4)
22 continue;
23 if (jet_pt < 40.)
24 continue;
25
26 bool jetIsLep = false;
27 for (unsigned int eidx = 0; eidx < cms2.els_p4().size(); eidx++) {
28
29 if (cms2.els_p4().at(eidx).pt() < elpt)
30 continue;
31 if (fabs(cms2.els_p4().at(eidx).eta()) > 2.4)
32 continue;
33 if (!samesign::isNumeratorLepton(11, eidx))
34 continue;
35
36 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.els_p4().at(eidx)) < 0.4) {
37 jetIsLep = true;
38 break;
39 }
40 }
41
42 if (jetIsLep)
43 continue;
44
45 for (unsigned int midx = 0; midx < cms2.mus_p4().size(); midx++) {
46
47 if (cms2.mus_p4().at(midx).pt() < mupt)
48 continue;
49 if (fabs(cms2.mus_p4().at(midx).eta()) > 2.4)
50 continue;
51 if (!samesign::isNumeratorLepton(13, midx))
52 continue;
53
54 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.mus_p4().at(midx)) < 0.4) {
55 jetIsLep = true;
56 break;
57 }
58 }
59
60 if (jetIsLep)
61 continue;
62
63 goodJets_.push_back(jidx);
64 }
65
66 return goodJets_;
67 }
68
69 std::vector<unsigned int> efftools::getRecoJets (struct GenParticleStruct& struct1, struct GenParticleStruct& struct2) {
70
71 const float mupt = 20.;
72 const float elpt = 20.;
73
74 std::vector<unsigned int> goodJets_;
75 for (unsigned int jidx = 0; jidx < cms2.pfjets_p4().size(); jidx++) {
76
77 float jet_eta = cms2.pfjets_p4().at(jidx).eta();
78 float jet_pt = cms2.pfjets_p4().at(jidx).pt() * cms2.pfjets_corL1FastL2L3().at(jidx);
79
80 if (fabs(jet_eta) > 2.4)
81 continue;
82 if (jet_pt < 40.)
83 continue;
84
85 std::pair<LorentzVector, int> recoLep1 = getRecoLepton(struct1);
86 if (recoLep1.second > -1) {
87 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), recoLep1.first) < 0.4)
88 continue;
89 }
90
91 std::pair<LorentzVector, int> recoLep2 = getRecoLepton(struct2);
92 if (recoLep2.second > -1) {
93 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), recoLep2.first) < 0.4)
94 continue;
95 }
96
97 bool jetIsLep = false;
98 for (unsigned int eidx = 0; eidx < cms2.els_p4().size(); eidx++) {
99
100 if (cms2.els_p4().at(eidx).pt() < elpt)
101 continue;
102 if (fabs(cms2.els_p4().at(eidx).eta()) > 2.4)
103 continue;
104 if (!samesign::isNumeratorLepton(11, eidx))
105 continue;
106
107 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.els_p4().at(eidx)) < 0.4) {
108 jetIsLep = true;
109 break;
110 }
111 }
112
113 if (jetIsLep)
114 continue;
115
116 for (unsigned int midx = 0; midx < cms2.mus_p4().size(); midx++) {
117
118 if (cms2.mus_p4().at(midx).pt() < mupt)
119 continue;
120 if (fabs(cms2.mus_p4().at(midx).eta()) > 2.4)
121 continue;
122 if (!samesign::isNumeratorLepton(13, midx))
123 continue;
124
125 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.mus_p4().at(midx)) < 0.4) {
126 jetIsLep = true;
127 break;
128 }
129 }
130
131 if (jetIsLep)
132 continue;
133
134 goodJets_.push_back(jidx);
135 }
136
137 return goodJets_;
138 }
139
140 std::vector<unsigned int> efftools::getRecoJets (unsigned int hidx) {
141
142 const float mupt = 20.;
143 const float elpt = 20.;
144
145 std::vector<unsigned int> goodJets_;
146 for (unsigned int jidx = 0; jidx < cms2.pfjets_p4().size(); jidx++) {
147
148 float jet_eta = cms2.pfjets_p4().at(jidx).eta();
149 float jet_pt = cms2.pfjets_p4().at(jidx).pt() * cms2.pfjets_corL1FastL2L3().at(jidx);
150
151 if (fabs(jet_eta) > 2.4)
152 continue;
153 if (jet_pt < 40.)
154 continue;
155
156 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.hyp_lt_p4().at(hidx)) < 0.4)
157 continue;
158 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.hyp_ll_p4().at(hidx)) < 0.4)
159 continue;
160
161 bool jetIsLep = false;
162 for (unsigned int eidx = 0; eidx < cms2.els_p4().size(); eidx++) {
163
164 if (cms2.els_p4().at(eidx).pt() < elpt)
165 continue;
166 if (fabs(cms2.els_p4().at(eidx).eta()) > 2.4)
167 continue;
168 if (!samesign::isNumeratorLepton(11, eidx))
169 continue;
170
171 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.els_p4().at(eidx)) < 0.4) {
172 jetIsLep = true;
173 break;
174 }
175 }
176
177 if (jetIsLep)
178 continue;
179
180 for (unsigned int midx = 0; midx < cms2.mus_p4().size(); midx++) {
181
182 if (cms2.mus_p4().at(midx).pt() < mupt)
183 continue;
184 if (fabs(cms2.mus_p4().at(midx).eta()) > 2.4)
185 continue;
186 if (!samesign::isNumeratorLepton(13, midx))
187 continue;
188
189 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), cms2.mus_p4().at(midx)) < 0.4) {
190 jetIsLep = true;
191 break;
192 }
193 }
194
195 if (jetIsLep)
196 continue;
197
198 goodJets_.push_back(jidx);
199 }
200
201 return goodJets_;
202 }
203
204 float efftools::getHT (std::vector<unsigned int>& jets) {
205
206 float ht_ = 0;
207 for (unsigned int jidx = 0; jidx < jets.size(); jidx++)
208 ht_ += cms2.pfjets_p4().at(jets.at(jidx)).pt() * cms2.pfjets_corL1FastL2L3().at(jets.at(jidx));
209
210 return ht_;
211 }
212
213 std::pair<LorentzVector, int> efftools::getRecoLepton (struct GenParticleStruct& struct1) {
214
215 const float mupt = 20.;
216 const float elpt = 20.;
217
218 int id = abs(struct1.id_) == 15 ? struct1.did_ : struct1.id_;
219 int idx = abs(struct1.id_) == 15 ? struct1.didx_ : struct1.idx_;
220
221 if (abs(id) == 11) {
222 for (unsigned int eidx = 0; eidx < cms2.els_p4().size(); eidx++) {
223
224 if (fabs(cms2.els_p4().at(eidx).eta()) > 2.4)
225 continue;
226 if (cms2.els_p4().at(eidx).pt() < elpt)
227 continue;
228 if (leptonIsFromW(eidx, 11, true) < 1)
229 continue;
230 if (mc3idx_eormu(11, eidx, 0.3, elpt) == idx)
231 return std::make_pair(cms2.els_p4().at(eidx), eidx);
232 }
233 }
234
235 if (abs(id) == 13) {
236 for (unsigned int midx = 0; midx < cms2.mus_p4().size(); midx++) {
237
238 if (fabs(cms2.mus_p4().at(midx).eta()) > 2.4)
239 continue;
240 if (cms2.mus_p4().at(midx).pt() < mupt)
241 continue;
242 if (leptonIsFromW(midx, 13, true) < 1)
243 continue;
244 if (mc3idx_eormu(13, midx, 0.3, mupt) == idx)
245 return std::make_pair(cms2.mus_p4().at(midx), midx);
246 }
247 }
248
249 return std::make_pair(LorentzVector(0., 0., 0., 0.), -1);
250 }
251
252 std::vector<std::pair<GenParticleStruct, GenParticleStruct> > efftools::makeGenHyps (float eta_cut, bool removeLeptonsOverlappingWithPartons) {
253
254 //
255 // loop over gen particles
256 //
257 std::vector<GenParticleStruct> gleps;
258 for (unsigned int gidx = 0; gidx < cms2.genps_p4().size(); gidx++) {
259
260 if (cms2.genps_status().at(gidx) != 3)
261 continue;
262
263 int id = cms2.genps_id().at(gidx);
264 float pt = cms2.genps_p4().at(gidx).pt();
265 float eta = cms2.genps_p4().at(gidx).eta();
266
267 if (fabs(eta) > eta_cut)
268 continue;
269
270 if (abs(id) == 11 || abs(id) == 13) {
271 if (removeLeptonsOverlappingWithPartons && leptonOverlapsWithParton(cms2.genps_p4().at(gidx)))
272 continue;
273
274 gleps.push_back(GenParticleStruct(id, gidx, pt, eta, -999, 999, -999., -999.));
275 }
276
277 //
278 // need to add protection here to not get more than one lepton from the same tau
279 //
280
281 if (abs(id) == 15) {
282
283 GenParticleStruct tmp_gp = GenParticleStruct(0, 0, 0, 0, 0, 0, 0, 0);
284 for (unsigned int didx = 0; didx < cms2.genps_lepdaughter_id().at(gidx).size(); didx++) {
285
286 int did = cms2.genps_lepdaughter_id().at(gidx).at(didx);
287 if (abs(did) != 11 && abs(did) != 13)
288 continue;
289
290 if (fabs(cms2.genps_lepdaughter_p4().at(gidx).at(didx).eta()) > eta_cut)
291 continue;
292
293 if (removeLeptonsOverlappingWithPartons && leptonOverlapsWithParton(cms2.genps_lepdaughter_p4().at(gidx).at(didx)))
294 continue;
295
296 float dpt = cms2.genps_lepdaughter_p4().at(gidx).at(didx).pt();
297 float deta = cms2.genps_lepdaughter_p4().at(gidx).at(didx).eta();
298 if (dpt > tmp_gp.dpt_)
299 tmp_gp = GenParticleStruct(id, gidx, pt, eta, did, didx, dpt, deta);
300 }
301
302 if (tmp_gp.dpt_ > 0.1)
303 gleps.push_back(tmp_gp);
304
305 } // end tau block
306
307 } // end loop over gen particles
308
309
310 //
311 // now make gen hyps
312 //
313 std::vector<std::pair<GenParticleStruct, GenParticleStruct> > glepPairs;
314 for (unsigned int idx1 = 0; idx1 < gleps.size(); idx1++) {
315 for (unsigned int idx2 = idx1 + 1; idx2 < gleps.size(); idx2++) {
316 glepPairs.push_back(make_pair(gleps.at(idx1), gleps.at(idx2)));
317 }
318 } // loop to make gen pairs
319
320 return glepPairs;
321 }
322
323 std::pair<GenParticleStruct, GenParticleStruct> efftools::getGenHyp (float pt1_cut, float pt2_cut, enum HypType hypType) {
324
325 float min_pt = std::min(pt1_cut, pt2_cut);
326 float max_pt = std::max(pt1_cut, pt2_cut);
327
328 std::vector<std::pair<GenParticleStruct, GenParticleStruct> > glepPairs = makeGenHyps();
329 unsigned int npairs = glepPairs.size();
330
331 GenParticleStruct gp = GenParticleStruct(0, -1, 0., 0., 0, -1, 0., 0.);
332 std::pair<GenParticleStruct, GenParticleStruct> good_gen_hyp = std::make_pair(gp, gp);
333 enum DileptonHypType good_hyp_type = DILEPTON_ALL;
334 for (unsigned int idx = 0; idx < npairs; idx++) {
335
336 GenParticleStruct gp1 = glepPairs.at(idx).first;
337 GenParticleStruct gp2 = glepPairs.at(idx).second;
338
339 float gpt1 = (abs(gp1.id_) == 15) ? gp1.dpt_ : gp1.pt_;
340 int gid1 = (abs(gp1.id_) == 15) ? gp1.did_ : gp1.id_;
341 float gpt2 = (abs(gp2.id_) == 15) ? gp2.dpt_ : gp2.pt_;
342 int gid2 = (abs(gp2.id_) == 15) ? gp2.did_ : gp2.id_;
343
344 // require SS
345 if (gid1 * gid2 < 0 && hypType == DILEPTON_SS)
346 continue;
347 else if (gid1 * gid2 > 0 && hypType == DILEPTON_OS)
348 continue;
349
350 if (std::min(gpt1, gpt2) < min_pt)
351 continue;
352
353 if (std::max(gpt1, gpt2) < max_pt)
354 continue;
355
356 enum DileptonHypType tmp_type = getHypType(gid1, gid2);
357 if (tmp_type == DILEPTON_ALL)
358 continue;
359
360 if (good_hyp_type == DILEPTON_ALL) {
361 good_hyp_type = tmp_type;
362 good_gen_hyp = glepPairs.at(idx);
363 }
364 else if (tmp_type < good_hyp_type) {
365 good_hyp_type = tmp_type;
366 good_gen_hyp = glepPairs.at(idx);
367 }
368 else if (tmp_type == good_hyp_type) {
369
370 float ggpt1 = (abs(good_gen_hyp.first.id_) == 15) ? good_gen_hyp.first.dpt_ : good_gen_hyp.first.pt_;
371 float ggpt2 = (abs(good_gen_hyp.second.id_) == 15) ? good_gen_hyp.second.dpt_ : good_gen_hyp.second.pt_;
372
373 if ( (gpt1+gpt2) > (ggpt1+ggpt2) ) {
374 good_hyp_type = tmp_type;
375 good_gen_hyp = glepPairs.at(idx);
376 }
377 }
378 } // end loop over gen hyps
379
380 return good_gen_hyp;
381 }
382
383 int efftools::getRecoJet (LorentzVector p4) {
384
385 float tmp_dr = 99.;
386 int b_index = -1;
387 for (unsigned int jidx = 0; jidx < cms2.pfjets_p4().size(); jidx++) {
388
389 if (fabs(cms2.pfjets_p4().at(jidx).eta()) > 2.4)
390 continue;
391
392 float jet_pt = cms2.pfjets_p4().at(jidx).pt() * cms2.pfjets_corL1FastL2L3().at(jidx);
393 if (jet_pt < 40.)
394 continue;
395 if (!passesPFJetID(jidx))
396 continue;
397
398 if (ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), p4) < tmp_dr) {
399 tmp_dr = ROOT::Math::VectorUtil::DeltaR(cms2.pfjets_p4().at(jidx), p4);
400 b_index = jidx;
401 }
402 }
403
404 return b_index;
405 }
406
407 std::vector<LorentzVector> efftools::getGenBjets (float pt_cut, float eta_cut) {
408
409 std::vector<LorentzVector> tmp_vec;
410 tmp_vec.clear();
411 for (unsigned int gidx = 0; gidx < cms2.genps_status().size(); gidx++) {
412 if (cms2.genps_status().at(gidx) != 3)
413 continue;
414 if (abs(cms2.genps_id().at(gidx)) != 5)
415 continue;
416 if (fabs(cms2.genps_p4().at(gidx).eta()) > eta_cut)
417 continue;
418 if (cms2.genps_p4().at(gidx).pt() < pt_cut)
419 continue;
420
421 tmp_vec.push_back(cms2.genps_p4().at(gidx));
422 }
423
424 return tmp_vec;
425 }
426
427 std::vector<LorentzVector> efftools::getGenJets (float pt_cut, float eta_cut) {
428
429 std::vector<LorentzVector> tmp_vec;
430 tmp_vec.clear();
431 for (unsigned int gidx = 0; gidx < cms2.genps_status().size(); gidx++) {
432 if (cms2.genps_status().at(gidx) != 3)
433 continue;
434 if ((abs(cms2.genps_id().at(gidx)) < 1 || abs(cms2.genps_id().at(gidx)) > 5) && abs(cms2.genps_id().at(gidx)) != 21)
435 continue;
436 if (fabs(cms2.genps_p4().at(gidx).eta()) > eta_cut)
437 continue;
438 if (cms2.genps_p4().at(gidx).pt() < pt_cut)
439 continue;
440
441 tmp_vec.push_back(cms2.genps_p4().at(gidx));
442 }
443
444 return tmp_vec;
445 }
446
447 float efftools::getGenHT(float pt_cut, float eta_cut) {
448
449 std::vector<LorentzVector> tmp_vec = getGenJets(pt_cut, eta_cut);
450
451 float tmp_ht = 0.;
452 for(unsigned int idx = 0; idx < tmp_vec.size(); idx++) {
453 tmp_ht += tmp_vec.at(idx).pt();
454 }
455
456 return tmp_ht;
457 }
458
459 bool efftools::leptonOverlapsWithParton(LorentzVector p4, float parton_pt, float dr) {
460
461 for (unsigned int idx = 0; idx < cms2.genps_p4().size(); idx++) {
462
463 if (cms2.genps_status().at(idx) != 3)
464 continue;
465 if ((abs(cms2.genps_id().at(idx)) < 1 || abs(cms2.genps_id().at(idx)) > 5) && abs(cms2.genps_id().at(idx)) != 21) // require a quark or gluon
466 continue;
467 if (cms2.genps_p4().at(idx).pt() < parton_pt)
468 continue;
469 if (ROOT::Math::VectorUtil::DeltaR(cms2.genps_p4().at(idx), p4) < dr)
470 return true;
471 }
472
473 return false;
474 }
475
476 enum DileptonHypType efftools::getHypType (int id1, int id2) {
477
478 if (abs(id1) != 11 && abs(id1) != 13)
479 return DILEPTON_ALL;
480 if (abs(id2) != 11 && abs(id2) != 13)
481 return DILEPTON_ALL;
482
483 if (abs(id1) == 11 && abs(id2) == 11)
484 return DILEPTON_EE;
485 else if (abs(id1) == 13 && abs(id2) == 13)
486 return DILEPTON_MUMU;
487 else
488 return DILEPTON_EMU;
489 }
490
491
492 int efftools::getGenParton (LorentzVector p4) {
493
494 float tmp_dr = 99.;
495 int index = -1;
496 for (unsigned int jidx = 0; jidx < cms2.genps_p4().size(); jidx++) {
497
498 if (fabs(cms2.genps_p4().at(jidx).eta()) > 2.4)
499 continue;
500 if (cms2.genps_p4().at(jidx).pt() < 40.)
501 continue;
502
503 if (ROOT::Math::VectorUtil::DeltaR(cms2.genps_p4().at(jidx), p4) < tmp_dr) {
504 tmp_dr = ROOT::Math::VectorUtil::DeltaR(cms2.genps_p4().at(jidx), p4);
505 index = jidx;
506 }
507 }
508
509 return index;
510 }