SFrameAnalysis/src/Cleaner.cxx

#include "../include/Cleaner.h"


Cleaner::Cleaner( BaseCycleContainer* input){

  bcc = input;

}

Cleaner::Cleaner(){

  ObjectHandler* objs = ObjectHandler::Instance();
  bcc = objs->GetBaseCycleContainer();

}

void Cleaner::resetEventCalc(){

  EventCalc* evc = EventCalc::Instance();
  evc->Reset();

}

void Cleaner::JetEnergyResolutionShifter(E_SystShift syst_shift){

  LorentzVector met(0,0,0,0);
  if(bcc->met) {
    met.SetPt(bcc->met->pt());
    met.SetPhi(bcc->met->phi());
  }

  for(unsigned int i=0; i<bcc->jets->size(); ++i){
    Jet jet = bcc->jets->at(i);
    float gen_pt = jet.genjet_pt();
    //ignore unmatched jets (which have zero vector) or jets with very low pt:
    if(gen_pt < 15.0) continue;

    met += jet.v4()*jet.JEC_factor_raw();
    
    float recopt = jet.pt();
    double factor = -10;
    double abseta = fabs(jet.eta());

    //numbers taken from https://twiki.cern.ch/twiki/bin/view/CMS/JetResolution
    if(syst_shift==e_Default){
      if(abseta < 0.5)
        factor = 0.052;
      else if(abseta >= 0.5 && abseta <1.1)
        factor = 0.057;
      else if(abseta >= 1.1 && abseta <1.7)
        factor = 0.096;
      else if(abseta >= 1.7 && abseta <2.3)
        factor = 0.134;
      else if(abseta >= 2.3)
        factor = 0.288;
    }
    else if(syst_shift==e_Up){
      if(abseta < 0.5)
        factor = 0.11;
      else if(abseta >= 0.5 && abseta <1.1)
        factor = 0.12;
      else if(abseta >= 1.1 && abseta <1.7)
        factor = 0.16;
      else if(abseta >= 1.7 && abseta <2.3)
        factor = 0.23;
      else if(abseta >= 2.3)
        factor = 0.49;
    }
    else if(syst_shift==e_Down){
      if(abseta < 0.5)
        factor = -0.01;
      else if(abseta >= 0.5 && abseta <1.1)
        factor = 0.00;
      else if(abseta >= 1.1 && abseta <1.7)
        factor = 0.04;
      else if(abseta >= 1.7 && abseta <2.3)
        factor = 0.03;
      else if(abseta >= 2.3)
        factor = 0.09;
    }

    float deltapt = (recopt - gen_pt) * factor;
    float ptscale = std::max(0.0f, (recopt + deltapt) / recopt);
    jet.set_pt(jet.pt() * ptscale);
    bcc->jets->at(i).set_pt( jet.pt());

    //propagate JER shifts to MET
    met -= jet.v4()*jet.JEC_factor_raw();
  }
  
  //store changed MET, flip phi if new MET is negative
  if(bcc->met){
    bcc->met->set_pt(met.Pt());
    bcc->met->set_phi(met.Phi());
  }

  sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
  resetEventCalc();
}


