CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HLTScoutingPFProducer.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: HLTrigger/JetMET
4 // Class: HLTScoutingPFProducer
5 //
11 //
12 // Original Author: Dustin James Anderson
13 // Created: Fri, 12 Jun 2015 15:49:20 GMT
14 //
15 //
16 
17 // system include files
18 #include <memory>
19 
20 // user include files
26 
35 
39 
41 
43 public:
44  explicit HLTScoutingPFProducer(const edm::ParameterSet &);
45  ~HLTScoutingPFProducer() override;
46 
47  static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
48 
49 private:
50  void produce(edm::StreamID sid, edm::Event &iEvent, edm::EventSetup const &setup) const final;
51 
58 
59  const double pfJetPtCut;
60  const double pfJetEtaCut;
61  const double pfCandidatePtCut;
62 
63  const bool doJetTags;
64  const bool doCandidates;
65  const bool doMet;
66 };
67 
68 //
69 // constructors and destructor
70 //
72  : pfJetCollection_(consumes<reco::PFJetCollection>(iConfig.getParameter<edm::InputTag>("pfJetCollection"))),
73  pfJetTagCollection_(consumes<reco::JetTagCollection>(iConfig.getParameter<edm::InputTag>("pfJetTagCollection"))),
75  consumes<reco::PFCandidateCollection>(iConfig.getParameter<edm::InputTag>("pfCandidateCollection"))),
76  vertexCollection_(consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertexCollection"))),
77  metCollection_(consumes<reco::PFMETCollection>(iConfig.getParameter<edm::InputTag>("metCollection"))),
78  rho_(consumes<double>(iConfig.getParameter<edm::InputTag>("rho"))),
79  pfJetPtCut(iConfig.getParameter<double>("pfJetPtCut")),
80  pfJetEtaCut(iConfig.getParameter<double>("pfJetEtaCut")),
81  pfCandidatePtCut(iConfig.getParameter<double>("pfCandidatePtCut")),
82  doJetTags(iConfig.getParameter<bool>("doJetTags")),
83  doCandidates(iConfig.getParameter<bool>("doCandidates")),
84  doMet(iConfig.getParameter<bool>("doMet")) {
85  //register products
86  produces<ScoutingPFJetCollection>();
87  produces<ScoutingParticleCollection>();
88  produces<double>("rho");
89  produces<double>("pfMetPt");
90  produces<double>("pfMetPhi");
91 }
92 
94 
95 // ------------ method called to produce the data ------------
97  using namespace edm;
98 
99  //get vertices
101  std::unique_ptr<ScoutingVertexCollection> outVertices(new ScoutingVertexCollection());
102  if (iEvent.getByToken(vertexCollection_, vertexCollection)) {
103  for (auto &vtx : *vertexCollection) {
104  outVertices->emplace_back(vtx.x(),
105  vtx.y(),
106  vtx.z(),
107  vtx.zError(),
108  vtx.xError(),
109  vtx.yError(),
110  vtx.tracksSize(),
111  vtx.chi2(),
112  vtx.ndof(),
113  vtx.isValid());
114  }
115  }
116 
117  //get rho
119  std::unique_ptr<double> outRho(new double(-999));
120  if (iEvent.getByToken(rho_, rho)) {
121  outRho.reset(new double(*rho));
122  }
123 
124  //get MET
126  std::unique_ptr<double> outMetPt(new double(-999));
127  std::unique_ptr<double> outMetPhi(new double(-999));
128  if (doMet && iEvent.getByToken(metCollection_, metCollection)) {
129  outMetPt.reset(new double(metCollection->front().pt()));
130  outMetPhi.reset(new double(metCollection->front().phi()));
131  }
132 
133  //get PF candidates
135  std::unique_ptr<ScoutingParticleCollection> outPFCandidates(new ScoutingParticleCollection());
136  if (doCandidates && iEvent.getByToken(pfCandidateCollection_, pfCandidateCollection)) {
137  for (auto &cand : *pfCandidateCollection) {
138  if (cand.pt() > pfCandidatePtCut) {
139  int vertex_index = -1;
140  int index_counter = 0;
141  double dr2 = 0.0001;
142  for (auto &vtx : *outVertices) {
143  double tmp_dr2 = pow(vtx.x() - cand.vx(), 2) + pow(vtx.y() - cand.vy(), 2) + pow(vtx.z() - cand.vz(), 2);
144  if (tmp_dr2 < dr2) {
145  dr2 = tmp_dr2;
146  vertex_index = index_counter;
147  }
148  if (dr2 == 0.0)
149  break;
150  ++index_counter;
151  }
152  outPFCandidates->emplace_back(cand.pt(), cand.eta(), cand.phi(), cand.mass(), cand.pdgId(), vertex_index);
153  }
154  }
155  }
156 
157  //get PF jets
159  std::unique_ptr<ScoutingPFJetCollection> outPFJets(new ScoutingPFJetCollection());
160  if (iEvent.getByToken(pfJetCollection_, pfJetCollection)) {
161  //get PF jet tags
163  bool haveJetTags = false;
164  if (doJetTags && iEvent.getByToken(pfJetTagCollection_, pfJetTagCollection)) {
165  haveJetTags = true;
166  }
167 
168  for (auto &jet : *pfJetCollection) {
169  if (jet.pt() < pfJetPtCut || fabs(jet.eta()) > pfJetEtaCut)
170  continue;
171  //find the jet tag corresponding to the jet
172  float tagValue = -20;
173  float minDR2 = 0.01;
174  if (haveJetTags) {
175  for (auto &tag : *pfJetTagCollection) {
176  float dR2 = reco::deltaR2(jet, *(tag.first));
177  if (dR2 < minDR2) {
178  minDR2 = dR2;
179  tagValue = tag.second;
180  }
181  }
182  }
183  //get the PF constituents of the jet
184  std::vector<int> candIndices;
185  if (doCandidates) {
186  for (auto &cand : jet.getPFConstituents()) {
187  if (cand->pt() > pfCandidatePtCut) {
188  //search for the candidate in the collection
189  float minDR2 = 0.0001;
190  int matchIndex = -1;
191  int outIndex = 0;
192  for (auto &outCand : *outPFCandidates) {
193  float dR2 = pow(cand->eta() - outCand.eta(), 2) + pow(cand->phi() - outCand.phi(), 2);
194  if (dR2 < minDR2) {
195  minDR2 = dR2;
196  matchIndex = outIndex;
197  }
198  if (minDR2 == 0) {
199  break;
200  }
201  outIndex++;
202  }
203  candIndices.push_back(matchIndex);
204  }
205  }
206  }
207  outPFJets->emplace_back(jet.pt(),
208  jet.eta(),
209  jet.phi(),
210  jet.mass(),
211  jet.jetArea(),
212  jet.chargedHadronEnergy(),
213  jet.neutralHadronEnergy(),
214  jet.photonEnergy(),
215  jet.electronEnergy(),
216  jet.muonEnergy(),
217  jet.HFHadronEnergy(),
218  jet.HFEMEnergy(),
219  jet.chargedHadronMultiplicity(),
220  jet.neutralHadronMultiplicity(),
221  jet.photonMultiplicity(),
222  jet.electronMultiplicity(),
223  jet.muonMultiplicity(),
224  jet.HFHadronMultiplicity(),
225  jet.HFEMMultiplicity(),
226  jet.hoEnergy(),
227  tagValue,
228  0.0,
229  candIndices);
230  }
231  }
232 
233  //put output
234  iEvent.put(std::move(outPFCandidates));
235  iEvent.put(std::move(outPFJets));
236  iEvent.put(std::move(outRho), "rho");
237  iEvent.put(std::move(outMetPt), "pfMetPt");
238  iEvent.put(std::move(outMetPhi), "pfMetPhi");
239 }
240 
241 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
244  desc.add<edm::InputTag>("pfJetCollection", edm::InputTag("hltAK4PFJets"));
245  desc.add<edm::InputTag>("pfJetTagCollection", edm::InputTag("hltDeepCombinedSecondaryVertexBJetTagsPF"));
246  desc.add<edm::InputTag>("pfCandidateCollection", edm::InputTag("hltParticleFlow"));
247  desc.add<edm::InputTag>("vertexCollection", edm::InputTag("hltPixelVertices"));
248  desc.add<edm::InputTag>("metCollection", edm::InputTag("hltPFMETProducer"));
249  desc.add<edm::InputTag>("rho", edm::InputTag("hltFixedGridRhoFastjetAll"));
250  desc.add<double>("pfJetPtCut", 20.0);
251  desc.add<double>("pfJetEtaCut", 3.0);
252  desc.add<double>("pfCandidatePtCut", 0.6);
253  desc.add<bool>("doJetTags", true);
254  desc.add<bool>("doCandidates", true);
255  desc.add<bool>("doMet", true);
256  descriptions.add("hltScoutingPFProducer", desc);
257 }
258 
259 // declare this class as a framework plugin
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:131
const edm::EDGetTokenT< reco::PFMETCollection > metCollection_
~HLTScoutingPFProducer() override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
JetFloatAssociation::Container JetTagCollection
Definition: JetTag.h:17
const edm::EDGetTokenT< reco::PFCandidateCollection > pfCandidateCollection_
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
const edm::EDGetTokenT< reco::VertexCollection > vertexCollection_
void produce(edm::StreamID sid, edm::Event &iEvent, edm::EventSetup const &setup) const final
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
std::vector< ScoutingVertex > ScoutingVertexCollection
std::vector< ScoutingParticle > ScoutingParticleCollection
const edm::EDGetTokenT< reco::PFJetCollection > pfJetCollection_
HLTScoutingPFProducer(const edm::ParameterSet &)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
const edm::EDGetTokenT< reco::JetTagCollection > pfJetTagCollection_
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16
const edm::EDGetTokenT< double > rho_
void add(std::string const &label, ParameterSetDescription const &psetDescription)
std::vector< PFJet > PFJetCollection
collection of PFJet objects
fixed size matrix
HLT enums.
std::vector< ScoutingPFJet > ScoutingPFJetCollection
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:30
def move(src, dest)
Definition: eostools.py:511
Collection of PF MET.