CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
MuonMCClassifier Class Reference

#include <MuonAnalysis/MuonAssociators/src/MuonMCClassifier.cc>

Inheritance diagram for MuonMCClassifier:
edm::EDProducer edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

 MuonMCClassifier (const edm::ParameterSet &)
 
 ~MuonMCClassifier ()
 
- Public Member Functions inherited from edm::EDProducer
 EDProducer ()
 
ModuleDescription const & moduleDescription () const
 
virtual ~EDProducer ()
 
- Public Member Functions inherited from edm::ProducerBase
 ProducerBase ()
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription
const &)> 
registrationCallback () const
 used by the fwk to register list of products More...
 
virtual ~ProducerBase ()
 
- Public Member Functions inherited from edm::EDConsumerBase
 EDConsumerBase ()
 
ProductHolderIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
std::vector
< ProductHolderIndexAndSkipBit >
const & 
itemsToGetFromEvent () const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesDependentUpon (const std::string &iProcessName, std::vector< const char * > &oModuleLabels) const
 
bool registeredToConsume (ProductHolderIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)
 
virtual ~EDConsumerBase ()
 

Private Member Functions

int convertAndPush (const TrackingParticle &tp, reco::GenParticleCollection &out, const TrackingParticleRef &momRef, const edm::Handle< reco::GenParticleCollection > &genParticles) const
 
int flavour (int pdgId) const
 Returns the flavour given a pdg id code. More...
 
TrackingParticleRef getTpMother (TrackingParticleRef tp)
 
virtual void produce (edm::Event &, const edm::EventSetup &) override
 
template<typename T >
void writeValueMap (edm::Event &iEvent, const edm::Handle< edm::View< reco::Muon > > &handle, const std::vector< T > &values, const std::string &label) const
 Write a ValueMap<int> in the event. More...
 

Private Attributes

std::string associatorLabel_
 The Associations. More...
 
double decayAbsZ_
 
double decayRho_
 Cylinder to use to decide if a decay is early or late. More...
 
edm::InputTag genParticles_
 
edm::EDGetTokenT
< reco::GenParticleCollection
genParticlesToken_
 
bool hasMuonCut_
 
bool linkToGenParticles_
 Create a link to the generator level particles. More...
 
StringCutObjectSelector
< pat::Muon
muonCut_
 
edm::EDGetTokenT< edm::View
< reco::Muon > > 
muonsToken_
 The RECO objects. More...
 
edm::EDGetTokenT
< TrackingParticleCollection
trackingParticlesToken_
 The TrackingParticle objects. More...
 
MuonAssociatorByHits::MuonTrackType trackType_
 Track to use. More...
 

Additional Inherited Members

- Public Types inherited from edm::EDProducer
typedef EDProducer ModuleType
 
- Public Types inherited from edm::ProducerBase
typedef
ProductRegistryHelper::TypeLabelList 
TypeLabelList
 
- Static Public Member Functions inherited from edm::EDProducer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

CLASSIFICATION: For each RECO Muon, match to SIM particle, and then:

In any case, if the TP is not preferentially matched back to the same RECO muon, label as Ghost (flip the classification)

FLAVOUR:

Definition at line 69 of file MuonMCClassifier.cc.

Constructor & Destructor Documentation

MuonMCClassifier::MuonMCClassifier ( const edm::ParameterSet iConfig)
explicit

Definition at line 128 of file MuonMCClassifier.cc.

References edm::hlt::Exception, genParticles_, genParticlesToken_, edm::ParameterSet::getParameter(), MuonAssociatorByHits::GlobalTk, MuonAssociatorByHits::InnerTk, linkToGenParticles_, MuonAssociatorByHits::OuterTk, MuonAssociatorByHits::Segments, AlCaHLTBitMon_QueryRunRegistry::string, and trackType_.