void Cleaner::JetLeptonSubtractor(FactorizedJetCorrector *corrector){

  for(unsigned int i=0; i<bcc->jets->size(); ++i){

    LorentzVector jet_v4_raw = bcc->jets->at(i).v4()*bcc->jets->at(i).JEC_factor_raw();

    //subtract lepton momenta from raw jet momentum
     
    double ele_energy =  bcc->jets->at(i).chargedEmEnergyFraction()*jet_v4_raw.E();
    double mu_energy = bcc->jets->at(i).muonEnergyFraction()*jet_v4_raw.E();

    if(bcc->electrons){
      for(unsigned int j=0; j<bcc->electrons->size(); ++j){
        if(bcc->jets->at(i).deltaR(bcc->electrons->at(j))<0.5){
//       if(jet_v4_raw.pt() >= bcc->electrons->at(j).pt()){
          jet_v4_raw -= bcc->electrons->at(j).v4();
//        bcc->jets->at(i).set_electronMultiplicity(bcc->jets->at(i).electronMultiplicity()-1);
//        ele_energy -= bcc->electrons->at(j).energy();
          //}
          //else{
          //std::cout << "Subtracting electron with pt="<< bcc->electrons->at(j).pt() <<" from jet with pt(corrected)= " << bcc->jets->at(i).pt() << ", pt(raw)=" << jet_v4_raw.pt()  << std::endl;
          //jet_v4_raw -= jet_v4_raw;
          //}
        }
      }
    }
    if(bcc->muons){
      for(unsigned int j=0; j<bcc->muons->size(); ++j){
        if(bcc->jets->at(i).deltaR(bcc->muons->at(j))<0.5){
          //if(jet_v4_raw.pt()>= bcc->muons->at(j).pt()){
          jet_v4_raw -= bcc->muons->at(j).v4();
//        bcc->jets->at(i).set_muonMultiplicity(bcc->jets->at(i).muonMultiplicity()-1);
//        mu_energy -= bcc->muons->at(j).energy();
          //}
          //else{
          //std::cout << "Subtracting muon with pt="<< bcc->muons->at(j).pt() <<" from jet with pt(corrected)= " << bcc->jets->at(i).pt() << ", pt(raw)=" << jet_v4_raw.pt()  << std::endl;
          //jet_v4_raw -= jet_v4_raw;   
          //}
        }
      }
    }
    if(ele_energy<=jet_v4_raw.E())
      bcc->jets->at(i).set_chargedEmEnergyFraction(ele_energy/jet_v4_raw.E());
    if(mu_energy<=jet_v4_raw.E())
     bcc->jets->at(i).set_muonEnergyFraction(mu_energy/jet_v4_raw.E());
    

    //apply jet energy corrections to modified raw momentum
    corrector->setJetPt(jet_v4_raw.Pt());
    corrector->setJetEta(jet_v4_raw.Eta());
    corrector->setJetE(jet_v4_raw.E());  
    corrector->setJetA(bcc->jets->at(i).jetArea());
    corrector->setRho(bcc->rho);
    corrector->setNPV(bcc->pvs->size());        

    float correctionfactor = corrector->getCorrection();
    
    LorentzVector jet_v4_corrected = jet_v4_raw *correctionfactor;

    bcc->jets->at(i).set_v4(jet_v4_corrected);
    bcc->jets->at(i).set_JEC_factor_raw(1./correctionfactor);

  }

  sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
  resetEventCalc();
}

//tight ele ID from https://twiki.cern.ch/twiki/bin/view/CMS/EgammaCutBasedIdentification
bool Cleaner::eleID(Electron ele){

  bool pass=false;
  
  if(fabs(ele.supercluster_eta())<1.4442){
    if(ele.dEtaIn()<0.004 && ele.dPhiIn()<0.03 && ele.sigmaIEtaIEta()<0.01 && ele.HoverE()<0.12 && fabs(1./ele.EcalEnergy()-1./ele.v4().P())<0.05) pass=true; 
  }
  else if( fabs(ele.supercluster_eta())>1.5660){
    if(ele.dEtaIn()<0.005 && ele.dPhiIn()<0.02 && ele.sigmaIEtaIEta()<0.03 && ele.HoverE()<0.10 && fabs(1./ele.EcalEnergy()-1./ele.v4().P())<0.05) pass=true; 
  }

  return pass;

}

//pf ID has already been applied when using goodPatJets
bool Cleaner::pfID(Jet jet){

  if(jet.numberOfDaughters()>1 
     && jet.neutralHadronEnergyFraction()<0.99
     && jet.neutralEmEnergyFraction()<0.99){
    
    if(fabs(jet.eta())>=2.4)
      return true;
    
    if(jet.chargedEmEnergyFraction()<0.99
       && jet.chargedHadronEnergyFraction()>0
       && jet.chargedMultiplicity()>0)
      return true;

  }
  std::cout << "Bloeder Jet, pt,eta= " << jet.pt() << ", " <<jet.eta() << ";  " << jet.numberOfDaughters() << " " << jet.neutralHadronEnergyFraction() << " " <<jet.neutralEmEnergyFraction() << " " << jet.chargedEmEnergyFraction() << " " << jet.chargedHadronEnergyFraction() << " " << jet.chargedMultiplicity() <<std::endl;
  return false;
     
}


