ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/UHHAnalysis/SFrameAnalysis/src/Cleaner.cxx
Revision: 1.33
Committed: Tue Jun 25 16:42:03 2013 UTC (11 years, 10 months ago) by rkogler
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.32: +6 -0 lines
Log Message:
update for new data format and minor fixes

File Contents

# User Rev Content
1 peiffer 1.1 #include "../include/Cleaner.h"
2    
3 bazterra 1.23 Cleaner::Cleaner( BaseCycleContainer* input)
4     {
5 peiffer 1.1
6 bazterra 1.23 bcc = input;
7 rkogler 1.28 m_jec_unc = NULL;
8     m_jecvar = e_Default;
9     m_jervar = e_Default;
10 peiffer 1.1
11     }
12    
13 bazterra 1.23 Cleaner::Cleaner()
14     {
15 peiffer 1.9
16 peiffer 1.32 EventCalc* calc = EventCalc::Instance();
17     bcc = calc->GetBaseCycleContainer();
18 rkogler 1.28 m_jec_unc = NULL;
19     m_jecvar = e_Default;
20     m_jervar = e_Default;
21 peiffer 1.9
22     }
23    
24 bazterra 1.23 void Cleaner::resetEventCalc()
25     {
26 peiffer 1.9
27 bazterra 1.23 //reset everything in EventCalc except the event weight
28     EventCalc* evc = EventCalc::Instance();
29     double weight = evc->GetWeight();
30     evc->Reset();
31     evc->ProduceWeight(weight);
32 peiffer 1.9
33     }
34    
35 rkogler 1.28 void Cleaner::JetEnergyResolutionShifter(bool sort)
36 bazterra 1.23 {
37 peiffer 1.2
38 bazterra 1.23 LorentzVector met(0,0,0,0);
39     if(bcc->met) {
40     met = bcc->met->v4();
41     }
42    
43     for(unsigned int i=0; i<bcc->jets->size(); ++i) {
44 peiffer 1.2
45 rkogler 1.28 //std::cout << "ResolutionShifter start: Jet " << i << ", pt = " << bcc->jets->at(i).pt() << " correction factor = " << 1./bcc->jets->at(i).JEC_factor_raw() << std::endl;
46    
47 bazterra 1.23 float genpt = bcc->jets->at(i).genjet_pt();
48     //ignore unmatched jets (which have zero vector) or jets with very low pt:
49     if(genpt < 15.0) {
50     // std::cout << "1.0 | " << bcc->jets->at(i).pt() << " | " << bcc->jets->at(i).eta() << " | " << genpt << std::endl;
51     continue;
52     }
53 peiffer 1.15
54 bazterra 1.23 LorentzVector jet_v4 = bcc->jets->at(i).v4();
55    
56     LorentzVector jet_v4_raw = jet_v4*bcc->jets->at(i).JEC_factor_raw();
57 peiffer 1.16
58 rkogler 1.28 double recopt = jet_v4.pt();
59     double factor = 0.0;
60     double uperr = 0.0;
61     double downerr = 0.0;
62     double abseta = fabs(jet_v4.eta());
63 peiffer 1.2
64 bazterra 1.23 //numbers taken from https://twiki.cern.ch/twiki/bin/view/CMS/JetResolution
65 rkogler 1.28 if(m_jervar==e_Default) {
66 bazterra 1.23 if(abseta < 0.5)
67     factor = 0.052;
68     else if(abseta >= 0.5 && abseta <1.1)
69     factor = 0.057;
70     else if(abseta >= 1.1 && abseta <1.7)
71     factor = 0.096;
72     else if(abseta >= 1.7 && abseta <2.3)
73     factor = 0.134;
74     else if(abseta >= 2.3)
75     factor = 0.288;
76     else
77     factor = 0.288;
78 rkogler 1.28 } else if(m_jervar==e_Up) {
79 bazterra 1.23 if(abseta < 0.5)
80     factor = 0.115;
81     else if(abseta >= 0.5 && abseta <1.1)
82     factor = 0.114;
83     else if(abseta >= 1.1 && abseta <1.7)
84     factor = 0.161;
85     else if(abseta >= 1.7 && abseta <2.3)
86     factor = 0.228;
87     else if(abseta >= 2.3)
88     factor = 0.488;
89 rkogler 1.28 } else if(m_jervar==e_Down) {
90 bazterra 1.23 if(abseta < 0.5)
91     factor = -0.01;
92     else if(abseta >= 0.5 && abseta <1.1)
93     factor = 0.00;
94     else if(abseta >= 1.1 && abseta <1.7)
95     factor = 0.032;
96     else if(abseta >= 1.7 && abseta <2.3)
97     factor = 0.042;
98     else if(abseta >= 2.3)
99     factor = 0.089;
100     }
101 peiffer 1.2
102 rkogler 1.28 double ptscale = std::max(0.0, 1 + factor * (recopt - genpt) / recopt);
103 peiffer 1.16
104 bazterra 1.23 jet_v4*=ptscale;
105 peiffer 1.16
106 bazterra 1.23 bcc->jets->at(i).set_v4(jet_v4);
107 rkogler 1.28 //std::cout << "ResolutionShifter: Jet " << i << ", pt = " << jet_v4.pt() << std::endl;
108 bazterra 1.23 //propagate JER shifts to MET
109     //if(jet_v4.Pt()>25) {
110     met += jet_v4_raw;
111     jet_v4_raw*=ptscale;
112     met -= jet_v4_raw;
113     //}
114 peiffer 1.16
115 bazterra 1.23 //std::cout << ptscale << " | " << jet_v4.Pt() << " | "<< jet_v4.eta() << " | " << genpt << std::endl;
116 peiffer 1.16
117 rkogler 1.28 //std::cout << "ResolutionShifter end: Jet " << i << ", pt = " << bcc->jets->at(i).pt() << " correction factor = " << 1./bcc->jets->at(i).JEC_factor_raw() << std::endl;
118    
119 bazterra 1.23 }
120 peiffer 1.6
121 bazterra 1.23 //store changed MET
122     if(bcc->met) {
123     bcc->met->set_pt(met.Pt());
124     bcc->met->set_phi(met.Phi());
125     }
126    
127     if(sort) std::sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
128     resetEventCalc();
129 peiffer 1.2 }
130    
131 peiffer 1.1
132 bazterra 1.23 void Cleaner::JetLeptonSubtractor(FactorizedJetCorrector *corrector, bool sort)
133     {
134    
135     //std::cout<< "event: " <<bcc->event <<std::endl;
136     //std::cout<< "ID| pt raw | PAT pt | area | rho | corr | old cor| JER ptscale | pt new | eta new | ptgen " <<std::endl;
137    
138 rkogler 1.33 if (!corrector){
139     std::cerr << "JetLeptonSubtractor: called without a valid JetCorrector! Please correct the error." << std::endl;
140     std::cerr << "Hint: have you supplied correct JEC files in your xml steering?" << std::endl;
141     exit(EXIT_FAILURE);
142     }
143    
144 peiffer 1.7
145 bazterra 1.23 for(unsigned int i=0; i<bcc->jets->size(); ++i) {
146 rkogler 1.28
147 bazterra 1.23 LorentzVector jet_v4_raw = bcc->jets->at(i).v4()*bcc->jets->at(i).JEC_factor_raw();
148 peiffer 1.7
149 bazterra 1.23 //subtract lepton momenta from raw jet momentum
150 peiffer 1.7
151 bazterra 1.23 //std::cout << i+1 << " | " << jet_v4_raw.Pt() << std::endl;
152     double ele_energy = bcc->jets->at(i).chargedEmEnergyFraction()*jet_v4_raw.E();
153     double mu_energy = bcc->jets->at(i).muonEnergyFraction()*jet_v4_raw.E();
154 peiffer 1.11
155 bazterra 1.23 if(bcc->electrons) {
156     for(unsigned int j=0; j<bcc->electrons->size(); ++j) {
157     if(bcc->jets->at(i).deltaR(bcc->electrons->at(j))<0.5) {
158 peiffer 1.13 // if(jet_v4_raw.pt() >= bcc->electrons->at(j).pt()){
159 bazterra 1.23 jet_v4_raw -= bcc->electrons->at(j).v4();
160 peiffer 1.13 // bcc->jets->at(i).set_electronMultiplicity(bcc->jets->at(i).electronMultiplicity()-1);
161     // ele_energy -= bcc->electrons->at(j).energy();
162 bazterra 1.23 //}
163     //else{
164     //std::cout << bcc->jets->at(i).pt()<< " -ele pt:"<< bcc->electrons->at(j).pt() << std::endl;
165     //jet_v4_raw -= jet_v4_raw;
166     //}
167     }
168     }
169     }
170     if(bcc->muons) {
171     for(unsigned int j=0; j<bcc->muons->size(); ++j) {
172     if(bcc->jets->at(i).deltaR(bcc->muons->at(j))<0.5) {
173     //if(jet_v4_raw.pt()>= bcc->muons->at(j).pt()){
174     jet_v4_raw -= bcc->muons->at(j).v4();
175 peiffer 1.13 // bcc->jets->at(i).set_muonMultiplicity(bcc->jets->at(i).muonMultiplicity()-1);
176     // mu_energy -= bcc->muons->at(j).energy();
177 bazterra 1.23 //}
178     //else{
179     //std::cout << bcc->jets->at(i).pt() << " -mu pt:"<< bcc->muons->at(j).pt() << std::endl;
180     //jet_v4_raw -= jet_v4_raw;
181     //}
182     }
183     }
184     }
185     if(ele_energy<=jet_v4_raw.E())
186     bcc->jets->at(i).set_chargedEmEnergyFraction(ele_energy/jet_v4_raw.E());
187     if(mu_energy<=jet_v4_raw.E())
188     bcc->jets->at(i).set_muonEnergyFraction(mu_energy/jet_v4_raw.E());
189    
190     //apply jet energy corrections to modified raw momentum
191     corrector->setJetPt(jet_v4_raw.Pt());
192     corrector->setJetEta(jet_v4_raw.Eta());
193     corrector->setJetE(jet_v4_raw.E());
194     corrector->setJetA(bcc->jets->at(i).jetArea());
195     corrector->setRho(bcc->rho);
196     corrector->setNPV(bcc->pvs->size());
197    
198 rkogler 1.28 double correctionfactor = corrector->getCorrection();
199 bazterra 1.23
200     LorentzVector jet_v4_corrected = jet_v4_raw *correctionfactor;
201    
202 rkogler 1.28 if (m_jecvar != e_Default){
203     if (m_jec_unc==NULL){
204     std::cerr << "JEC variation should be applied, but JEC uncertainty object is NULL! Abort." << std::endl;
205     exit(EXIT_FAILURE);
206     }
207 rkogler 1.29 // ignore jets with very low pt or high eta, avoiding a crash from the JESUncertainty tool
208     double pt = jet_v4_corrected.Pt();
209     double eta = jet_v4_corrected.Eta();
210     if (pt<5. || fabs(eta)>5.) continue;
211    
212     m_jec_unc->setJetEta(eta);
213     m_jec_unc->setJetPt(pt);
214    
215 rkogler 1.28 double unc = 0.;
216     if (m_jecvar == e_Up){
217     unc = m_jec_unc->getUncertainty(1);
218     correctionfactor *= (1 + fabs(unc));
219     } else if (m_jecvar == e_Down){
220     unc = m_jec_unc->getUncertainty(-1);
221     correctionfactor *= (1 - fabs(unc));
222     }
223     jet_v4_corrected = jet_v4_raw * correctionfactor;
224     }
225    
226 bazterra 1.23 //std::cout << i+1 << " | " << jet_v4_raw.Pt() << " | " << bcc->jets->at(i).pt() << " | " << bcc->jets->at(i).jetArea()<< " | " << bcc->rho << " | " << correctionfactor << " | " << 1./bcc->jets->at(i).JEC_factor_raw() << " | " << std::endl;
227    
228     bcc->jets->at(i).set_v4(jet_v4_corrected);
229     bcc->jets->at(i).set_JEC_factor_raw(1./correctionfactor);
230     }
231    
232     if(sort) std::sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
233     resetEventCalc();
234     }
235    
236     void Cleaner::JetRecorrector(FactorizedJetCorrector *corrector, bool sort)
237     {
238    
239     for(unsigned int i=0; i<bcc->jets->size(); ++i) {
240    
241     LorentzVector jet_v4_raw = bcc->jets->at(i).v4()*bcc->jets->at(i).JEC_factor_raw();
242     corrector->setJetPt(jet_v4_raw.Pt());
243     corrector->setJetEta(jet_v4_raw.Eta());
244     corrector->setJetE(jet_v4_raw.E());
245     corrector->setJetA(bcc->jets->at(i).jetArea());
246     corrector->setRho(bcc->rho);
247     corrector->setNPV(bcc->pvs->size());
248    
249     float correctionfactor = corrector->getCorrection();
250    
251     LorentzVector jet_v4_corrected = jet_v4_raw *correctionfactor;
252    
253 rkogler 1.28 if (m_jecvar != e_Default){
254     if (m_jec_unc==NULL){
255     std::cerr << "JEC variation should be applied, but JEC uncertainty object is NULL! Abort." << std::endl;
256     exit(EXIT_FAILURE);
257     }
258     m_jec_unc->setJetEta(jet_v4_corrected.Eta());
259     m_jec_unc->setJetPt(jet_v4_corrected.Pt());
260     double unc = 0.;
261     if (m_jecvar == e_Up){
262     unc = m_jec_unc->getUncertainty(1);
263     correctionfactor *= (1 + fabs(unc));
264     } else if (m_jecvar == e_Down){
265     unc = m_jec_unc->getUncertainty(-1);
266     correctionfactor *= (1 - fabs(unc));
267     }
268     jet_v4_corrected = jet_v4_raw * correctionfactor;
269     }
270    
271    
272 bazterra 1.23 bcc->jets->at(i).set_v4(jet_v4_corrected);
273     bcc->jets->at(i).set_JEC_factor_raw(1./correctionfactor);
274     }
275    
276     if(sort) std::sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
277     resetEventCalc();
278     }
279    
280 bazterra 1.27 // Old electron id
281 peiffer 1.30 /*
282     bool Cleaner::passElectronId(BaseCycleContainer * bcc, unsigned int index)
283 bazterra 1.26 {
284     Electron ele = bcc->electrons->at(index);
285     if(fabs(ele.supercluster_eta())<1.4442 || fabs(ele.supercluster_eta())>1.5660) {
286     if(bcc->pvs->size()>0) {
287     if(fabs(ele.gsfTrack_dxy_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y()))<0.02) {
288     if(fabs(ele.gsfTrack_dz_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y(), bcc->pvs->at(0).z()))<0.1) {
289     if(ele.passconversionveto()) {
290     if(ele.mvaTrigV0()>0.0) {
291     if(ele.eleID(Electron::e_Tight)) {
292     return true;
293     }
294     }
295     }
296     }
297     }
298     }
299     }
300     return false;
301 bazterra 1.27 }*/
302    
303     // Base on the new top simple cut recomendation for 53x
304 peiffer 1.30 /*
305     bool Cleaner::passElectronId(BaseCycleContainer * bcc, unsigned int index)
306 bazterra 1.27 {
307     Electron ele = bcc->electrons->at(index);
308     if(fabs(ele.supercluster_eta())<1.4442 || fabs(ele.supercluster_eta())>1.5660) {
309     if(bcc->pvs->size()>0) {
310     if(fabs(ele.gsfTrack_dxy_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y()))<0.02) {
311     if(fabs(ele.gsfTrack_dz_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y(), bcc->pvs->at(0).z()))<0.1) {
312 peiffer 1.30 if(ele.passconversionveto()) {
313     if(ele.eleID(Electron::e_Tight)) {
314     return true;
315     }
316     }
317 bazterra 1.27 }
318     }
319     }
320     }
321     return false;
322 peiffer 1.30 }
323     */
324     // Base on the new top mva recomendation for 53x
325    
326 bazterra 1.27
327     bool Cleaner::passElectronId(BaseCycleContainer * bcc, unsigned int index)
328     {
329     Electron ele = bcc->electrons->at(index);
330     if(fabs(ele.supercluster_eta())<1.4442 || fabs(ele.supercluster_eta())>1.5660) {
331     if(bcc->pvs->size()>0) {
332     if(fabs(ele.gsfTrack_dxy_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y()))<0.02) {
333     if(ele.passconversionveto()) {
334     if(ele.mvaTrigV0()>0.9) {
335     if (ele.gsfTrack_trackerExpectedHitsInner_numberOfLostHits() <= 0) {
336     return true;
337     }
338     }
339     }
340     }
341     }
342     }
343     return false;
344 bazterra 1.26 }
345    
346 peiffer 1.30
347 bazterra 1.26 void Cleaner::ElectronCleaner_noID_noIso(double ptmin, double etamax)
348 bazterra 1.23 {
349     std::vector<Electron> good_eles;
350     for(unsigned int i=0; i<bcc->electrons->size(); ++i) {
351     Electron ele = bcc->electrons->at(i);
352     if(ele.pt()>ptmin) {
353     if(fabs(ele.eta())<etamax) {
354     good_eles.push_back(ele);
355     }
356     }
357     }
358    
359     bcc->electrons->clear();
360    
361     for(unsigned int i=0; i<good_eles.size(); ++i) {
362     bcc->electrons->push_back(good_eles[i]);
363     }
364     sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
365     resetEventCalc();
366     }
367    
368 bazterra 1.26 void Cleaner::ElectronCleaner_noIso(double ptmin, double etamax, bool reverseID)
369     {
370     ElectronCleaner_noID_noIso(ptmin, etamax);
371     std::vector<Electron> good_eles;
372     for(unsigned int i=0; i<bcc->electrons->size(); ++i) {
373     Electron ele = bcc->electrons->at(i);
374     if(!reverseID && passElectronId(bcc, i))
375     good_eles.push_back(ele);
376     else if (reverseID && !passElectronId(bcc, i))
377     good_eles.push_back(ele);
378     }
379     bcc->electrons->clear();
380    
381     for(unsigned int i=0; i<good_eles.size(); ++i) {
382     bcc->electrons->push_back(good_eles[i]);
383     }
384     sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
385     resetEventCalc();
386     }
387    
388    
389     void Cleaner::ElectronCleaner(double ptmin, double etamax, double relisomax, bool reverseID, bool reverseIso)
390 bazterra 1.23 {
391 bazterra 1.26 ElectronCleaner_noIso(ptmin, etamax, reverseID);
392 bazterra 1.23 std::vector<Electron> good_eles;
393     for(unsigned int i=0; i<bcc->electrons->size(); ++i) {
394     Electron ele = bcc->electrons->at(i);
395 bazterra 1.26 if(!reverseIso && ele.relIsorho(bcc->rho)<relisomax)
396     good_eles.push_back(ele);
397     else if (reverseIso && ele.relIsorho(bcc->rho)>=relisomax)
398 bazterra 1.23 good_eles.push_back(ele);
399     }
400     bcc->electrons->clear();
401    
402     for(unsigned int i=0; i<good_eles.size(); ++i) {
403     bcc->electrons->push_back(good_eles[i]);
404     }
405     sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
406     resetEventCalc();
407     }
408    
409 bazterra 1.26
410 bazterra 1.23 void Cleaner::MuonCleaner_noID_noIso(double ptmin, double etamax)
411     {
412     std::vector<Muon> good_mus;
413     for(unsigned int i=0; i<bcc->muons->size(); ++i) {
414     Muon mu = bcc->muons->at(i);
415     if(mu.pt()>ptmin) {
416     if(fabs(mu.eta())<etamax) {
417     good_mus.push_back(mu);
418     }
419     }
420     }
421     bcc->muons->clear();
422    
423     for(unsigned int i=0; i<good_mus.size(); ++i) {
424     bcc->muons->push_back(good_mus[i]);
425     }
426     sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
427     resetEventCalc();
428     }
429    
430     void Cleaner::MuonCleaner_noIso(double ptmin, double etamax)
431     {
432    
433     MuonCleaner_noID_noIso(ptmin, etamax);
434     std::vector<Muon> good_mus;
435     for(unsigned int i=0; i<bcc->muons->size(); ++i) {
436     Muon mu = bcc->muons->at(i);
437     if(mu.isGlobalMuon()) {
438     if(mu.isPFMuon()) {
439     if(mu.globalTrack_chi2()/mu.globalTrack_ndof()<10) {
440     if(mu.globalTrack_numberOfValidMuonHits()>0) {
441     if(mu.innerTrack_trackerLayersWithMeasurement()>5) {
442 peiffer 1.24 if(mu.dB()<0.2) {
443 bazterra 1.23 if(fabs(mu.vertex_z()-bcc->pvs->at(0).z())<0.5) {
444     if(mu.innerTrack_numberOfValidPixelHits()>0) {
445     if(mu.numberOfMatchedStations()>1) {
446     good_mus.push_back(mu);
447     }
448     }
449     }
450     }
451     }
452     }
453     }
454     }
455     }
456     }
457     bcc->muons->clear();
458    
459     for(unsigned int i=0; i<good_mus.size(); ++i) {
460     bcc->muons->push_back(good_mus[i]);
461     }
462     sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
463     resetEventCalc();
464     }
465    
466     void Cleaner::MuonCleaner(double ptmin, double etamax, double relisomax)
467     {
468    
469     MuonCleaner_noIso(ptmin, etamax);
470     std::vector<Muon> good_mus;
471     for(unsigned int i=0; i<bcc->muons->size(); ++i) {
472     Muon mu = bcc->muons->at(i);
473     if(mu.relIso()<relisomax) {
474     good_mus.push_back(mu);
475     }
476     }
477     bcc->muons->clear();
478    
479     for(unsigned int i=0; i<good_mus.size(); ++i) {
480     bcc->muons->push_back(good_mus[i]);
481     }
482     sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
483     resetEventCalc();
484     }
485    
486     void Cleaner::TauCleaner(double ptmin, double etamax)
487     {
488    
489     std::vector<Tau> good_taus;
490     for(unsigned int i=0; i<bcc->taus->size(); ++i) {
491     Tau tau = bcc->taus->at(i);
492     if(tau.pt()>ptmin) {
493     if(fabs(tau.eta())<etamax) {
494     if(bcc->taus->at(i).decayModeFinding()) {
495     if(bcc->taus->at(i).byMediumCombinedIsolationDeltaBetaCorr()) {
496 peiffer 1.31 if(bcc->taus->at(i).againstElectronTightMVA3()) {
497     if(bcc->taus->at(i).againstMuonTight2()) {
498 bazterra 1.23 good_taus.push_back(tau);
499     }
500     }
501     }
502     }
503     }
504     }
505     }
506    
507     bcc->taus->clear();
508    
509     for(unsigned int i=0; i<good_taus.size(); ++i) {
510     bcc->taus->push_back(good_taus[i]);
511     }
512     sort(bcc->taus->begin(), bcc->taus->end(), HigherPt());
513     resetEventCalc();
514     }
515    
516    
517     void Cleaner::JetCleaner(double ptmin, double etamax, bool doPFID)
518     {
519    
520     std::vector<Jet> good_jets;
521     for(unsigned int i=0; i<bcc->jets->size(); ++i) {
522     Jet jet = bcc->jets->at(i);
523     if(jet.pt()>ptmin) {
524     if(fabs(jet.eta())<etamax) {
525     if(!doPFID || jet.pfID()) {
526     good_jets.push_back(jet);
527     }
528     }
529     }
530     }
531    
532     bcc->jets->clear();
533    
534     for(unsigned int i=0; i<good_jets.size(); ++i) {
535     bcc->jets->push_back(good_jets[i]);
536     }
537     sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
538     resetEventCalc();
539     }
540    
541     void Cleaner::TopJetCleaner(double ptmin, double etamax, bool doPFID)
542     {
543    
544     std::vector<TopJet> good_topjets;
545     for(unsigned int i=0; i<bcc->topjets->size(); ++i) {
546     TopJet topjet = bcc->topjets->at(i);
547     if(topjet.pt()>ptmin) {
548     if(fabs(topjet.eta())<etamax) {
549     if(!doPFID || topjet.pfID()) {
550     good_topjets.push_back(topjet);
551     }
552     }
553     }
554     }
555     bcc->topjets->clear();
556    
557     for(unsigned int i=0; i<good_topjets.size(); ++i) {
558     bcc->topjets->push_back(good_topjets[i]);
559     }
560     sort(bcc->topjets->begin(), bcc->topjets->end(), HigherPt());
561     resetEventCalc();
562     }
563    
564    
565     void Cleaner::PrimaryVertexCleaner(int ndofmax, double zmax, double rhomax)
566     {
567    
568     std::vector<PrimaryVertex> good_pvs;
569     for(unsigned int i=0; i<bcc->pvs->size(); ++i) {
570     PrimaryVertex pv = bcc->pvs->at(i);
571     if(pv.ndof()>=ndofmax && fabs(pv.z())<zmax && pv.rho()<=rhomax) {
572     good_pvs.push_back(pv);
573     }
574     }
575     bcc->pvs->clear();
576    
577     for(unsigned int i=0; i<good_pvs.size(); ++i) {
578     bcc->pvs->push_back(good_pvs[i]);
579     }
580     resetEventCalc();
581 peiffer 1.7 }
582 peiffer 1.30
583     void Cleaner::METPhiCorrector()
584     {
585    
586     float met_x = bcc->met->v4().x();
587     float met_y = bcc->met->v4().y();
588    
589     met_x += 0.2661+0.3217*bcc->pvs->size();
590     met_y += -0.2252-0.1747*bcc->pvs->size();
591    
592     float met = sqrt(met_x*met_x + met_y*met_y);
593     float phi = atan(met_y/met_x);
594    
595     bcc->met->set_phi(phi);
596     bcc->met->set_pt(met);
597    
598     resetEventCalc();
599     }