CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PFCandidateChecker.cc
Go to the documentation of this file.
1 
18 
19 namespace edm {
20  class EventSetup;
21  class Run;
22 } // namespace edm
23 
25 public:
26  explicit PFCandidateChecker(const edm::ParameterSet&);
27 
28  ~PFCandidateChecker() override;
29 
30  void analyze(const edm::Event&, const edm::EventSetup&) override;
31 
32  void beginRun(const edm::Run& r, const edm::EventSetup& c) override;
33 
34 private:
35  void printJets(const reco::PFJetCollection& pfJetsReco, const reco::PFJetCollection& pfJetsReReco) const;
36 
37  void printMet(const reco::PFCandidateCollection& pfReco, const reco::PFCandidateCollection& pfReReco) const;
38 
39  void printElementsInBlocks(const reco::PFCandidate& cand, std::ostream& out = std::cout) const;
40 
46 
48  double deltaEMax_;
49  double deltaEtaMax_;
50  double deltaPhiMax_;
51 
53  bool verbose_;
54 
57 
59  bool rankByPt_;
60 
62  unsigned entry_;
63 
64  static bool greaterPt(const reco::PFCandidate& a, const reco::PFCandidate& b) { return (a.pt() > b.pt()); }
65 };
66 
68 
69 using namespace std;
70 using namespace edm;
71 using namespace reco;
72 
74  inputTagPFCandidatesReco_ = iConfig.getParameter<InputTag>("pfCandidatesReco");
75 
76  inputTagPFCandidatesReReco_ = iConfig.getParameter<InputTag>("pfCandidatesReReco");
77 
78  inputTagPFJetsReco_ = iConfig.getParameter<InputTag>("pfJetsReco");
79 
80  inputTagPFJetsReReco_ = iConfig.getParameter<InputTag>("pfJetsReReco");
81 
82  deltaEMax_ = iConfig.getParameter<double>("deltaEMax");
83 
84  deltaEtaMax_ = iConfig.getParameter<double>("deltaEtaMax");
85 
86  deltaPhiMax_ = iConfig.getParameter<double>("deltaPhiMax");
87 
88  verbose_ = iConfig.getUntrackedParameter<bool>("verbose", false);
89 
90  printBlocks_ = iConfig.getUntrackedParameter<bool>("printBlocks", false);
91 
92  rankByPt_ = iConfig.getUntrackedParameter<bool>("rankByPt", false);
93 
94  entry_ = 0;
95 
96  LogDebug("PFCandidateChecker") << " input collections : " << inputTagPFCandidatesReco_ << " "
97  << inputTagPFCandidatesReReco_;
98 }
99 
101 
103 
104 void PFCandidateChecker::analyze(const Event& iEvent, const EventSetup& iSetup) {
105  LogDebug("PFCandidateChecker") << "START event: " << iEvent.id().event() << " in run " << iEvent.id().run() << endl;
106 
107  // get PFCandidates
108 
109  Handle<PFCandidateCollection> pfCandidatesReco;
110  iEvent.getByLabel(inputTagPFCandidatesReco_, pfCandidatesReco);
111 
112  Handle<PFCandidateCollection> pfCandidatesReReco;
113  iEvent.getByLabel(inputTagPFCandidatesReReco_, pfCandidatesReReco);
114 
115  Handle<PFJetCollection> pfJetsReco;
116  iEvent.getByLabel(inputTagPFJetsReco_, pfJetsReco);
117 
118  Handle<PFJetCollection> pfJetsReReco;
119  iEvent.getByLabel(inputTagPFJetsReReco_, pfJetsReReco);
120 
121  reco::PFCandidateCollection pfReco, pfReReco;
122 
123  // to sort, one needs to copy
124  if (rankByPt_) {
125  pfReco = *pfCandidatesReco;
126  pfReReco = *pfCandidatesReReco;
127  sort(pfReco.begin(), pfReco.end(), greaterPt);
128  sort(pfReReco.begin(), pfReReco.end(), greaterPt);
129  }
130 
131  unsigned minSize = pfReco.size() < pfReReco.size() ? pfReco.size() : pfReReco.size();
132  bool differentCand = false;
133  bool differentSize = pfReco.size() != pfReReco.size();
134  if (differentSize)
135  std::cout << "+++WARNING+++ PFCandidate size changed for entry " << entry_ << " !" << endl
136  << " - RECO size : " << pfReco.size() << endl
137  << " - Re-RECO size : " << pfReReco.size() << endl;
138 
139  unsigned npr = 0;
140  for (unsigned i = 0; i < minSize; i++) {
141  const reco::PFCandidate& candReco = (rankByPt_) ? pfReco[i] : (*pfCandidatesReco)[i];
142  const reco::PFCandidate& candReReco = (rankByPt_) ? pfReReco[i] : (*pfCandidatesReReco)[i];
143 
144  double deltaE = (candReReco.energy() - candReco.energy()) / (candReReco.energy() + candReco.energy());
145  double deltaEta = candReReco.eta() - candReco.eta();
146  double deltaPhi = candReReco.phi() - candReco.phi();
147  if (fabs(deltaE) > deltaEMax_ || fabs(deltaEta) > deltaEtaMax_ || fabs(deltaPhi) > deltaPhiMax_) {
148  differentCand = true;
149  std::cout << "+++WARNING+++ PFCandidate " << i << " changed for entry " << entry_ << " ! " << std::endl
150  << " - RECO : " << candReco << std::endl
151  << " - Re-RECO : " << candReReco << std::endl
152  << " DeltaE = : " << deltaE << std::endl
153  << " DeltaEta = : " << deltaEta << std::endl
154  << " DeltaPhi = : " << deltaPhi << std::endl
155  << std::endl;
156  if (printBlocks_) {
157  std::cout << "Elements in Block for RECO: " << std::endl;
158  printElementsInBlocks(candReco);
159  std::cout << "Elements in Block for Re-RECO: " << std::endl;
160  printElementsInBlocks(candReReco);
161  }
162  if (++npr == 5)
163  break;
164  }
165  }
166 
167  if (differentSize || differentCand) {
168  printJets(*pfJetsReco, *pfJetsReReco);
169  printMet(pfReco, pfReReco);
170  }
171 
172  ++entry_;
173  LogDebug("PFCandidateChecker") << "STOP event: " << iEvent.id().event() << " in run " << iEvent.id().run()
174  << std::endl;
175 }
176 
178  double metX = 0.;
179  double metY = 0.;
180  for (unsigned i = 0; i < pfReco.size(); i++) {
181  metX += pfReco[i].px();
182  metY += pfReco[i].py();
183  }
184  double met = std::sqrt(metX * metX + metY * metY);
185  std::cout << "MET RECO = " << metX << " " << metY << " " << met << std::endl;
186 
187  metX = 0.;
188  metY = 0.;
189  for (unsigned i = 0; i < pfReReco.size(); i++) {
190  metX += pfReReco[i].px();
191  metY += pfReReco[i].py();
192  }
193  met = std::sqrt(metX * metX + metY * metY);
194  std::cout << "MET Re-RECO = " << metX << " " << metY << " " << met << std::endl;
195 }
196 
197 void PFCandidateChecker::printJets(const PFJetCollection& pfJetsReco, const PFJetCollection& pfJetsReReco) const {
198  bool differentSize = pfJetsReco.size() != pfJetsReReco.size();
199  if (differentSize)
200  std::cout << "+++WARNING+++ PFJet size changed for entry " << entry_ << " !" << endl
201  << " - RECO size : " << pfJetsReco.size() << endl
202  << " - Re-RECO size : " << pfJetsReReco.size() << endl;
203  unsigned minSize = pfJetsReco.size() < pfJetsReReco.size() ? pfJetsReco.size() : pfJetsReReco.size();
204  unsigned npr = 0;
205  for (unsigned i = 0; i < minSize; ++i) {
206  const reco::PFJet& candReco = pfJetsReco[i];
207  const reco::PFJet& candReReco = pfJetsReReco[i];
208  if (candReco.et() < 20. && candReReco.et() < 20.)
209  break;
210  double deltaE = (candReReco.et() - candReco.et()) / (candReReco.et() + candReco.et());
211  double deltaEta = candReReco.eta() - candReco.eta();
212  double deltaPhi = candReReco.phi() - candReco.phi();
213  if (fabs(deltaE) > deltaEMax_ || fabs(deltaEta) > deltaEtaMax_ || fabs(deltaPhi) > deltaPhiMax_) {
214  std::cout << "+++WARNING+++ PFJet " << i << " changed for entry " << entry_ << " ! " << std::endl
215  << " - RECO : " << candReco.et() << " " << candReco.eta() << " " << candReco.phi() << std::endl
216  << " - Re-RECO : " << candReReco.et() << " " << candReReco.eta() << " " << candReReco.phi()
217  << std::endl
218  << " DeltaE = : " << deltaE << std::endl
219  << " DeltaEta = : " << deltaEta << std::endl
220  << " DeltaPhi = : " << deltaPhi << std::endl
221  << std::endl;
222  if (++npr == 5)
223  break;
224  } else {
225  std::cout << "Jet " << i << " " << candReco.et() << std::endl;
226  }
227  }
228 }
229 
231  if (!out)
232  return;
233 
234  PFBlockRef firstRef;
235 
236  assert(!cand.elementsInBlocks().empty());
237  for (unsigned i = 0; i < cand.elementsInBlocks().size(); i++) {
238  PFBlockRef blockRef = cand.elementsInBlocks()[i].first;
239 
240  if (blockRef.isNull()) {
241  cerr << "ERROR! no block ref!";
242  continue;
243  }
244 
245  if (!i) {
246  out << (*blockRef);
247  firstRef = blockRef;
248  } else if (blockRef != firstRef) {
249  cerr << "WARNING! This PFCandidate is not made from a single block" << endl;
250  }
251 
252  out << "\t" << cand.elementsInBlocks()[i].second << endl;
253  }
254 }
#define LogDebug(id)
RunNumber_t run() const
Definition: EventID.h:38
T getParameter(std::string const &) const
EventNumber_t event() const
Definition: EventID.h:40
void analyze(const edm::Event &, const edm::EventSetup &) override
T getUntrackedParameter(std::string const &, T const &) const
void beginRun(const edm::Run &r, const edm::EventSetup &c) override
static bool greaterPt(const reco::PFCandidate &a, const reco::PFCandidate &b)
double eta() const final
momentum pseudorapidity
unsigned entry_
Counter.
double pt() const final
transverse momentum
double deltaEMax_
Cuts for comparison.
Jets made from PFObjects.
Definition: PFJet.h:20
example_stream void analyze(const edm::Event &, const edm::EventSetup &) override
static const double deltaEta
Definition: CaloConstants.h:8
Checks what a re-reco changes in PFCandidates.
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
double et() const final
transverse energy
T sqrt(T t)
Definition: SSEVec.h:19
bool verbose_
verbose ?
double energy() const final
energy
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:488
bool isNull() const
Checks for null.
Definition: Ref.h:235
edm::InputTag inputTagPFJetsReReco_
bool rankByPt_
rank the candidates by Pt
edm::InputTag inputTagPFJetsReco_
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
double b
Definition: hdecay.h:118
void printElementsInBlocks(const reco::PFCandidate &cand, std::ostream &out=std::cout) const
edm::EventID id() const
Definition: EventBase.h:59
std::vector< PFJet > PFJetCollection
collection of PFJet objects
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:40
fixed size matrix
HLT enums.
edm::InputTag inputTagPFCandidatesReco_
PFCandidates in which we&#39;ll look for pile up particles.
double a
Definition: hdecay.h:119
bool printBlocks_
print the blocks associated to a given candidate ?
PFCandidateChecker(const edm::ParameterSet &)
const ElementsInBlocks & elementsInBlocks() const
Definition: PFCandidate.cc:636
double phi() const final
momentum azimuthal angle
edm::InputTag inputTagPFCandidatesReReco_
Definition: Run.h:45
void printMet(const reco::PFCandidateCollection &pfReco, const reco::PFCandidateCollection &pfReReco) const
void printJets(const reco::PFJetCollection &pfJetsReco, const reco::PFJetCollection &pfJetsReReco) const