void Cleaner::ElectronCleaner_noIso(double ptmin, double etamax){

  std::vector<Electron> good_eles;
  for(unsigned int i=0; i<bcc->electrons->size(); ++i){
    Electron ele = bcc->electrons->at(i);
    if(ele.pt()>ptmin){
      if(fabs(ele.eta())<etamax){
//      if(fabs(ele.supercluster_eta())<1.4442 || fabs(ele.supercluster_eta())>1.5660){
//        if(bcc->pvs->size()>0){
//          if(ele.gsfTrack_dxy_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y())<0.02){
//            if(ele.passconversionveto()){
//              if(ele.mvaTrigV0()>0.0){
//                if(eleID(ele)){
                    good_eles.push_back(ele);
                  }
//              }
//            }
//          }
//        }
//      }
//       }
    }
  }

  bcc->electrons->clear();

  for(unsigned int i=0; i<good_eles.size(); ++i){
    bcc->electrons->push_back(good_eles[i]);
  }
  sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
}


void Cleaner::ElectronCleaner(double ptmin, double etamax, double relisomax){

  ElectronCleaner_noIso(ptmin, etamax);
  std::vector<Electron> good_eles;
  for(unsigned int i=0; i<bcc->electrons->size(); ++i){
    Electron ele = bcc->electrons->at(i);
    if(ele.relIso()<relisomax){
      good_eles.push_back(ele);
    }
  }
  bcc->electrons->clear();

  for(unsigned int i=0; i<good_eles.size(); ++i){
    bcc->electrons->push_back(good_eles[i]);
  }
  sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
  resetEventCalc();
}


void Cleaner::MuonCleaner_noIso(double ptmin, double etamax){

  std::vector<Muon> good_mus;
  for(unsigned int i=0; i<bcc->muons->size(); ++i){
    Muon mu = bcc->muons->at(i);
    if(mu.pt()>ptmin){
      if(fabs(mu.eta())<etamax){
//      if(mu.isGlobalMuon()){
//        if(mu.globalTrack_chi2()/mu.globalTrack_ndof()<10){
//          if(mu.innerTrack_trackerLayersWithMeasurement()>5){
//            if(mu.dB()<0.02){
//              if(fabs(mu.vertex_z()-bcc->pvs->at(0).z())<1){
//                if(mu.innerTrack_numberOfValidPixelHits()>0){
//                  if(mu.numberOfMatchedStations()>1){
                      good_mus.push_back(mu);
 //                 }
//                }
//              }
//            }
//          }
//        }
//      }
      }
    }
  }
  bcc->muons->clear();

  for(unsigned int i=0; i<good_mus.size(); ++i){
    bcc->muons->push_back(good_mus[i]);
  }
  sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
  resetEventCalc();
}

void Cleaner::MuonCleaner(double ptmin, double etamax, double relisomax){

  MuonCleaner_noIso(ptmin, etamax);

  std::vector<Muon> good_mus;
  for(unsigned int i=0; i<bcc->muons->size(); ++i){
    Muon mu = bcc->muons->at(i);
    if(mu.relIso()<relisomax){
      good_mus.push_back(mu);
    }
  }
  bcc->muons->clear();

  for(unsigned int i=0; i<good_mus.size(); ++i){
    bcc->muons->push_back(good_mus[i]);
  }
  sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
}

void Cleaner::TauCleaner(double ptmin, double etamax){

  std::vector<Tau> good_taus;
  for(unsigned int i=0; i<bcc->taus->size(); ++i){
    Tau tau = bcc->taus->at(i);
    if(tau.pt()>ptmin){
      if(fabs(tau.eta())<etamax){
        if(bcc->taus->at(i).decayModeFinding()){
          if(bcc->taus->at(i).byMediumCombinedIsolationDeltaBetaCorr()){
            if(bcc->taus->at(i).againstElectronTight()){
              if(bcc->taus->at(i).againstMuonTight()){
                good_taus.push_back(tau);
              }
            }
          }
        }
      }
    }
  }
  
  bcc->taus->clear();

  for(unsigned int i=0; i<good_taus.size(); ++i){
    bcc->taus->push_back(good_taus[i]);
  }
  sort(bcc->taus->begin(), bcc->taus->end(), HigherPt());
  resetEventCalc();
}


