1 |
|
2 |
#include "PhysicsTools/PatUtils/interface/CaloJetSelector.h"
|
3 |
#include "DataFormats/Math/interface/deltaR.h"
|
4 |
|
5 |
using pat::CaloJetSelector;
|
6 |
|
7 |
//______________________________________________________________________________
|
8 |
const pat::ParticleStatus
|
9 |
CaloJetSelector::filter( //const unsigned int& index,
|
10 |
// const edm::View<reco::CaloJet>& Jets
|
11 |
const reco::CaloJet& Jet
|
12 |
) const
|
13 |
{
|
14 |
ParticleStatus result = GOOD;
|
15 |
|
16 |
///Retrieve information
|
17 |
///Pt Jet
|
18 |
if (Jet.p4().Pt()<config_.Ptmin) result = BAD;
|
19 |
|
20 |
///eta region
|
21 |
double eta = fabs(Jet.p4().Eta());
|
22 |
if (eta>config_.Etamax) result = BAD;
|
23 |
|
24 |
///electromagnetic fraction
|
25 |
double EMF = Jet.emEnergyFraction();
|
26 |
if (EMF<config_.EMFmin ||
|
27 |
EMF>config_.EMFmax ) result = BAD;
|
28 |
|
29 |
///(EMCalEnergyFraction + HadCalEnergyFraction) / (EMCalEnergyFraction - HadCalEnergyFraction)
|
30 |
double HadF = Jet.emEnergyFraction();
|
31 |
double EMvsHadF = 0.;
|
32 |
if (EMF-HadF!=0.) EMvsHadF = (EMF+HadF)/(EMF-HadF);
|
33 |
if (EMvsHadF<config_.EMvsHadFmin ||
|
34 |
EMvsHadF>config_.EMvsHadFmax ) result = BAD;
|
35 |
|
36 |
//ratio Energy over Momentum (both from calorimeter)
|
37 |
//Useful? Both come from a lorentz-vector
|
38 |
//double EoverP = 0.;
|
39 |
//if (Jet.p4().P()!=0.) EoverP = Jet.energy() / Jet.p4().P();
|
40 |
//if (EoverP > config_.EoverPmax) result = BAD;
|
41 |
|
42 |
///n90: number of towers containing 90% of the jet's energy
|
43 |
double n90 = Jet.n90();
|
44 |
if (n90<config_.N90min ||
|
45 |
n90>config_.N90max ) result = BAD;
|
46 |
|
47 |
///Tower Number
|
48 |
std::vector<CaloTowerPtr> jetTowers = Jet.getCaloConstituents();
|
49 |
if (jetTowers.size()<config_.NCaloTowersmin ||
|
50 |
jetTowers.size()>config_.NCaloTowersmax ) result = BAD;
|
51 |
|
52 |
//calculate tower related variables:
|
53 |
double MaxEnergyTower = 0., SumTowPt=0., SumTowPtR=0.;
|
54 |
for(std::vector<CaloTowerPtr>::const_iterator tow = jetTowers.begin(),
|
55 |
towend = jetTowers.end(); tow != towend; ++tow){
|
56 |
|
57 |
SumTowPt += (*tow)->et();
|
58 |
SumTowPtR += (*tow)->et()*deltaR( Jet.p4().Eta(), Jet.p4().Phi(),
|
59 |
(*tow)->eta(), (*tow)->phi() );
|
60 |
if ( (*tow)->et() > MaxEnergyTower )
|
61 |
MaxEnergyTower = (*tow)->et();
|
62 |
}
|
63 |
|
64 |
///Highest Et Tower / Et Jet
|
65 |
double EtTowerMaxOverEtJet = 0.;
|
66 |
if (Jet.p4().Et()!=0.) EtTowerMaxOverEtJet = MaxEnergyTower /Jet.p4().Et();
|
67 |
if (EtTowerMaxOverEtJet < config_.HighestTowerOverJetmin ||
|
68 |
EtTowerMaxOverEtJet > config_.HighestTowerOverJetmax ) result = BAD;
|
69 |
|
70 |
///Sum(E Twr * DeltaR(Twr-Jet)) / Sum(E Twr)
|
71 |
double RWidth = 0.;
|
72 |
if (SumTowPt!=0.) RWidth = SumTowPtR /SumTowPt;
|
73 |
if (RWidth < config_.RWidthmin ||
|
74 |
RWidth > config_.RWidthmax ) result = BAD;
|
75 |
|
76 |
///Pt Jet / Towers Area
|
77 |
double PtJetoverArea = 0.;
|
78 |
if (Jet.towersArea() !=0.) PtJetoverArea = Jet.p4().Pt() / Jet.towersArea();
|
79 |
if (PtJetoverArea < config_.PtJetOverArea_min ||
|
80 |
PtJetoverArea > config_.PtJetOverArea_max ) result = BAD;
|
81 |
|
82 |
///Highest Et Tower / Towers Area
|
83 |
double PtToweroverArea = 0.;
|
84 |
if (Jet.towersArea() !=0.) PtToweroverArea = MaxEnergyTower / Jet.towersArea();
|
85 |
if (PtToweroverArea < config_.PtTowerOverArea_min ||
|
86 |
PtToweroverArea > config_.PtTowerOverArea_max ) result = BAD;
|
87 |
|
88 |
return result;
|
89 |
}
|