CMS 3D CMS Logo

PFBlockBasedIsolation.cc
Go to the documentation of this file.
2 #include <cmath>
4 
19 
20 //--------------------------------------------------------------------------------------------------
21 
23  // Default Constructor.
24 }
25 
26 //--------------------------------------------------------------------------------------------------
28 
29 void PFBlockBasedIsolation::setup(const edm::ParameterSet& conf) { coneSize_ = conf.getParameter<double>("coneSize"); }
30 
31 std::vector<reco::PFCandidateRef> PFBlockBasedIsolation::calculate(
33  const reco::PFCandidateRef pfEGCand,
34  const edm::Handle<reco::PFCandidateCollection> pfCandidateHandle) {
35  std::vector<reco::PFCandidateRef> myVec;
36 
37  math::XYZVector candidateMomentum(p4.px(), p4.py(), p4.pz());
38  math::XYZVector candidateDirection = candidateMomentum.Unit();
39 
40  const reco::PFCandidate::ElementsInBlocks& theElementsInpfEGcand = (*pfEGCand).elementsInBlocks();
41  reco::PFCandidate::ElementsInBlocks::const_iterator ieg = theElementsInpfEGcand.begin();
42  const reco::PFBlockRef egblock = ieg->first;
43 
44  unsigned nObj = pfCandidateHandle->size();
45  for (unsigned int lCand = 0; lCand < nObj; lCand++) {
46  reco::PFCandidateRef pfCandRef(reco::PFCandidateRef(pfCandidateHandle, lCand));
47 
48  float dR = 0.0;
49  if (coneSize_ < 10.0) {
50  dR = deltaR(candidateDirection.Eta(), candidateDirection.Phi(), pfCandRef->eta(), pfCandRef->phi());
51  if (dR > coneSize_)
52  continue;
53  }
54 
55  const reco::PFCandidate::ElementsInBlocks& theElementsInPFcand = pfCandRef->elementsInBlocks();
56 
57  bool elementFound = false;
58  for (reco::PFCandidate::ElementsInBlocks::const_iterator ipf = theElementsInPFcand.begin();
59  ipf < theElementsInPFcand.end();
60  ++ipf) {
61  if (ipf->first == egblock && !elementFound) {
62  for (ieg = theElementsInpfEGcand.begin(); ieg < theElementsInpfEGcand.end(); ++ieg) {
63  if (ipf->second == ieg->second && !elementFound) {
64  if (elementPassesCleaning(pfCandRef, pfEGCand)) {
65  myVec.push_back(pfCandRef);
66  elementFound = true;
67  }
68  }
69  }
70  }
71  }
72  }
73 
74  return myVec;
75 }
76 
78  const reco::PFCandidateRef& pfEGCand) {
79  if (pfCand->particleId() == reco::PFCandidate::h)
80  return passesCleaningChargedHadron(pfCand, pfEGCand);
81  else if (pfCand->particleId() == reco::PFCandidate::h0)
82  return passesCleaningNeutralHadron(pfCand, pfEGCand);
83  else if (pfCand->particleId() == reco::PFCandidate::gamma)
84  return passesCleaningPhoton(pfCand, pfEGCand);
85  else
86  return true; //doesnt really matter here as if its not a photon,neutral or charged it wont be included in isolation
87 }
88 
89 //currently the record of which candidates came from the charged hadron is acceptable, no further cleaning is needed
91  const reco::PFCandidateRef& pfEGCand) {
92  return true;
93 }
94 
95 //neutral hadrons are not part of the PF E/gamma reco, therefore they cant currently come from an electron/photon and so should be rejected
96 //but we still think there may be some useful info here and given we can easily
97 //fix this at AOD level, we will auto accept them for now and clean later
99  const reco::PFCandidateRef& pfEGCand) {
100  return true;
101 }
102 
103 //the highest et ECAL element of the photon must match to the electron superclusters or one of its sub clusters
105  const reco::PFCandidateRef& pfEGCand) {
106  bool passesCleaning = false;
107  const reco::PFBlockElementCluster* ecalClusWithMaxEt = getHighestEtECALCluster(*pfCand);
108  if (ecalClusWithMaxEt) {
109  if (ecalClusWithMaxEt->superClusterRef().isNonnull() &&
110  ecalClusWithMaxEt->superClusterRef()->seed()->seed() ==
111  pfEGCand->superClusterRef()
112  ->seed()
113  ->seed()) { //being sure to match, some concerned about different collections, shouldnt be but to be safe
114  passesCleaning = true;
115  } else {
116  for (auto cluster : pfEGCand->superClusterRef()->clusters()) {
117  //the PF clusters there are in two different collections so cant reference match
118  //but we can match on the seed id, no clusters can share a seed so if the seeds are
119  //equal, it must be the same cluster
120  if (ecalClusWithMaxEt->clusterRef()->seed() == cluster->seed()) {
121  passesCleaning = true;
122  }
123  } //end of loop over clusters
124  }
125  } //end of null check for highest ecal cluster
126  return passesCleaning;
127 }
128 
130  float maxECALEt = -1;
131  const reco::PFBlockElement* maxEtECALCluster = nullptr;
132  const reco::PFCandidate::ElementsInBlocks& elementsInPFCand = pfCand.elementsInBlocks();
133  for (auto& elemIndx : elementsInPFCand) {
134  const reco::PFBlockElement* elem =
135  elemIndx.second < elemIndx.first->elements().size() ? &elemIndx.first->elements()[elemIndx.second] : nullptr;
136  if (elem && elem->type() == reco::PFBlockElement::ECAL && elem->clusterRef()->pt() > maxECALEt) {
137  maxECALEt = elem->clusterRef()->pt();
138  maxEtECALCluster = elem;
139  }
140  }
141  return dynamic_cast<const reco::PFBlockElementCluster*>(maxEtECALCluster);
142 }
PFBlockBasedIsolation::getHighestEtECALCluster
const reco::PFBlockElementCluster * getHighestEtECALCluster(const reco::PFCandidate &pfCand)
Definition: PFBlockBasedIsolation.cc:129
PFBlockElementCluster.h
PFCandidate.h
math::XYZTLorentzVectorD
ROOT::Math::LorentzVector< ROOT::Math::PxPyPzE4D< double > > XYZTLorentzVectorD
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:14
reco::PFCandidate::h
Definition: PFCandidate.h:45
PFBlockBasedIsolation::calculate
std::vector< reco::PFCandidateRef > calculate(math::XYZTLorentzVectorD p4, const reco::PFCandidateRef pfEGCand, const edm::Handle< reco::PFCandidateCollection > pfCandidateHandle)
Definition: PFBlockBasedIsolation.cc:31
EcalClusterLazyTools.h
PFBlockBasedIsolation::~PFBlockBasedIsolation
~PFBlockBasedIsolation()
Definition: PFBlockBasedIsolation.cc:27
PFBlockBasedIsolation::passesCleaningPhoton
bool passesCleaningPhoton(const reco::PFCandidateRef &pfCand, const reco::PFCandidateRef &pfEGCand)
Definition: PFBlockBasedIsolation.cc:104
reco::PFCandidate::elementsInBlocks
const ElementsInBlocks & elementsInBlocks() const
Definition: PFCandidate.cc:636
PFBlockBasedIsolation::passesCleaningNeutralHadron
bool passesCleaningNeutralHadron(const reco::PFCandidateRef &pfCand, const reco::PFCandidateRef &pfEGCand)
Definition: PFBlockBasedIsolation.cc:98
edm::Handle
Definition: AssociativeIterator.h:50
PFBlockBasedIsolation::setup
void setup(const edm::ParameterSet &conf)
Definition: PFBlockBasedIsolation.cc:29
edm::Ref< PFCandidateCollection >
deltaR.h
Track.h
TrackFwd.h
PFCluster.h
reco::PFBlockElementCluster::superClusterRef
const SuperClusterRef & superClusterRef() const
Definition: PFBlockElementCluster.h:30
PbPb_ZMuSkimMuonDPG_cff.deltaR
deltaR
Definition: PbPb_ZMuSkimMuonDPG_cff.py:63
reco::PFCandidate::ElementsInBlocks
std::vector< ElementInBlock > ElementsInBlocks
Definition: PFCandidate.h:378
Vertex.h
TransientTrackBuilder.h
edm::ParameterSet
Definition: ParameterSet.h:36
edm::Ref::isNonnull
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
math::XYZVector
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
reco::PFBlockElement::ECAL
Definition: PFBlockElement.h:35
RefToPtr.h
p4
double p4[4]
Definition: TauolaWrapper.h:92
reco::PFCandidate::gamma
Definition: PFCandidate.h:48
PFBlockBasedIsolation.h
PFBlockBasedIsolation::coneSize_
double coneSize_
Definition: PFBlockBasedIsolation.h:55
reco::PFBlockElement
Abstract base class for a PFBlock element (track, cluster...)
Definition: PFBlockElement.h:26
PFBlockFwd.h
PFBlockBasedIsolation::elementPassesCleaning
bool elementPassesCleaning(const reco::PFCandidateRef &pfCand, const reco::PFCandidateRef &pfEGCand)
Definition: PFBlockBasedIsolation.cc:77
PFBlock.h
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
IPTools.h
reco::PFCandidate::h0
Definition: PFCandidate.h:49
reco::PFBlockElementCluster
Cluster Element.
Definition: PFBlockElementCluster.h:16
SuperCluster.h
reco::PFBlockElementCluster::clusterRef
const PFClusterRef & clusterRef() const override
Definition: PFBlockElementCluster.h:29
reco::PFBlockElement::type
Type type() const
Definition: PFBlockElement.h:69
PFBlockBasedIsolation::passesCleaningChargedHadron
bool passesCleaningChargedHadron(const reco::PFCandidateRef &pfCand, const reco::PFCandidateRef &pfEGCand)
Definition: PFBlockBasedIsolation.cc:90
reco::PFCandidate
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:40
HGC3DClusterGenMatchSelector_cfi.dR
dR
Definition: HGC3DClusterGenMatchSelector_cfi.py:7
PFBlockBasedIsolation::PFBlockBasedIsolation
PFBlockBasedIsolation()
Definition: PFBlockBasedIsolation.cc:22
reco::PFBlockElement::clusterRef
virtual const PFClusterRef & clusterRef() const
Definition: PFBlockElement.h:90