void Cleaner::JetCleaner(double ptmin, double etamax, bool doPFID){

  std::vector<Jet> good_jets;
  for(unsigned int i=0; i<bcc->jets->size(); ++i){
    Jet jet = bcc->jets->at(i);
    if(jet.pt()>ptmin){
      if(fabs(jet.eta())<etamax){
        if(!doPFID || pfID(jet)){
          good_jets.push_back(jet);
        }
      }
    }
  }

  bcc->jets->clear();

  for(unsigned int i=0; i<good_jets.size(); ++i){
    bcc->jets->push_back(good_jets[i]);
  }
  sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
  resetEventCalc();
}

void Cleaner::TopJetCleaner(double ptmin, double etamax, bool doPFID){

  std::vector<TopJet> good_topjets;
  for(unsigned int i=0; i<bcc->topjets->size(); ++i){
    TopJet topjet = bcc->topjets->at(i);
    if(topjet.pt()>ptmin){
      if(fabs(topjet.eta())<etamax){
        if(!doPFID || pfID(topjet)){
          good_topjets.push_back(topjet);
        }
      }
    }
  }
  bcc->topjets->clear();

  for(unsigned int i=0; i<good_topjets.size(); ++i){
    bcc->topjets->push_back(good_topjets[i]);
  }
  sort(bcc->topjets->begin(), bcc->topjets->end(), HigherPt());
  resetEventCalc();
}