128  :
130  hasMuonCut_(iConfig.existsAs<std::string>("muonPreselection")),
131  muonCut_(hasMuonCut_ ? iConfig.getParameter<std::string>("muonPreselection") : ""),
132  trackingParticlesToken_(consumes<TrackingParticleCollection>(iConfig.getParameter<edm::InputTag>("trackingParticles"))),
133  associatorLabel_(iConfig.getParameter< std::string >("associatorLabel")),
134  decayRho_(iConfig.getParameter<double>("decayRho")),
135  decayAbsZ_(iConfig.getParameter<double>("decayAbsZ")),
136  linkToGenParticles_(iConfig.getParameter<bool>("linkToGenParticles")),
137  genParticles_(linkToGenParticles_ ? iConfig.getParameter<edm::InputTag>("genParticles") : edm::InputTag("NONE"))
138 {
139  std::string trackType = iConfig.getParameter< std::string >("trackType");
140  if (trackType == "inner") trackType_ = MuonAssociatorByHits::InnerTk;
141  else if (trackType == "outer") trackType_ = MuonAssociatorByHits::OuterTk;
142  else if (trackType == "global") trackType_ = MuonAssociatorByHits::GlobalTk;
143  else if (trackType == "segments") trackType_ = MuonAssociatorByHits::Segments;
144  else throw cms::Exception("Configuration") << "Track type '" << trackType << "' not supported.\n";
145  if (linkToGenParticles_) {
146  genParticlesToken_ = consumes<reco::GenParticleCollection>(genParticles_);
147  }
148 
149  produces<edm::ValueMap<int> >();
150  produces<edm::ValueMap<int> >("ext");
151  produces<edm::ValueMap<int> >("flav");
152  produces<edm::ValueMap<int> >("hitsPdgId");
153  produces<edm::ValueMap<int> >("momPdgId");
154  produces<edm::ValueMap<int> >("momFlav");
155  produces<edm::ValueMap<int> >("momStatus");
156  produces<edm::ValueMap<int> >("gmomPdgId");
157  produces<edm::ValueMap<int> >("gmomFlav");
158  produces<edm::ValueMap<int> >("hmomFlav"); // heaviest mother flavour
159  produces<edm::ValueMap<int> >("tpId");
160  produces<edm::ValueMap<float> >("prodRho");
161  produces<edm::ValueMap<float> >("prodZ");
162  produces<edm::ValueMap<float> >("momRho");
163  produces<edm::ValueMap<float> >("momZ");
164  produces<edm::ValueMap<float> >("tpAssoQuality");
165  if (linkToGenParticles_) {
166  produces<reco::GenParticleCollection>("secondaries");
167  produces<edm::Association<reco::GenParticleCollection> >("toPrimaries");
168  produces<edm::Association<reco::GenParticleCollection> >("toSecondaries");
169  }
170 }
edm::InputTag genParticles_
bool linkToGenParticles_
Create a link to the generator level particles.
T getParameter(std::string const &) const
bool existsAs(std::string const &parameterName, bool trackiness=true) const
checks if a parameter exists as a given type
Definition: ParameterSet.h:184
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
std::string associatorLabel_
The Associations.
edm::EDGetTokenT< reco::GenParticleCollection > genParticlesToken_
MuonAssociatorByHits::MuonTrackType trackType_
Track to use.
double decayRho_
Cylinder to use to decide if a decay is early or late.
edm::EDGetTokenT< edm::View< reco::Muon > > muonsToken_
The RECO objects.
edm::EDGetTokenT< TrackingParticleCollection > trackingParticlesToken_
The TrackingParticle objects.
StringCutObjectSelector< pat::Muon > muonCut_
MuonMCClassifier::~MuonMCClassifier ( )

Definition at line 172 of file MuonMCClassifier.cc.

173 {
174 }

Member Function Documentation

int MuonMCClassifier::convertAndPush ( const TrackingParticle tp,
reco::GenParticleCollection out,
const TrackingParticleRef momRef,
const edm::Handle< reco::GenParticleCollection > &  genParticles 
) const
private

Convert TrackingParticle into GenParticle, save into output collection, if mother is primary set reference to it, return index in output collection

Definition at line 484 of file MuonMCClassifier.cc.

References TrackingParticle::charge(), edm::hlt::Exception, genParticles_, edm::HandleBase::id(), edm::Ref< C, T, F >::id(), edm::Ref< C, T, F >::isNonnull(), TrackingParticle::p4(), TrackingParticle::pdgId(), TrackingParticle::status(), and TrackingParticle::vertex().

Referenced by produce().

487  {
488  out.push_back(reco::GenParticle(tp.charge(), tp.p4(), tp.vertex(), tp.pdgId(), tp.status(), true));
489  if (simMom.isNonnull() && !simMom->genParticles().empty()) {
490  if (genParticles.id() != simMom->genParticles().id()) {
491  throw cms::Exception("Configuration") << "Product ID mismatch between the genParticle collection (" << genParticles_ << ", id " << genParticles.id() << ") and the references in the TrackingParticles (id " << simMom->genParticles().id() << ")\n";
492  }
493  out.back().addMother(simMom->genParticles()[0]);
494  }
495  return out.size()-1;
496 }
edm::InputTag genParticles_
ProductID id() const
Definition: HandleBase.cc:15
int pdgId() const
PDG ID.
int status() const
Status word.
float charge() const
Gives charge in unit of quark charge (should be 3 time the abaove)
tuple out
Definition: dbtoconf.py:99
const LorentzVector & p4() const
Four-momentum Lorentz vector. Note this is taken from the first SimTrack only.
Point vertex() const
ProductIndex id() const
Definition: ProductID.h:38
int MuonMCClassifier::flavour ( int  pdgId) const
private