void Cleaner::PrimaryVertexCleaner(int ndofmax, double zmax, double rhomax){

  std::vector<PrimaryVertex> good_pvs;
  for(unsigned int i=0; i<bcc->pvs->size(); ++i){
    PrimaryVertex pv = bcc->pvs->at(i);
    if(pv.ndof()>=ndofmax && fabs(pv.z())<zmax && pv.rho()<=rhomax){
      good_pvs.push_back(pv);
    } 
  }
  bcc->pvs->clear();

  for(unsigned int i=0; i<good_pvs.size(); ++i){
    bcc->pvs->push_back(good_pvs[i]);
  }
  resetEventCalc();
}
Revision:	1.14
Committed:	Mon Jun 18 16:23:45 2012 UTC (12 years, 10 months ago) by peiffer
Content type:	text/plain
Branch:	MAIN
Changes since 1.13:	+0 -15 lines
Log Message:	sync update
#	User	Rev	Content
1	peiffer	1.1	#include "../include/Cleaner.h"
2
3
4	peiffer	1.7	Cleaner::Cleaner( BaseCycleContainer* input){
5	peiffer	1.1
6	peiffer	1.7	bcc = input;
7	peiffer	1.1
8			}
9
10	peiffer	1.9	Cleaner::Cleaner(){
11
12			ObjectHandler* objs = ObjectHandler::Instance();
13			bcc = objs->GetBaseCycleContainer();
14
15			}
16
17			void Cleaner::resetEventCalc(){
18
19			EventCalc* evc = EventCalc::Instance();
20			evc->Reset();
21
22			}
23
24	peiffer	1.13	void Cleaner::JetEnergyResolutionShifter(E_SystShift syst_shift){
25	peiffer	1.2
26	peiffer	1.13	LorentzVector met(0,0,0,0);
27			if(bcc->met) {
28			met.SetPt(bcc->met->pt());
29			met.SetPhi(bcc->met->phi());
30			}
31	peiffer	1.2
32			for(unsigned int i=0; i<bcc->jets->size(); ++i){
33			Jet jet = bcc->jets->at(i);
34	peiffer	1.5	float gen_pt = jet.genjet_pt();
35	peiffer	1.2	//ignore unmatched jets (which have zero vector) or jets with very low pt:
36			if(gen_pt < 15.0) continue;
37
38	peiffer	1.13	met += jet.v4()*jet.JEC_factor_raw();
39	peiffer	1.2
40	peiffer	1.5	float recopt = jet.pt();
41	peiffer	1.2	double factor = -10;
42	peiffer	1.5	double abseta = fabs(jet.eta());
43	peiffer	1.2
44			//numbers taken from https://twiki.cern.ch/twiki/bin/view/CMS/JetResolution
45	peiffer	1.13	if(syst_shift==e_Default){
46	peiffer	1.2	if(abseta < 0.5)
47			factor = 0.052;
48			else if(abseta >= 0.5 && abseta <1.1)
49			factor = 0.057;
50			else if(abseta >= 1.1 && abseta <1.7)
51			factor = 0.096;
52			else if(abseta >= 1.7 && abseta <2.3)
53			factor = 0.134;
54			else if(abseta >= 2.3)
55			factor = 0.288;
56			}
57	peiffer	1.13	else if(syst_shift==e_Up){
58	peiffer	1.2	if(abseta < 0.5)
59			factor = 0.11;
60			else if(abseta >= 0.5 && abseta <1.1)
61			factor = 0.12;
62			else if(abseta >= 1.1 && abseta <1.7)
63			factor = 0.16;
64			else if(abseta >= 1.7 && abseta <2.3)
65			factor = 0.23;
66			else if(abseta >= 2.3)
67			factor = 0.49;
68			}
69	peiffer	1.13	else if(syst_shift==e_Down){
70	peiffer	1.2	if(abseta < 0.5)
71			factor = -0.01;
72			else if(abseta >= 0.5 && abseta <1.1)
73			factor = 0.00;
74			else if(abseta >= 1.1 && abseta <1.7)
75			factor = 0.04;
76			else if(abseta >= 1.7 && abseta <2.3)
77			factor = 0.03;
78			else if(abseta >= 2.3)
79			factor = 0.09;
80			}
81
82			float deltapt = (recopt - gen_pt) * factor;
83			float ptscale = std::max(0.0f, (recopt + deltapt) / recopt);
84	peiffer	1.5	jet.set_pt(jet.pt() * ptscale);
85			bcc->jets->at(i).set_pt( jet.pt());
86	peiffer	1.2
87			//propagate JER shifts to MET
88	peiffer	1.13	met -= jet.v4()*jet.JEC_factor_raw();
89	peiffer	1.2	}
90
91			//store changed MET, flip phi if new MET is negative
92	peiffer	1.6	if(bcc->met){
93	peiffer	1.13	bcc->met->set_pt(met.Pt());
94			bcc->met->set_phi(met.Phi());
95	peiffer	1.2	}
96	peiffer	1.6
97	peiffer	1.5	sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
98	peiffer	1.9	resetEventCalc();
99	peiffer	1.2	}
100
101	peiffer	1.1
102	peiffer	1.7	void Cleaner::JetLeptonSubtractor(FactorizedJetCorrector *corrector){
103
104			for(unsigned int i=0; i<bcc->jets->size(); ++i){
105
106			LorentzVector jet_v4_raw = bcc->jets->at(i).v4()*bcc->jets->at(i).JEC_factor_raw();
107
108			//subtract lepton momenta from raw jet momentum
109	peiffer	1.13
110	peiffer	1.11	double ele_energy = bcc->jets->at(i).chargedEmEnergyFraction()*jet_v4_raw.E();
111			double mu_energy = bcc->jets->at(i).muonEnergyFraction()*jet_v4_raw.E();
112
113	peiffer	1.7	if(bcc->electrons){
114			for(unsigned int j=0; j<bcc->electrons->size(); ++j){
115			if(bcc->jets->at(i).deltaR(bcc->electrons->at(j))<0.5){
116	peiffer	1.13	// if(jet_v4_raw.pt() >= bcc->electrons->at(j).pt()){
117			jet_v4_raw -= bcc->electrons->at(j).v4();
118			// bcc->jets->at(i).set_electronMultiplicity(bcc->jets->at(i).electronMultiplicity()-1);
119			// ele_energy -= bcc->electrons->at(j).energy();
120			//}
121			//else{
122			//std::cout << "Subtracting electron with pt="<< bcc->electrons->at(j).pt() <<" from jet with pt(corrected)= " << bcc->jets->at(i).pt() << ", pt(raw)=" << jet_v4_raw.pt() << std::endl;
123			//jet_v4_raw -= jet_v4_raw;
124			//}
125	peiffer	1.7	}
126			}
127			}
128			if(bcc->muons){
129			for(unsigned int j=0; j<bcc->muons->size(); ++j){
130	peiffer	1.10	if(bcc->jets->at(i).deltaR(bcc->muons->at(j))<0.5){
131	peiffer	1.13	//if(jet_v4_raw.pt()>= bcc->muons->at(j).pt()){
132	peiffer	1.10	jet_v4_raw -= bcc->muons->at(j).v4();
133	peiffer	1.13	// bcc->jets->at(i).set_muonMultiplicity(bcc->jets->at(i).muonMultiplicity()-1);
134			// mu_energy -= bcc->muons->at(j).energy();
135			//}
136			//else{
137			//std::cout << "Subtracting muon with pt="<< bcc->muons->at(j).pt() <<" from jet with pt(corrected)= " << bcc->jets->at(i).pt() << ", pt(raw)=" << jet_v4_raw.pt() << std::endl;
138			//jet_v4_raw -= jet_v4_raw;
139			//}
140	peiffer	1.10	}
141	peiffer	1.7	}
142			}
143	peiffer	1.11	if(ele_energy<=jet_v4_raw.E())
144			bcc->jets->at(i).set_chargedEmEnergyFraction(ele_energy/jet_v4_raw.E());
145			if(mu_energy<=jet_v4_raw.E())
146	peiffer	1.13	bcc->jets->at(i).set_muonEnergyFraction(mu_energy/jet_v4_raw.E());
147	peiffer	1.11
148	peiffer	1.7
149			//apply jet energy corrections to modified raw momentum
150			corrector->setJetPt(jet_v4_raw.Pt());
151			corrector->setJetEta(jet_v4_raw.Eta());
152			corrector->setJetE(jet_v4_raw.E());
153			corrector->setJetA(bcc->jets->at(i).jetArea());
154			corrector->setRho(bcc->rho);
155	peiffer	1.13	corrector->setNPV(bcc->pvs->size());
156	peiffer	1.7
157			float correctionfactor = corrector->getCorrection();
158	peiffer	1.13
159	peiffer	1.7	LorentzVector jet_v4_corrected = jet_v4_raw *correctionfactor;
160
161			bcc->jets->at(i).set_v4(jet_v4_corrected);
162			bcc->jets->at(i).set_JEC_factor_raw(1./correctionfactor);
163
164			}
165
166	peiffer	1.8	sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
167	peiffer	1.9	resetEventCalc();
168	peiffer	1.7	}
169
170	peiffer	1.1	//tight ele ID from https://twiki.cern.ch/twiki/bin/view/CMS/EgammaCutBasedIdentification
171			bool Cleaner::eleID(Electron ele){
172
173			bool pass=false;
174
175	peiffer	1.5	if(fabs(ele.supercluster_eta())<1.4442){
176			if(ele.dEtaIn()<0.004 && ele.dPhiIn()<0.03 && ele.sigmaIEtaIEta()<0.01 && ele.HoverE()<0.12 && fabs(1./ele.EcalEnergy()-1./ele.v4().P())<0.05) pass=true;
177	peiffer	1.1	}
178	peiffer	1.5	else if( fabs(ele.supercluster_eta())>1.5660){
179			if(ele.dEtaIn()<0.005 && ele.dPhiIn()<0.02 && ele.sigmaIEtaIEta()<0.03 && ele.HoverE()<0.10 && fabs(1./ele.EcalEnergy()-1./ele.v4().P())<0.05) pass=true;
180	peiffer	1.1	}
181
182			return pass;
183
184			}
185
186	peiffer	1.2	//pf ID has already been applied when using goodPatJets
187			bool Cleaner::pfID(Jet jet){
188
189	peiffer	1.5	if(jet.numberOfDaughters()>1
190			&& jet.neutralHadronEnergyFraction()<0.99
191			&& jet.neutralEmEnergyFraction()<0.99){
192	peiffer	1.2
193	peiffer	1.5	if(fabs(jet.eta())>=2.4)
194	peiffer	1.2	return true;
195
196	peiffer	1.5	if(jet.chargedEmEnergyFraction()<0.99
197			&& jet.chargedHadronEnergyFraction()>0
198			&& jet.chargedMultiplicity()>0)
199	peiffer	1.2	return true;
200
201			}
202	peiffer	1.13	std::cout << "Bloeder Jet, pt,eta= " << jet.pt() << ", " <<jet.eta() << "; " << jet.numberOfDaughters() << " " << jet.neutralHadronEnergyFraction() << " " <<jet.neutralEmEnergyFraction() << " " << jet.chargedEmEnergyFraction() << " " << jet.chargedHadronEnergyFraction() << " " << jet.chargedMultiplicity() <<std::endl;
203	peiffer	1.2	return false;
204
205			}
206
207	peiffer	1.1
208	peiffer	1.7	void Cleaner::ElectronCleaner_noIso(double ptmin, double etamax){
209	peiffer	1.1
210	peiffer	1.3	std::vector<Electron> good_eles;
211	peiffer	1.1	for(unsigned int i=0; i<bcc->electrons->size(); ++i){
212			Electron ele = bcc->electrons->at(i);
213	peiffer	1.5	if(ele.pt()>ptmin){
214			if(fabs(ele.eta())<etamax){
215	peiffer	1.13	// if(fabs(ele.supercluster_eta())<1.4442 \|\| fabs(ele.supercluster_eta())>1.5660){
216			// if(bcc->pvs->size()>0){
217			// if(ele.gsfTrack_dxy_vertex(bcc->pvs->at(0).x(), bcc->pvs->at(0).y())<0.02){
218			// if(ele.passconversionveto()){
219			// if(ele.mvaTrigV0()>0.0){
220			// if(eleID(ele)){
221	peiffer	1.3	good_eles.push_back(ele);
222	peiffer	1.13	}
223			// }
224			// }
225			// }
226			// }
227			// }
228			// }
229	peiffer	1.1	}
230			}
231
232	peiffer	1.3	bcc->electrons->clear();
233
234			for(unsigned int i=0; i<good_eles.size(); ++i){
235			bcc->electrons->push_back(good_eles[i]);
236			}
237	peiffer	1.5	sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
238	peiffer	1.1	}
239
240	rkogler	1.12
241	peiffer	1.7	void Cleaner::ElectronCleaner(double ptmin, double etamax, double relisomax){
242	peiffer	1.1
243	peiffer	1.7	ElectronCleaner_noIso(ptmin, etamax);
244			std::vector<Electron> good_eles;
245			for(unsigned int i=0; i<bcc->electrons->size(); ++i){
246			Electron ele = bcc->electrons->at(i);
247			if(ele.relIso()<relisomax){
248			good_eles.push_back(ele);
249			}
250			}
251			bcc->electrons->clear();
252
253			for(unsigned int i=0; i<good_eles.size(); ++i){
254			bcc->electrons->push_back(good_eles[i]);
255			}
256			sort(bcc->electrons->begin(), bcc->electrons->end(), HigherPt());
257	peiffer	1.9	resetEventCalc();
258	peiffer	1.7	}
259
260
261			void Cleaner::MuonCleaner_noIso(double ptmin, double etamax){
262	peiffer	1.3
263			std::vector<Muon> good_mus;
264	peiffer	1.1	for(unsigned int i=0; i<bcc->muons->size(); ++i){
265			Muon mu = bcc->muons->at(i);
266	peiffer	1.5	if(mu.pt()>ptmin){
267			if(fabs(mu.eta())<etamax){
268	peiffer	1.13	// if(mu.isGlobalMuon()){
269			// if(mu.globalTrack_chi2()/mu.globalTrack_ndof()<10){
270			// if(mu.innerTrack_trackerLayersWithMeasurement()>5){
271			// if(mu.dB()<0.02){
272			// if(fabs(mu.vertex_z()-bcc->pvs->at(0).z())<1){
273			// if(mu.innerTrack_numberOfValidPixelHits()>0){
274			// if(mu.numberOfMatchedStations()>1){
275	peiffer	1.7	good_mus.push_back(mu);
276	peiffer	1.13	// }
277			// }
278			// }
279			// }
280			// }
281			// }
282			// }
283	peiffer	1.1	}
284			}
285			}
286	peiffer	1.3	bcc->muons->clear();
287
288			for(unsigned int i=0; i<good_mus.size(); ++i){
289			bcc->muons->push_back(good_mus[i]);
290			}
291	peiffer	1.5	sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
292	peiffer	1.9	resetEventCalc();
293	peiffer	1.1	}
294	peiffer	1.2
295	peiffer	1.7	void Cleaner::MuonCleaner(double ptmin, double etamax, double relisomax){
296
297			MuonCleaner_noIso(ptmin, etamax);
298
299			std::vector<Muon> good_mus;
300			for(unsigned int i=0; i<bcc->muons->size(); ++i){
301			Muon mu = bcc->muons->at(i);
302			if(mu.relIso()<relisomax){
303			good_mus.push_back(mu);
304			}
305			}
306			bcc->muons->clear();
307
308			for(unsigned int i=0; i<good_mus.size(); ++i){
309			bcc->muons->push_back(good_mus[i]);
310			}
311			sort(bcc->muons->begin(), bcc->muons->end(), HigherPt());
312			}
313
314	peiffer	1.4	void Cleaner::TauCleaner(double ptmin, double etamax){
315
316			std::vector<Tau> good_taus;
317			for(unsigned int i=0; i<bcc->taus->size(); ++i){
318			Tau tau = bcc->taus->at(i);
319	peiffer	1.5	if(tau.pt()>ptmin){
320			if(fabs(tau.eta())<etamax){
321			if(bcc->taus->at(i).decayModeFinding()){
322			if(bcc->taus->at(i).byMediumCombinedIsolationDeltaBetaCorr()){
323			if(bcc->taus->at(i).againstElectronTight()){
324			if(bcc->taus->at(i).againstMuonTight()){
325	peiffer	1.4	good_taus.push_back(tau);
326			}
327			}
328			}
329			}
330			}
331			}
332			}
333
334			bcc->taus->clear();
335
336			for(unsigned int i=0; i<good_taus.size(); ++i){
337			bcc->taus->push_back(good_taus[i]);
338			}
339	peiffer	1.5	sort(bcc->taus->begin(), bcc->taus->end(), HigherPt());
340	peiffer	1.9	resetEventCalc();
341	peiffer	1.4	}
342
343	peiffer	1.2
344	peiffer	1.3	void Cleaner::JetCleaner(double ptmin, double etamax, bool doPFID){
345
346			std::vector<Jet> good_jets;
347	peiffer	1.2	for(unsigned int i=0; i<bcc->jets->size(); ++i){
348			Jet jet = bcc->jets->at(i);
349	peiffer	1.5	if(jet.pt()>ptmin){
350			if(fabs(jet.eta())<etamax){
351	peiffer	1.2	if(!doPFID \|\| pfID(jet)){
352	peiffer	1.3	good_jets.push_back(jet);
353	peiffer	1.2	}
354			}
355			}
356			}
357
358	peiffer	1.3	bcc->jets->clear();
359
360			for(unsigned int i=0; i<good_jets.size(); ++i){
361			bcc->jets->push_back(good_jets[i]);
362			}
363	peiffer	1.5	sort(bcc->jets->begin(), bcc->jets->end(), HigherPt());
364	peiffer	1.9	resetEventCalc();
365	peiffer	1.2	}
366
367	peiffer	1.3	void Cleaner::TopJetCleaner(double ptmin, double etamax, bool doPFID){
368	peiffer	1.2
369	peiffer	1.3	std::vector<TopJet> good_topjets;
370	peiffer	1.2	for(unsigned int i=0; i<bcc->topjets->size(); ++i){
371			TopJet topjet = bcc->topjets->at(i);
372	peiffer	1.5	if(topjet.pt()>ptmin){
373			if(fabs(topjet.eta())<etamax){
374	peiffer	1.2	if(!doPFID \|\| pfID(topjet)){
375	peiffer	1.3	good_topjets.push_back(topjet);
376	peiffer	1.2	}
377			}
378			}
379			}
380	peiffer	1.3	bcc->topjets->clear();
381	peiffer	1.2
382	peiffer	1.3	for(unsigned int i=0; i<good_topjets.size(); ++i){
383			bcc->topjets->push_back(good_topjets[i]);
384			}
385	peiffer	1.5	sort(bcc->topjets->begin(), bcc->topjets->end(), HigherPt());
386	peiffer	1.9	resetEventCalc();
387	peiffer	1.2	}
388	peiffer	1.7
389
390			void Cleaner::PrimaryVertexCleaner(int ndofmax, double zmax, double rhomax){
391
392			std::vector<PrimaryVertex> good_pvs;
393			for(unsigned int i=0; i<bcc->pvs->size(); ++i){
394			PrimaryVertex pv = bcc->pvs->at(i);
395			if(pv.ndof()>=ndofmax && fabs(pv.z())<zmax && pv.rho()<=rhomax){
396			good_pvs.push_back(pv);
397			}
398			}
399			bcc->pvs->clear();
400
401			for(unsigned int i=0; i<good_pvs.size(); ++i){
402			bcc->pvs->push_back(good_pvs[i]);
403			}
404	peiffer	1.9	resetEventCalc();
405	peiffer	1.7	}