Returns the flavour given a pdg id code.

Definition at line 468 of file MuonMCClassifier.cc.

References funct::abs().

Referenced by cuy.graphElement::__init__(), and produce().

468  {
469  int flav = abs(pdgId);
470  // for quarks, leptons and bosons except gluons, take their pdgId
471  // muons and taus have themselves as flavour
472  if (flav <= 37 && flav != 21) return flav;
473  // look for barions
474  int bflav = ((flav / 1000) % 10);
475  if (bflav != 0) return bflav;
476  // look for mesons
477  int mflav = ((flav / 100) % 10);
478  if (mflav != 0) return mflav;
479  return 0;
480 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
TrackingParticleRef MuonMCClassifier::getTpMother ( TrackingParticleRef  tp)
inlineprivate

Definition at line 111 of file MuonMCClassifier.cc.

References edm::Ref< C, T, F >::isNonnull().

Referenced by produce().

111  {
112  if (tp.isNonnull() && tp->parentVertex().isNonnull() && !tp->parentVertex()->sourceTracks().empty()) {
113  return tp->parentVertex()->sourceTracks()[0];
114  } else {
115  return TrackingParticleRef();
116  }
117  }
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:250
edm::Ref< TrackingParticleCollection > TrackingParticleRef
void MuonMCClassifier::produce ( edm::Event iEvent,
const edm::EventSetup iSetup 
)
overrideprivatevirtual

Implements edm::EDProducer.

Definition at line 177 of file MuonMCClassifier.cc.

References funct::abs(), MuonAssociatorByHits::associateMuons(), associatorLabel_, edm::RefToBaseVector< T >::begin(), convertAndPush(), decayAbsZ_, decayRho_, edm::RefToBaseVector< T >::end(), edm::hlt::Exception, edm::helper::Filler< Map >::fill(), spr::find(), flavour(), EgammaValidation_cff::genp, genParticleCandidates2GenParticles_cfi::genParticles, genParticles_, genParticlesToken_, edm::EventSetup::get(), edm::Event::getByToken(), getTpMother(), MuonAssociatorByHits::GlobalTk, hasMuonCut_, i, edm::HandleBase::id(), edm::Ref< C, T, F >::id(), edm::helper::Filler< Map >::insert(), reco::Muon::isGlobalMuon(), edm::Ref< C, T, F >::isNonnull(), edm::Ref< C, T, F >::key(), linkToGenParticles_, match(), RPCpg::mu, muonCut_, patZpeak::muons, muonsToken_, n, MuonAssociatorByHits::OuterTk, edm::ESHandle< class >::product(), edm::RefToBaseVector< T >::push_back(), edm::RefVector< C, T, F >::push_back(), edm::Event::put(), trackingParticlesFastSim_cfi::trackingParticles, trackingParticlesToken_, trackType_, and writeValueMap().

178 {
179  edm::LogVerbatim("MuonMCClassifier") <<"\n sono in MuonMCClassifier !";
180 
182  iEvent.getByToken(muonsToken_, muons);
183 
185  iEvent.getByToken(trackingParticlesToken_, trackingParticles);
186 
188  if (linkToGenParticles_) {
189  iEvent.getByToken(genParticlesToken_, genParticles);
190  }
191 
192  edm::ESHandle<TrackAssociatorBase> associatorBase;
193  iSetup.get<TrackAssociatorRecord>().get(associatorLabel_, associatorBase);
194  const MuonAssociatorByHits * assoByHits = dynamic_cast<const MuonAssociatorByHits *>(associatorBase.product());
195  if (assoByHits == 0) throw cms::Exception("Configuration") << "The Track Associator with label '" << associatorLabel_ << "' is not a MuonAssociatorByHits.\n";
196 
199  edm::LogVerbatim("MuonMCClassifier") <<"\n ***************************************************************** ";
200  edm::LogVerbatim("MuonMCClassifier") << " RECO MUON association, type: "<< trackType_;
201  edm::LogVerbatim("MuonMCClassifier") << " ***************************************************************** \n";
202 
204  if (!hasMuonCut_) {
205  // all muons
206  for (size_t i = 0, n = muons->size(); i < n; ++i) {
207  edm::RefToBase<reco::Muon> rmu = muons->refAt(i);
208  selMuons.push_back(rmu);
209  }
210  } else {
211  // filter, fill refvectors, associate
212  // I pass through pat::Muon so that I can access muon id selectors
213  for (size_t i = 0, n = muons->size(); i < n; ++i) {
214  edm::RefToBase<reco::Muon> rmu = muons->refAt(i);
215  if (muonCut_(pat::Muon(rmu))) selMuons.push_back(rmu);
216  }
217  }
218 
220  for (size_t i = 0, n = trackingParticles->size(); i < n; ++i) {
221  allTPs.push_back(TrackingParticleRef(trackingParticles,i));
222  }
223 
224  assoByHits->associateMuons(recSimColl, simRecColl, selMuons, trackType_, allTPs, &iEvent, &iSetup);
225 
226  // for global muons without hits on muon detectors, look at the linked standalone muon
230  edm::LogVerbatim("MuonMCClassifier") <<"\n ***************************************************************** ";
231  edm::LogVerbatim("MuonMCClassifier") << " STANDALONE (UpdAtVtx) MUON association ";
232  edm::LogVerbatim("MuonMCClassifier") << " ***************************************************************** \n";
233  assoByHits->associateMuons(UpdSTA_recSimColl, UpdSTA_simRecColl, selMuons, MuonAssociatorByHits::OuterTk,
234  allTPs, &iEvent, &iSetup);
235  }
236 
237  typedef MuonAssociatorByHits::MuonToSimCollection::const_iterator r2s_it;
238  typedef MuonAssociatorByHits::SimToMuonCollection::const_iterator s2r_it;
239 
240  size_t nmu = muons->size();
241  edm::LogVerbatim("MuonMCClassifier") <<"\n There are "<<nmu<<" reco::Muons.";
242 
243  std::vector<int> classif(nmu, 0), ext(nmu, 0);
244  std::vector<int> hitsPdgId(nmu, 0), momPdgId(nmu, 0), gmomPdgId(nmu, 0), momStatus(nmu, 0);
245  std::vector<int> flav(nmu, 0), momFlav(nmu, 0), gmomFlav(nmu, 0), hmomFlav(nmu, 0);
246  std::vector<int> tpId(nmu, -1);
247  std::vector<float> prodRho(nmu, 0.0), prodZ(nmu, 0.0), momRho(nmu, 0.0), momZ(nmu, 0.0);
248  std::vector<float> tpAssoQuality(nmu, -1);
249 
250  std::auto_ptr<reco::GenParticleCollection> secondaries; // output collection of secondary muons
251  std::map<TrackingParticleRef, int> tpToSecondaries; // map from tp to (index+1) in output collection
252  std::vector<int> muToPrimary(nmu, -1), muToSecondary(nmu, -1); // map from input into (index) in output, -1 for null
253  if (linkToGenParticles_) secondaries.reset(new reco::GenParticleCollection());
254 
255  for(size_t i = 0; i < nmu; ++i) {
256  edm::LogVerbatim("MuonMCClassifier") <<"\n reco::Muons # "<<i;
257  edm::RefToBase<reco::Muon> mu = muons->refAt(i);
258  if (hasMuonCut_ && (std::find(selMuons.begin(), selMuons.end(), mu) == selMuons.end()) ) {
259  edm::LogVerbatim("MuonMCClassifier") <<"\t muon didn't pass the selection. classified as -99 and skipped";
260  classif[i] = -99; continue;
261  }
262 
264  edm::RefToBase<reco::Muon> muMatchBack;
265  r2s_it match = recSimColl.find(mu);
266  s2r_it matchback;
267  if (match != recSimColl.end()) {
268  edm::LogVerbatim("MuonMCClassifier") <<"\t RtS matched Ok...";
269  // match->second is vector, front is first element, first is the ref (second would be the quality)
270  tp = match->second.front().first;
271  tpId[i] = tp.isNonnull() ? tp.key() : -1; // we check, even if null refs should not appear here at all
272  tpAssoQuality[i] = match->second.front().second;
273  s2r_it matchback = simRecColl.find(tp);
274  if (matchback != simRecColl.end()) {
275  muMatchBack = matchback->second.front().first;
276  } else {
277  edm::LogWarning("MuonMCClassifier") << "\n***WARNING: This I do NOT understand: why no match back? *** \n";
278  }
279  } else if ((trackType_ == MuonAssociatorByHits::GlobalTk) &&
280  mu->isGlobalMuon()) {
281  // perform a second attempt, matching with the standalone muon
282  r2s_it matchSta = UpdSTA_recSimColl.find(mu);
283  if (matchSta != UpdSTA_recSimColl.end()) {
284  edm::LogVerbatim("MuonMCClassifier") <<"\t RtS matched Ok... from the UpdSTA_recSimColl ";
285  tp = matchSta->second.front().first;
286  tpId[i] = tp.isNonnull() ? tp.key() : -1; // we check, even if null refs should not appear here at all
287  tpAssoQuality[i] = matchSta->second.front().second;
288  s2r_it matchback = UpdSTA_simRecColl.find(tp);
289  if (matchback != UpdSTA_simRecColl.end()) {
290  muMatchBack = matchback->second.front().first;
291  } else {
292  edm::LogWarning("MuonMCClassifier") << "\n***WARNING: This I do NOT understand: why no match back in UpdSTA? *** \n";
293  }
294  }
295  }
296  if (tp.isNonnull()) {
297  bool isGhost = muMatchBack != mu;
298  if (isGhost) edm::LogVerbatim("MuonMCClassifier") <<"\t This seems a GHOST ! classif[i] will be < 0";
299 
300  hitsPdgId[i] = tp->pdgId();
301  prodRho[i] = tp->vertex().Rho();
302  prodZ[i] = tp->vertex().Z();
303  edm::LogVerbatim("MuonMCClassifier") <<"\t TP pdgId = "<<hitsPdgId[i] << ", vertex rho = " << prodRho[i] << ", z = " << prodZ[i];
304 
305  // Try to extract mother and grand mother of this muon.
306  // Unfortunately, SIM and GEN histories require diffent code :-(
307  if (!tp->genParticles().empty()) { // Muon is in GEN
308  reco::GenParticleRef genp = tp->genParticles()[0];
309  reco::GenParticleRef genMom = genp->numberOfMothers() > 0 ? genp->motherRef() : reco::GenParticleRef();
310  if (genMom.isNonnull()) {
311  momPdgId[i] = genMom->pdgId();
312  momStatus[i] = genMom->status();
313  momRho[i] = genMom->vertex().Rho(); momZ[i] = genMom->vz();
314  edm::LogVerbatim("MuonMCClassifier") << "\t Particle pdgId = "<<hitsPdgId[i] << " produced at rho = " << prodRho[i] << ", z = " << prodZ[i] << ", has GEN mother pdgId = " << momPdgId[i];
315  reco::GenParticleRef genGMom = genMom->numberOfMothers() > 0 ? genMom->motherRef() : reco::GenParticleRef();
316  if (genGMom.isNonnull()) {
317  gmomPdgId[i] = genGMom->pdgId();
318  edm::LogVerbatim("MuonMCClassifier") << "\t\t mother prod. vertex rho = " << momRho[i] << ", z = " << momZ[i] << ", grand-mom pdgId = " << gmomPdgId[i];
319  }
320  // in this case, we might want to know the heaviest mom flavour
321  for (reco::GenParticleRef nMom = genMom;
322  nMom.isNonnull() && abs(nMom->pdgId()) >= 100; // stop when we're no longer looking at hadrons or mesons
323  nMom = nMom->numberOfMothers() > 0 ? nMom->motherRef() : reco::GenParticleRef()) {
324  int flav = flavour(nMom->pdgId());
325  if (hmomFlav[i] < flav) hmomFlav[i] = flav;
326  edm::LogVerbatim("MuonMCClassifier") << "\t\t backtracking flavour: mom pdgId = "<<nMom->pdgId()<< ", flavour = " << flav << ", heaviest so far = " << hmomFlav[i];
327  }
328  }
329  } else { // Muon is in SIM Only
330  TrackingParticleRef simMom = getTpMother(tp);
331  if (simMom.isNonnull()) {
332  momPdgId[i] = simMom->pdgId();
333  momRho[i] = simMom->vertex().Rho();
334  momZ[i] = simMom->vertex().Z();
335  edm::LogVerbatim("MuonMCClassifier") << "\t Particle pdgId = "<<hitsPdgId[i] << " produced at rho = " << prodRho[i] << ", z = " << prodZ[i] <<
336  ", has SIM mother pdgId = " << momPdgId[i] << " produced at rho = " << simMom->vertex().Rho() << ", z = " << simMom->vertex().Z();
337  if (!simMom->genParticles().empty()) {
338  momStatus[i] = simMom->genParticles()[0]->status();
339  reco::GenParticleRef genGMom = (simMom->genParticles()[0]->numberOfMothers() > 0 ? simMom->genParticles()[0]->motherRef() : reco::GenParticleRef());
340  if (genGMom.isNonnull()) gmomPdgId[i] = genGMom->pdgId();
341  edm::LogVerbatim("MuonMCClassifier") << "\t\t SIM mother is in GEN (status " << momStatus[i] << "), grand-mom id = " << gmomPdgId[i];
342  } else {
343  momStatus[i] = -1;
344  TrackingParticleRef simGMom = getTpMother(simMom);
345  if (simGMom.isNonnull()) gmomPdgId[i] = simGMom->pdgId();
346  edm::LogVerbatim("MuonMCClassifier") << "\t\t SIM mother is in SIM only, grand-mom id = " << gmomPdgId[i];
347  }
348  } else {
349  edm::LogVerbatim("MuonMCClassifier") << "\t Particle pdgId = "<<hitsPdgId[i] << " produced at rho = " << prodRho[i] << ", z = " << prodZ[i] << ", has NO mother!";
350  }
351  }
352  momFlav[i] = flavour(momPdgId[i]);
353  gmomFlav[i] = flavour(gmomPdgId[i]);
354 
355  // Check first IF this is a muon at all
356  if (abs(tp->pdgId()) != 13) {
357  classif[i] = isGhost ? -1 : 1;
358  ext[i] = isGhost ? -1 : 1;
359  edm::LogVerbatim("MuonMCClassifier") <<"\t This is not a muon. Sorry. classif[i] = " << classif[i];
360  continue;
361  }
362 
363  // Is this SIM muon also a GEN muon, with a mother?
364  if (!tp->genParticles().empty() && (momPdgId[i] != 0)) {
365  if (abs(momPdgId[i]) < 100 && (abs(momPdgId[i]) != 15)) {
366  classif[i] = isGhost ? -4 : 4;
367  flav[i] = (abs(momPdgId[i]) == 15 ? 15 : 13);
368  edm::LogVerbatim("MuonMCClassifier") <<"\t This seems PRIMARY MUON ! classif[i] = " << classif[i];
369  ext[i] = 10;
370  } else if (momFlav[i] == 4 || momFlav[i] == 5 || momFlav[i] == 15) {
371  classif[i] = isGhost ? -3 : 3;
372  flav[i] = momFlav[i];
373  if (momFlav[i] == 15) ext[i] = 9; // tau->mu
374  else if (momFlav[i] == 5) ext[i] = 8; // b->mu
375  else if (hmomFlav[i] == 5) ext[i] = 7; // b->c->mu
376  else ext[i] = 6; // c->mu
377  edm::LogVerbatim("MuonMCClassifier") <<"\t This seems HEAVY FLAVOUR ! classif[i] = " << classif[i];
378  } else {
379  classif[i] = isGhost ? -2 : 2;
380  flav[i] = momFlav[i];
381  edm::LogVerbatim("MuonMCClassifier") <<"\t This seems LIGHT FLAVOUR ! classif[i] = " << classif[i];
382  }
383  } else {
384  classif[i] = isGhost ? -2 : 2;
385  flav[i] = momFlav[i];
386  edm::LogVerbatim("MuonMCClassifier") <<"\t This seems LIGHT FLAVOUR ! classif[i] = " << classif[i];
387  }
388  // extended classification
389  if (momPdgId[i] == 0) ext[i] = 2; // if it has no mom, it's not a primary particle so it won't be in ppMuX
390  else if (abs(momPdgId[i]) < 100) ext[i] = (momFlav[i] == 15 ? 9 : 10); // primary mu, or tau->mu
391  else if (momFlav[i] == 5) ext[i] = 8; // b->mu
392  else if (momFlav[i] == 4) ext[i] = (hmomFlav[i] == 5 ? 7 : 6); // b->c->mu and c->mu
393  else if (momStatus[i] != -1) { // primary light particle
394  int id = abs(momPdgId[i]);
395  if (id != /*pi+*/211 && id != /*K+*/321 && id != 130 /*K0L*/) ext[i] = 5; // other light particle, possibly short-lived
396  else if (prodRho[i] < decayRho_ && abs(prodZ[i]) < decayAbsZ_) ext[i] = 4; // decay a la ppMuX (primary pi/K within a cylinder)
397  else ext[i] = 3; // late decay that wouldn't be in ppMuX
398  } else ext[i] = 2; // decay of non-primary particle, would not be in ppMuX
399  if (isGhost) ext[i] = -ext[i];
400 
401  if (linkToGenParticles_ && abs(ext[i]) >= 2) {
402  // Link to the genParticle if possible, but not decays in flight (in ppMuX they're in GEN block, but they have wrong parameters)
403  if (!tp->genParticles().empty() && abs(ext[i]) >= 5) {
404  if (genParticles.id() != tp->genParticles().id()) {
405  throw cms::Exception("Configuration") << "Product ID mismatch between the genParticle collection (" << genParticles_ << ", id " << genParticles.id() << ") and the references in the TrackingParticles (id " << tp->genParticles().id() << ")\n";
406  }
407  muToPrimary[i] = tp->genParticles()[0].key();
408  } else {
409  // Don't put the same trackingParticle twice!
410  int &indexPlus1 = tpToSecondaries[tp]; // will create a 0 if the tp is not in the list already
411  if (indexPlus1 == 0) indexPlus1 = convertAndPush(*tp, *secondaries, getTpMother(tp), genParticles) + 1;
412  muToSecondary[i] = indexPlus1 - 1;
413  }
414  }
415  edm::LogVerbatim("MuonMCClassifier") <<"\t Extended classification code = " << ext[i];
416  }
417  }
418 
419  writeValueMap(iEvent, muons, classif, "");
420  writeValueMap(iEvent, muons, ext, "ext");
421  writeValueMap(iEvent, muons, flav, "flav");
422  writeValueMap(iEvent, muons, tpId, "tpId");
423  writeValueMap(iEvent, muons, hitsPdgId, "hitsPdgId");
424  writeValueMap(iEvent, muons, momPdgId, "momPdgId");
425  writeValueMap(iEvent, muons, momStatus, "momStatus");
426  writeValueMap(iEvent, muons, momFlav, "momFlav");
427  writeValueMap(iEvent, muons, gmomPdgId, "gmomPdgId");
428  writeValueMap(iEvent, muons, gmomFlav, "gmomFlav");
429  writeValueMap(iEvent, muons, hmomFlav, "hmomFlav");
430  writeValueMap(iEvent, muons, prodRho, "prodRho");
431  writeValueMap(iEvent, muons, prodZ, "prodZ");
432  writeValueMap(iEvent, muons, momRho, "momRho");
433  writeValueMap(iEvent, muons, momZ, "momZ");
434  writeValueMap(iEvent, muons, tpAssoQuality, "tpAssoQuality");
435 
436  if (linkToGenParticles_) {
437  edm::OrphanHandle<reco::GenParticleCollection> secHandle = iEvent.put(secondaries, "secondaries");
438  edm::RefProd<reco::GenParticleCollection> priRP(genParticles);
440  std::auto_ptr<edm::Association<reco::GenParticleCollection> > outPri(new edm::Association<reco::GenParticleCollection>(priRP));
441  std::auto_ptr<edm::Association<reco::GenParticleCollection> > outSec(new edm::Association<reco::GenParticleCollection>(secRP));
442  edm::Association<reco::GenParticleCollection>::Filler fillPri(*outPri), fillSec(*outSec);
443  fillPri.insert(muons, muToPrimary.begin(), muToPrimary.end());
444  fillSec.insert(muons, muToSecondary.begin(), muToSecondary.end());
445  fillPri.fill(); fillSec.fill();
446  iEvent.put(outPri, "toPrimaries");
447  iEvent.put(outSec, "toSecondaries");
448  }
449 }
edm::InputTag genParticles_
bool linkToGenParticles_
Create a link to the generator level particles.
std::vector< GenParticle > GenParticleCollection
collection of GenParticles
int i
Definition: DBlmapReader.cc:9
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:250
ProductID id() const
Definition: HandleBase.cc:15
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:446
void associateMuons(MuonToSimCollection &recoToSim, SimToMuonCollection &simToReco, const edm::RefToBaseVector< reco::Muon > &, MuonTrackType, const edm::RefVector< TrackingParticleCollection > &, const edm::Event *event=0, const edm::EventSetup *setup=0) const
edm::Ref< GenParticleCollection > GenParticleRef
persistent reference to a GenParticle
const_iterator end() const
bool isGlobalMuon() const
Definition: Muon.h:218
int convertAndPush(const TrackingParticle &tp, reco::GenParticleCollection &out, const TrackingParticleRef &momRef, const edm::Handle< reco::GenParticleCollection > &genParticles) const
key_type key() const
Accessor for product key.
Definition: Ref.h:266
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
ProductID id() const
Accessor for product ID.
Definition: Ref.h:256
std::map< TrackingParticleRef, std::vector< std::pair< edm::RefToBase< reco::Muon >, double > > > SimToMuonCollection
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:113
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
TrackingParticleRef getTpMother(TrackingParticleRef tp)
const int mu
Definition: Constants.h:22
std::string associatorLabel_
The Associations.
tuple genp
produce generated paricles in acceptance #
edm::EDGetTokenT< reco::GenParticleCollection > genParticlesToken_
MuonAssociatorByHits::MuonTrackType trackType_
Track to use.
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:86
void writeValueMap(edm::Event &iEvent, const edm::Handle< edm::View< reco::Muon > > &handle, const std::vector< T > &values, const std::string &label) const
Write a ValueMap&lt;int&gt; in the event.
const_iterator begin() const
double decayRho_
Cylinder to use to decide if a decay is early or late.
edm::EDGetTokenT< edm::View< reco::Muon > > muonsToken_
The RECO objects.
tuple muons
Definition: patZpeak.py:38
void push_back(const RefToBase< T > &)
void push_back(value_type const &ref)
Add a Ref&lt;C, T&gt; to the RefVector.
Definition: RefVector.h:64
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:6
int flavour(int pdgId) const
Returns the flavour given a pdg id code.
Analysis-level muon class.
Definition: Muon.h:50
std::map< edm::RefToBase< reco::Muon >, std::vector< std::pair< TrackingParticleRef, double > >, RefToBaseSort > MuonToSimCollection
edm::Ref< TrackingParticleCollection > TrackingParticleRef
edm::EDGetTokenT< TrackingParticleCollection > trackingParticlesToken_
The TrackingParticle objects.
StringCutObjectSelector< pat::Muon > muonCut_
template<typename T >
void MuonMCClassifier::writeValueMap ( edm::Event iEvent,
const edm::Handle< edm::View< reco::Muon > > &  handle,
const std::vector< T > &  values,
const std::string &  label 
) const
private

Write a ValueMap<int> in the event.

Definition at line 453 of file MuonMCClassifier.cc.

References edm::helper::Filler< Map >::fill(), patZpeak::handle, edm::helper::Filler< Map >::insert(), and edm::Event::put().

Referenced by produce().

457 {
458  using namespace edm;
459  using namespace std;
460  auto_ptr<ValueMap<T> > valMap(new ValueMap<T>());
461  typename edm::ValueMap<T>::Filler filler(*valMap);
462  filler.insert(handle, values.begin(), values.end());
463  filler.fill();
464  iEvent.put(valMap, label);
465 }
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:113

Member Data Documentation

std::string MuonMCClassifier::associatorLabel_
private

The Associations.

Definition at line 91 of file MuonMCClassifier.cc.

Referenced by produce().

double MuonMCClassifier::decayAbsZ_
private

Definition at line 94 of file MuonMCClassifier.cc.

Referenced by produce().

double MuonMCClassifier::decayRho_
private

Cylinder to use to decide if a decay is early or late.

Definition at line 94 of file MuonMCClassifier.cc.

Referenced by produce().

edm::InputTag MuonMCClassifier::genParticles_
private

Definition at line 98 of file MuonMCClassifier.cc.

Referenced by convertAndPush(), MuonMCClassifier(), and produce().

edm::EDGetTokenT<reco::GenParticleCollection> MuonMCClassifier::genParticlesToken_
private

Definition at line 99 of file MuonMCClassifier.cc.

Referenced by MuonMCClassifier(), and produce().

bool MuonMCClassifier::hasMuonCut_
private

A preselection cut for the muon. I pass through pat::Muon so that I can access muon id selectors

Definition at line 81 of file MuonMCClassifier.cc.

Referenced by produce().

bool MuonMCClassifier::linkToGenParticles_
private

Create a link to the generator level particles.

Definition at line 97 of file MuonMCClassifier.cc.

Referenced by MuonMCClassifier(), and produce().

StringCutObjectSelector<pat::Muon> MuonMCClassifier::muonCut_
private

Definition at line 82 of file MuonMCClassifier.cc.

Referenced by produce().

edm::EDGetTokenT<edm::View<reco::Muon> > MuonMCClassifier::muonsToken_
private

The RECO objects.

Definition at line 77 of file MuonMCClassifier.cc.

Referenced by produce().

edm::EDGetTokenT<TrackingParticleCollection> MuonMCClassifier::trackingParticlesToken_
private

The TrackingParticle objects.

Definition at line 88 of file MuonMCClassifier.cc.

Referenced by produce().

MuonAssociatorByHits::MuonTrackType MuonMCClassifier::trackType_
private

Track to use.

Definition at line 85 of file MuonMCClassifier.cc.

Referenced by MuonMCClassifier(), and produce().