CMS 3D CMS Logo

List of all members | Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes
GenHFHadronMatcher Class Reference

Finds the origin of each heavy flavour hadron and associated jets to it. More...

#include <GenHFHadronMatcher.cc>

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

Public Member Functions

 GenHFHadronMatcher (const edm::ParameterSet &)
 constructor initialising producer products and config parameters More...
 
 ~GenHFHadronMatcher ()
 
- Public Member Functions inherited from edm::EDProducer
 EDProducer ()
 
ModuleDescription const & moduleDescription () const
 
virtual ~EDProducer ()
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector< edm::ProductResolverIndex > const & indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription const &)> registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, std::unordered_multimap< std::string, edm::ProductResolverIndex > const &iIndicies, std::string const &moduleLabel)
 
virtual ~ProducerBase () noexcept(false)
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
virtual ~EDConsumerBase () noexcept(false)
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 description of the run-time parameters More...
 
- Static Public Member Functions inherited from edm::EDProducer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 

Private Member Functions

int analyzeMothers (const reco::Candidate *thisParticle, int &topDaughterQId, int &topBarDaughterQId, std::vector< const reco::Candidate * > &hadMothers, std::vector< std::vector< int > > &hadMothersIndices, std::set< const reco::Candidate * > *analyzedParticles, const int prevPartIndex)
 do a recursive search for the mother particles until the b-quark is found or the absolute mother is found More...
 
bool checkForLoop (std::vector< const reco::Candidate * > &particleChain, const reco::Candidate *particle)
 
std::vector< int > findHadronJets (const reco::GenParticleCollection *genParticles, const reco::JetFlavourInfoMatchingCollection *jetFlavourInfos, std::vector< int > &hadIndex, std::vector< reco::GenParticle > &hadMothersGenPart, std::vector< std::vector< int > > &hadMothersIndices, std::vector< int > &hadLeptonIndex, std::vector< int > &hadLeptonHadIndex, std::vector< int > &hadLeptonViaTau, std::vector< int > &hadFlavour, std::vector< int > &hadFromTopWeakDecay, std::vector< int > &hadBHadronId)
 identify the jets that contain b-hadrons More...
 
int findInMothers (int idx, std::vector< int > &mothChains, const std::vector< std::vector< int > > &hadMothersIndices, const std::vector< reco::GenParticle > &hadMothers, int status, int pdgId, bool pdgAbs, int stopId, int firstLast, bool verbose)
 helper function to find indices of particles with particular pdgId and status from the list of mothers of a given particle More...
 
bool fixExtraSameFlavours (const unsigned int hadId, const std::vector< int > &hadIndices, const std::vector< reco::GenParticle > &hadMothers, const std::vector< std::vector< int > > &hadMothersIndices, const std::vector< int > &isFromTopWeakDecay, const std::vector< std::vector< int > > &LastQuarkIds, const std::vector< std::vector< int > > &LastQuarkMotherIds, std::vector< int > &lastQuarkIndices, std::vector< int > &hadronFlavour, std::set< int > &checkedHadronIds, const int lastQuarkIndex)
 
int flavourSign (const int pdgId)
 Sign of the flavour (matter/antimatter) More...
 
std::string getParticleName (int id) const
 
bool hasHadronDaughter (const int flavour, const reco::Candidate *thisParticle)
 Check if the particle has bHadron among daughters. More...
 
int idInList (std::vector< const reco::Candidate * > particleList, const reco::Candidate *particle)
 Check if the cpecified particle is already in the list of particles. More...
 
int idInList (std::vector< int > list, const int value)
 
bool isBaryonPdgId (const int flavour, const int pdgId)
 Check the pdgId if it represents a baryon of particular flavour. More...
 
bool isHadron (const int flavour, const reco::Candidate *thisParticle)
 Check the pdgId of a given particle if it is a hadron. More...
 
bool isHadronPdgId (const int flavour, const int pdgId)
 Check the pdgId if it represents a hadron of particular flavour. More...
 
bool isMesonPdgId (const int flavour, const int pdgId)
 Check the pdgId if it represents a meson of particular flavour. More...
 
bool isNeutralPdg (int pdgId)
 Check whether a given pdgId represents neutral particle. More...
 
virtual void produce (edm::Event &, const edm::EventSetup &)
 
bool putMotherIndex (std::vector< std::vector< int > > &hadMothersIndices, int partIndex, int mothIndex)
 puts mother index to the list of mothers of particle, if it isn't there already More...
 

Private Attributes

int flavour_
 
std::string flavourStr_
 
edm::EDGetTokenT< reco::GenParticleCollectiongenParticlesToken_
 
edm::EDGetTokenT< reco::JetFlavourInfoMatchingCollectionjetFlavourInfosToken_
 
bool noBBbarResonances_
 
bool onlyJetClusteredHadrons_
 
edm::ESHandle< ParticleDataTablepdt_
 

Additional Inherited Members

- Public Types inherited from edm::EDProducer
typedef EDProducer ModuleType
 
- Public Types inherited from edm::ProducerBase
typedef ProductRegistryHelper::TypeLabelList TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- 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

Finds the origin of each heavy flavour hadron and associated jets to it.

Starting from each consituent of each jet, tracks back in chain to find heavy flavour hadrons. From each hadron traces back until finds the b quark and its mother. For each hadron identifies the jet to which it was injected as a ghost hadron.

The description of the run-time parameters can be found at fillDescriptions()

The description of the products can be found at GenHFHadronMatcher()

Definition at line 58 of file GenHFHadronMatcher.cc.

Constructor & Destructor Documentation

GenHFHadronMatcher::GenHFHadronMatcher ( const edm::ParameterSet cfg)
explicit

constructor initialising producer products and config parameters

Warning
Definition of b-hadron and anti-b-hadron: The term b-hadron and anti-b-hadron is in reference to the quark content and not to distinguish particles from anti-particles. Here a b-hadron contains a b-quark and an anti-b-hadron contains an anti-b-quark. For mesons this means an inversion with respect to the PDG definition, as mesons actually contain anti-b-quarks and anti-mesons contain b-quarks.

Definition at line 132 of file GenHFHadronMatcher.cc.

References funct::abs(), flavour_, flavourStr_, edm::ParameterSet::getParameter(), noBBbarResonances_, and onlyJetClusteredHadrons_.

132  :
133 genParticlesToken_(consumes<reco::GenParticleCollection>(cfg.getParameter<edm::InputTag>("genParticles"))),
134 jetFlavourInfosToken_(consumes<reco::JetFlavourInfoMatchingCollection>(cfg.getParameter<edm::InputTag>("jetFlavourInfos")))
135 {
136  flavour_ = cfg.getParameter<int> ( "flavour" );
137  noBBbarResonances_ = cfg.getParameter<bool> ( "noBBbarResonances" );
138  onlyJetClusteredHadrons_ = cfg.getParameter<bool> ( "onlyJetClusteredHadrons" );
139 
140  flavour_ = std::abs( flavour_ ); // Make flavour independent of sign given in configuration
141  if ( flavour_==5 ) {
142  flavourStr_="B";
143  } else if ( flavour_==4 ) {
144  flavourStr_="C";
145  } else {
146  edm::LogError ( "GenHFHadronMatcher" ) << "Flavour option must be 4 (c-jet) or 5 (b-jet), but is: " << flavour_ << ". Correct this!";
147  }
148 
149  // Hadron matching products
150  produces< std::vector<reco::GenParticle> > ( "gen"+flavourStr_+"HadPlusMothers" ); // All mothers in all decay chains above any hadron of specified flavour
151  produces< std::vector< std::vector<int> > > ( "gen"+flavourStr_+"HadPlusMothersIndices" ); // Indices of mothers of each hadMother
152  produces< std::vector<int> > ( "gen"+flavourStr_+"HadIndex" ); // Index of hadron in the vector of hadMothers
153  produces< std::vector<int> > ( "gen"+flavourStr_+"HadFlavour" ); // PdgId of the first non-b(c) quark mother with sign corresponding to hadron charge
154  produces< std::vector<int> > ( "gen"+flavourStr_+"HadJetIndex" ); // Index of genJet matched to each hadron by jet clustering algorithm
155  produces< std::vector<int> > ( "gen"+flavourStr_+"HadLeptonIndex" ); // Index of lepton found among the hadron decay products in the list of mothers
156  produces< std::vector<int> > ( "gen"+flavourStr_+"HadLeptonHadronIndex" ); // Index of hadron the lepton is associated to
157  produces< std::vector<int> > ( "gen"+flavourStr_+"HadLeptonViaTau" ); // Whether lepton comes directly from hadron or via tau decay
158  produces< std::vector<int> > ( "gen"+flavourStr_+"HadFromTopWeakDecay" ); // Tells whether the hadron appears in the chain after top decay
159  produces< std::vector<int> > ( "gen"+flavourStr_+"HadBHadronId" ); // Index of a b-hadron which the current hadron comes from (for c-hadrons)
160 }
T getParameter(std::string const &) const
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::EDGetTokenT< reco::JetFlavourInfoMatchingCollection > jetFlavourInfosToken_
edm::EDGetTokenT< reco::GenParticleCollection > genParticlesToken_
GenHFHadronMatcher::~GenHFHadronMatcher ( )

Definition at line 162 of file GenHFHadronMatcher.cc.

163 {
164 }

Member Function Documentation

int GenHFHadronMatcher::analyzeMothers ( const reco::Candidate thisParticle,
int &  topDaughterQId,
int &  topBarDaughterQId,
std::vector< const reco::Candidate * > &  hadMothers,
std::vector< std::vector< int > > &  hadMothersIndices,
std::set< const reco::Candidate * > *  analyzedParticles,
const int  prevPartIndex 
)
private

do a recursive search for the mother particles until the b-quark is found or the absolute mother is found

the treatment of b-bar resonances depends on the global parameter noBBbarResonances_

Parameters
[in]thisParticlecurrent particle from which starts the search of the hadron and all its mothers up to proton
[out]hadronthe last hadron in the decay chain [that decays weekly]
[out]leptonlepton found in the current decay chain
[out]topDaughterQIdId of the top quark daughter b(c) quark
[out]topBarDaughterQIdId of the antitop quark daughter b(c) quark
[out]hadMotherslist of all processed particles ending with proton
[out]hadMothersIndiceslist of i-vectors containing j-indices representing particles that are mothers of each i-particle from hadMothers
[out]analyzedParticleslist of particles analysed in the chain [used for infinite loop detection]
[out]prevPartIndexindex of the previous particle in the current chain [used for infinite loop detection]
Returns
index of hadron in the hadMothers list [-1 if no hadron found]

Definition at line 605 of file GenHFHadronMatcher.cc.

References funct::abs(), mps_fire::i, idInList(), diffTreeTool::index, reco::Candidate::mother(), reco::Candidate::numberOfMothers(), reco::Candidate::pdgId(), reco::Candidate::pt(), and putMotherIndex().

Referenced by findHadronJets().

606 {
607  // Getting the index of the particle which is a hadron in the first call
608  int hadronIndex=-1; // Index of the hadron that is returned by this function
609  int index = idInList( hadMothers, thisParticle );
610  if ( index<0 ) { // If hadron is not in the list of mothers yet
611  hadMothers.push_back ( thisParticle );
612  hadronIndex=hadMothers.size()-1;
613  } else { // If hadron is in the list of mothers already
614  hadronIndex=index;
615  }
616 
617  int partIndex = -1; // Index of particle being checked in the list of mothers
618  partIndex = idInList( hadMothers, thisParticle );
619 
620  // Checking whether this particle is already in the chain of analyzed particles in order to identify a loop
621  bool isLoop = false;
622  if ( !analyzedParticles ) {
623  analyzedParticles = new std::set<const reco::Candidate*>;
624  }
625  for ( unsigned int i=0; i<analyzedParticles->size(); i++ ) {
626  if ( analyzedParticles->count ( thisParticle ) <=0 ) {
627  continue;
628  }
629  isLoop = true;
630  break;
631  }
632 
633  // If a loop has been detected
634  if ( isLoop ) {
635  if ( prevPartIndex>=0 ) {
636  putMotherIndex ( hadMothersIndices, prevPartIndex, -1 ); // Setting mother index of previous particle to -1
637  }
638  return hadronIndex; // Stopping further processing of the current chain
639  }
640  analyzedParticles->insert ( thisParticle );
641 
642  // Putting the mothers to the list of mothers
643  for ( size_t iMother = 0; iMother < thisParticle->numberOfMothers(); ++iMother ) {
644  const reco::Candidate* mother = thisParticle->mother ( iMother );
645  int mothIndex = idInList( hadMothers, mother );
646  if ( mothIndex == partIndex && partIndex>=0 ) {
647  continue; // Skipping the mother that is its own daughter
648  }
649 
650  // If this mother isn't yet in the list and hadron or lepton is in the list
651  if ( mothIndex<0 ) {
652  hadMothers.push_back ( mother );
653  mothIndex=hadMothers.size()-1;
654  }
655  // If hadron has already been found in current chain and the mother isn't a duplicate of the particle being checked
656  if ( mothIndex!=partIndex && partIndex>=0 ) {
657  putMotherIndex ( hadMothersIndices, partIndex, mothIndex ); // Putting the index of mother for current particle
658  }
659  analyzeMothers ( mother, topDaughterQId, topBarDaughterQId, hadMothers, hadMothersIndices, analyzedParticles, partIndex );
660  // Setting the id of the particle which is a quark from the top decay
661  if(std::abs(mother->pdgId())==6) {
662  int& bId = mother->pdgId() < 0 ? topBarDaughterQId : topDaughterQId;
663  int thisFlav = std::abs(thisParticle->pdgId());
664  if( bId<0){
665  if(thisFlav <= 5) bId = partIndex;
666  } else {
667  int bIdFlav = std::abs(hadMothers.at(bId)->pdgId());
668  if( bIdFlav != 5 && thisFlav == 5) bId = partIndex;
669  else if( thisFlav == 5 && thisParticle->pt() > hadMothers.at(bId)->pt() ) bId = partIndex;
670  } // If daughter quark of the top not found yet
671  } // If the mother is a top quark and hadron has been found
672  } // End of loop over mothers
673 
674  analyzedParticles->erase ( thisParticle ); // Removing current particle from the current chain that is being analyzed
675 
676  if ( partIndex<0 ) {
677  return hadronIndex; // Safety check
678  }
679 
680  // Adding -1 to the list of mother indices for current particle if it has no mothers (for consistency between numbering of indices and mothers)
681  if ( (int)thisParticle->numberOfMothers() <=0) {
682  putMotherIndex ( hadMothersIndices, partIndex, -1 );
683  }
684 
685  return hadronIndex;
686 
687 }
int idInList(std::vector< const reco::Candidate * > particleList, const reco::Candidate *particle)
Check if the cpecified particle is already in the list of particles.
virtual const Candidate * mother(size_type i=0) const =0
return pointer to mother
int analyzeMothers(const reco::Candidate *thisParticle, int &topDaughterQId, int &topBarDaughterQId, std::vector< const reco::Candidate * > &hadMothers, std::vector< std::vector< int > > &hadMothersIndices, std::set< const reco::Candidate * > *analyzedParticles, const int prevPartIndex)
do a recursive search for the mother particles until the b-quark is found or the absolute mother is f...
virtual int pdgId() const =0
PDG identifier.
virtual size_type numberOfMothers() const =0
number of mothers (zero or one in most of but not all the cases)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool putMotherIndex(std::vector< std::vector< int > > &hadMothersIndices, int partIndex, int mothIndex)
puts mother index to the list of mothers of particle, if it isn&#39;t there already
virtual double pt() const =0
transverse momentum
bool GenHFHadronMatcher::checkForLoop ( std::vector< const reco::Candidate * > &  particleChain,
const reco::Candidate particle 
)
private
void GenHFHadronMatcher::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

description of the run-time parameters

Definition at line 172 of file GenHFHadronMatcher.cc.

References edm::ConfigurationDescriptions::add(), edm::ParameterSetDescription::add(), funct::false, and funct::true.

173 {
175  desc.add<edm::InputTag>("genParticles")->setComment( "Collection of GenParticle objects which contains all particles produced in the event" );
176  desc.add<edm::InputTag>("jetFlavourInfos")->setComment( "Output from the JetFlavour tool. Contains information about partons/hadrons/leptons associated to jets" );
177  desc.add<bool> ( "noBBbarResonances",true )->setComment ( "Whether resonances should not be treated as hadrons" );
178  desc.add<bool> ( "onlyJetClusteredHadrons",false )->setComment ( "Whether only hadrons that are matched to jets should be analysed. Runs x1000 faster in Sherpa" );
179  desc.add<int> ( "flavour",5 )->setComment ( "Flavour of weakly decaying hadron that should be analysed (4-c, 5-b)" );
180  descriptions.add ( "matchGenHFHadron",desc );
181 }
ParameterDescriptionBase * add(U const &iLabel, T const &value)
void add(std::string const &label, ParameterSetDescription const &psetDescription)
std::vector< int > GenHFHadronMatcher::findHadronJets ( const reco::GenParticleCollection genParticles,
const reco::JetFlavourInfoMatchingCollection jetFlavourInfos,
std::vector< int > &  hadIndex,
std::vector< reco::GenParticle > &  hadMothers,
std::vector< std::vector< int > > &  hadMothersIndices,
std::vector< int > &  hadLeptonIndex,
std::vector< int > &  hadLeptonHadIndex,
std::vector< int > &  hadLeptonViaTau,
std::vector< int > &  hadFlavour,
std::vector< int > &  hadFromTopWeakDecay,
std::vector< int > &  hadBHadronId 
)
private

identify the jets that contain b-hadrons

For each hadron all mothers from all levels and chains are analysed to find the quark or gluon from which the hadron has originated. This is done by searching through the generator particle decay tree starting from the hadron, performing checks for flavour and kinematic consistency. The hadrons that are not last in the decay chain (don't decay weakly) are skipped. Additionally for each hadron it is checked whether it comes from the top weak decay branch and whether it comes from some other b-hadron decay

b-bbar (c-cbar) resonances can either be considered as hadrons or not, depending on the configuration.

Parameters
[out]hadIndexvector of indices of found hadrons in hadMothers
[out]hadMothersvector of all mothers at all levels of each found hadron
[out]hadMothersIndicesconnection between each particle from hadMothers and its mothers
[out]hadLeptonIndexindex of lepton among the hadMothers
[out]hadLeptonHadIndexindex of hadron associated to each lepton
[out]hadLeptonViaTauwhether lepton is from direct hadron->lepton or hadron->tau->lepton
[out]hadFlavourflavour of each found hadron
[out]hadFromTopWeakDecaywhether each hadron contains the top quark in its decay chain [works only for B-Hadrons]
[out]hadBHadronIdfor each hadron - index of the ancestor b-hadron [-1 if hadron doesn't come from another b-hdaron]
Returns
vector of jet indices that were matched to each hadron [by the jet clustering algorithm]

Definition at line 252 of file GenHFHadronMatcher.cc.

References funct::abs(), analyzeMothers(), edm::AssociationVector< KeyRefProd, CVal, KeyRef, SizeType, KeyReferenceHelper >::begin(), edm::RefVector< C, T, F >::begin(), deltaR(), PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi::dR, edm::AssociationVector< KeyRefProd, CVal, KeyRef, SizeType, KeyReferenceHelper >::end(), edm::RefVector< C, T, F >::end(), funct::false, spr::find(), findInMothers(), fixExtraSameFlavours(), flavour_, flavourSign(), reco::JetFlavourInfo::getbHadrons(), reco::JetFlavourInfo::getcHadrons(), hasHadronDaughter(), mps_fire::i, createfilelist::int, isHadron(), reco::CompositeRefCandidateT< D >::mother(), reco::Candidate::mother(), onlyJetClusteredHadrons_, cosmictrackingParticleSelector_cfi::pdgId, reco::Candidate::pdgId(), and reco::LeafCandidate::pdgId().

Referenced by produce().

259 {
260  std::vector<int> hadJetIndex;
261  std::vector<const reco::Candidate*> hadMothersCand;
262 
263  int topDaughterQId = -1;
264  int topBarDaughterQId = -1;
265 
266  // Looping over all jets to get hadrons associated to them
267  for(reco::JetFlavourInfoMatchingCollection::const_iterator i_info = jetFlavourInfos->begin(); i_info != jetFlavourInfos->end(); ++i_info){
268  reco::JetFlavourInfo jetInfo = i_info->second;
269  const int jetIndex = i_info - jetFlavourInfos->begin();
270  // Looping over each hadron associated with the jet and finding its origin
271  const reco::GenParticleRefVector& hadronsInJet = flavour_==5 ? jetInfo.getbHadrons() : flavour_==4 ? jetInfo.getcHadrons() : reco::GenParticleRefVector();
272  for(reco::GenParticleRefVector::const_iterator hadron = hadronsInJet.begin(); hadron != hadronsInJet.end(); ++hadron) {
273  // Check that the hadron satisfies criteria configured in the module
274  if(!isHadron ( flavour_, (&**hadron) )) continue;
275  if(hasHadronDaughter ( flavour_, (reco::Candidate*)(&**hadron) )) continue;
276  // Scanning the chain starting from the hadron
277  int hadronIndex = analyzeMothers ( (reco::Candidate*)(&**hadron), topDaughterQId, topBarDaughterQId, hadMothersCand, hadMothersIndices, 0, -1 );
278  // Storing the index of the hadron to the list
279  hadIndex.push_back ( hadronIndex );
280  hadJetIndex.push_back ( jetIndex ); // Putting jet index to the result list
281  }
282  } // End of loop over jets
283 
284  // Access all hadrons which are not associated with jets, if requested
286  for(reco::GenParticleCollection::const_iterator i_particle = genParticles->begin(); i_particle != genParticles->end(); ++i_particle){
287  const reco::GenParticle* thisParticle = &*i_particle;
288  if(!isHadron(flavour_, thisParticle)) continue;
289  // Skipping the hadron if it was already found directly from jets
290  if(std::find(hadMothersCand.begin(), hadMothersCand.end(), thisParticle) != hadMothersCand.end()) continue;
291 
292  // Scanning the chain starting from the hadron
293  int hadronIndex = analyzeMothers ( thisParticle, topDaughterQId, topBarDaughterQId, hadMothersCand, hadMothersIndices, 0, -1 );
294  // Storing the index of the hadron to the list
295  hadIndex.push_back ( hadronIndex );
296  hadJetIndex.push_back ( -1 ); // Jet index undefined
297  }
298  }
299 
300  // Transfering Candidates to the list of processed particles for further analysis
301  for ( int i=0; i< (int)hadMothersCand.size(); i++ ) {
302  const reco::GenParticle* particle = dynamic_cast<const reco::GenParticle*>( hadMothersCand.at(i) );
303  hadMothers.push_back(*particle);
304  }
305 
306  // Adding leptons from hadron decays
307  for(reco::GenParticleCollection::const_iterator i_particle = genParticles->begin(); i_particle != genParticles->end(); ++i_particle){
308  const reco::GenParticle lepton = *i_particle;
309  const int pdg_abs = lepton.pdgId();
310  // Skipping if not a lepton: e/mu
311  if(pdg_abs != 11 && pdg_abs != 13) continue;
312  bool leptonViaTau = false;
313  const reco::Candidate* leptonMother = lepton.mother();
314  if(!leptonMother) continue;
315  // Taking next mother if direct mother is a tau
316  if(std::abs(leptonMother->pdgId()) == 15) {
317  leptonViaTau = 1;
318  leptonMother = leptonMother->mother();
319  }
320  // Skipping this lepton if its mother is not a proper hadron
321  if(!isHadron(flavour_, leptonMother)) continue;
322  // Finding the index of this hadron in the list of analysed particles
323  size_t leptonHadronParticleIndex = std::find(hadMothersCand.begin(), hadMothersCand.end(), leptonMother) - hadMothersCand.begin();
324  if(leptonHadronParticleIndex >= hadMothersCand.size()) continue;
325  // Finding the actual hadron index among those that were found
326  size_t leptonHadronIndex = std::find(hadIndex.begin(), hadIndex.end(), leptonHadronParticleIndex) - hadIndex.begin();
327  if(leptonHadronIndex >= hadIndex.size()) continue;
328  // Putting the lepton, its index and hadron index to the corresponding lists
329  hadMothers.push_back(lepton);
330  const int leptonIndex = hadMothersCand.size()-1;
331  hadLeptonIndex.push_back(leptonIndex);
332  hadLeptonViaTau.push_back(leptonViaTau);
333  hadLeptonHadIndex.push_back(leptonHadronIndex);
334  }
335 
336  // Checking mothers of hadrons in order to assign flags (where the hadron comes from)
337  unsigned int nHad = hadIndex.size();
338 
339  std::vector<std::vector<int> > LastQuarkMotherIds;
340  std::vector<std::vector<int> > LastQuarkIds;
341  std::vector<int> lastQuarkIndices(nHad, -1);
342 
343  // Looping over all hadrons
344  for ( unsigned int hadNum=0; hadNum<nHad; hadNum++ ) {
345  unsigned int hadIdx = hadIndex.at(hadNum); // Index of hadron in the hadMothers
346 
347  std::vector <int> FirstQuarkId;
348  std::vector <int> LastQuarkId;
349  std::vector <int> LastQuarkMotherId;
350 
351  const int hadronFlavourSign = flavourSign( hadMothers.at(hadIdx).pdgId() );
352 
353  // Searching only first quark in the chain with the same flavour as hadron
354  if(hadronFlavourSign != 0) {
355  findInMothers( hadIdx, FirstQuarkId, hadMothersIndices, hadMothers, 0, hadronFlavourSign*flavour_, false, -1, 1, false );
356  }
357  // Searching for quarks with both flavours since it is a bb/cc resonance
358  else {
359  findInMothers( hadIdx, FirstQuarkId, hadMothersIndices, hadMothers, 0, flavour_, false, -1, 1, false );
360  findInMothers( hadIdx, FirstQuarkId, hadMothersIndices, hadMothers, 0, -1*flavour_, false, -1, 1, false );
361  }
362 
363  // Finding last quark for each first quark
364  for ( unsigned int qId=0; qId<FirstQuarkId.size(); qId++ ) {
365  // Identifying the flavour of the first quark to find the last quark of the same flavour
366  const int quarkFlavourSign = flavourSign( hadMothers.at(FirstQuarkId.at(qId)).pdgId() );
367  // Finding last quark of the hadron starting from the first quark
368  findInMothers( FirstQuarkId.at(qId), LastQuarkId, hadMothersIndices, hadMothers, 0, quarkFlavourSign*flavour_, false, -1, 2, false );
369  } // End of loop over all first quarks of the hadron
370 
371 
372  // Setting initial flavour of the hadron
373  int hadronFlavour = 0;
374 
375  // Initialising pairs of last quark index and its distance from the hadron (to sort by distance)
376  std::vector<std::pair<double, int> > lastQuark_dR_id_pairs;
377 
378  // Finding the closest quark in dR
379  for ( unsigned int qId = 0; qId < LastQuarkId.size(); qId++ ) {
380  int qIdx = LastQuarkId.at(qId);
381  // Calculating the dR between hadron and quark
382  float dR = deltaR ( hadMothers.at(hadIdx).eta(),hadMothers.at(hadIdx).phi(),hadMothers.at(qIdx).eta(),hadMothers.at(qIdx).phi() );
383 
384  std::pair<double, int> dR_hadId_pair(dR,qIdx);
385  lastQuark_dR_id_pairs.push_back(dR_hadId_pair);
386  } // End of loop over all last quarks of the hadron
387 
388  std::sort(lastQuark_dR_id_pairs.begin(), lastQuark_dR_id_pairs.end());
389 
390  if(lastQuark_dR_id_pairs.size() > 1) {
391  double dRratio = (lastQuark_dR_id_pairs.at(1).first - lastQuark_dR_id_pairs.at(0).first)/lastQuark_dR_id_pairs.at(1).first;
392  int qIdx_closest = lastQuark_dR_id_pairs.at(0).second;
393  LastQuarkId.clear();
394  if(dRratio>0.5) LastQuarkId.push_back(qIdx_closest);
395  else for(std::pair<double, int> qIdDrPair : lastQuark_dR_id_pairs) LastQuarkId.push_back(qIdDrPair.second);
396  }
397  for(int qIdx : LastQuarkId) {
398  int qmIdx = hadMothersIndices.at(qIdx).at(0);
399  LastQuarkMotherId.push_back(qmIdx);
400  }
401 
402  if((int)LastQuarkId.size()>0) lastQuarkIndices.at(hadNum) = 0; // Setting the first quark in array as a candidate if it exists
403 
404  LastQuarkIds.push_back( LastQuarkId );
405 
406  LastQuarkMotherIds.push_back ( LastQuarkMotherId );
407 
408  if(LastQuarkMotherId.size()<1) {
409  hadronFlavour = 0;
410  } else {
411  int qIdx = LastQuarkId.at( lastQuarkIndices.at(hadNum) );
412  int qFlav = flavourSign( hadMothers.at(qIdx).pdgId() );
413  hadronFlavour = qFlav*std::abs( hadMothers.at( LastQuarkMotherId.at( lastQuarkIndices.at(hadNum) ) ).pdgId() );
414  }
415  hadFlavour.push_back(hadronFlavour); // Adding hadron flavour to the list of flavours
416 
417  // Checking whether hadron comes from the Top weak decay
418  int isFromTopWeakDecay = 1;
419  std::vector <int> checkedParticles;
420  if(hadFlavour.at(hadNum) != 0) {
421  int lastQIndex = LastQuarkId.at(lastQuarkIndices.at(hadNum));
422  bool fromTB = topDaughterQId >= 0 ? findInMothers( lastQIndex, checkedParticles, hadMothersIndices, hadMothers, -1, 0, false, topDaughterQId, 2, false ) >= 0 : false;
423  checkedParticles.clear();
424  bool fromTbarB = topBarDaughterQId >= 0 ? findInMothers( lastQIndex, checkedParticles, hadMothersIndices, hadMothers, -1, 0, false, topBarDaughterQId, 2, false) >= 0 : false;
425  checkedParticles.clear();
426  if(!fromTB && !fromTbarB) {
427  isFromTopWeakDecay = 0;
428  }
429  } else isFromTopWeakDecay = 2;
430  hadFromTopWeakDecay.push_back(isFromTopWeakDecay);
431  int bHadronMotherId = findInMothers( hadIdx, checkedParticles, hadMothersIndices, hadMothers, 0, 555555, true, -1, 1, false );
432  hadBHadronId.push_back(bHadronMotherId);
433 
434 
435  if(LastQuarkMotherId.size()>0) {
436  std::set<int> checkedHadronIds;
437  fixExtraSameFlavours(hadNum, hadIndex, hadMothers, hadMothersIndices, hadFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadFlavour, checkedHadronIds, 0);
438  }
439 
440  } // End of loop over all hadrons
441 
442  return hadJetIndex;
443 }
const GenParticleRefVector & getbHadrons() const
Return a vector of GenParticleRef&#39;s to b hadrons clustered inside the jet.
int flavourSign(const int pdgId)
Sign of the flavour (matter/antimatter)
const GenParticleRefVector & getcHadrons() const
Return a vector of GenParticleRef&#39;s to c hadrons clustered inside the jet.
const_iterator end() const
virtual const Candidate * mother(size_type i=0) const =0
return pointer to mother
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:253
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:248
int findInMothers(int idx, std::vector< int > &mothChains, const std::vector< std::vector< int > > &hadMothersIndices, const std::vector< reco::GenParticle > &hadMothers, int status, int pdgId, bool pdgAbs, int stopId, int firstLast, bool verbose)
helper function to find indices of particles with particular pdgId and status from the list of mother...
int analyzeMothers(const reco::Candidate *thisParticle, int &topDaughterQId, int &topBarDaughterQId, std::vector< const reco::Candidate * > &hadMothers, std::vector< std::vector< int > > &hadMothersIndices, std::set< const reco::Candidate * > *analyzedParticles, const int prevPartIndex)
do a recursive search for the mother particles until the b-quark is found or the absolute mother is f...
virtual int pdgId() const =0
PDG identifier.
virtual int pdgId() const final
PDG identifier.
Class storing the jet flavour information.
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
edm::RefVector< GenParticleCollection > GenParticleRefVector
vector of reference to GenParticle in the same collection
bool isHadron(const int flavour, const reco::Candidate *thisParticle)
Check the pdgId of a given particle if it is a hadron.
const_iterator begin() const
virtual const Candidate * mother(size_type=0) const
return mother at a given position, i = 0, ... numberOfMothers() - 1 (read only mode) ...
bool hasHadronDaughter(const int flavour, const reco::Candidate *thisParticle)
Check if the particle has bHadron among daughters.
bool fixExtraSameFlavours(const unsigned int hadId, const std::vector< int > &hadIndices, const std::vector< reco::GenParticle > &hadMothers, const std::vector< std::vector< int > > &hadMothersIndices, const std::vector< int > &isFromTopWeakDecay, const std::vector< std::vector< int > > &LastQuarkIds, const std::vector< std::vector< int > > &LastQuarkMotherIds, std::vector< int > &lastQuarkIndices, std::vector< int > &hadronFlavour, std::set< int > &checkedHadronIds, const int lastQuarkIndex)
int GenHFHadronMatcher::findInMothers ( int  idx,
std::vector< int > &  mothChains,
const std::vector< std::vector< int > > &  hadMothersIndices,
const std::vector< reco::GenParticle > &  hadMothers,
int  status,
int  pdgId,
bool  pdgAbs = false,
int  stopId = -1,
int  firstLast = 0,
bool  verbose = false 
)
private

helper function to find indices of particles with particular pdgId and status from the list of mothers of a given particle

Parameters
[in]idxindex of particle, mothers of which should be searched
[in]mothChainsvector of indices where the found mothers are stored
[in]hadMothersIndiceslist of indices pointing to mothers of each particle from list of mothers
[in]hadMothersvector of all hadron mother particles of all levels
[in]statusstatus of mother that is being looked for
[in]pdgIdpdgId of mother that is being looked for [flavour*111111 used to identify hadrons of particular flavour]
[in]pdgAbswhether the sign of pdgId should be taken into account
[in]stopIdid of the particle in the hadMothers array after which the checking should stop
[in]firstLastshould all(0), first(1) or last(2) occurances of the searched particle be stored
[in]verboseoption to print every particle that is processed during the search
Returns
index of the found particle in the hadMothers array [-1 if the specified particle not found]

Definition at line 752 of file GenHFHadronMatcher.cc.

References funct::abs(), funct::false, spr::find(), mps_fire::i, training_settings::idx, isHadronPdgId(), and funct::true.

Referenced by findHadronJets().

755 {
756  int foundStopId = -1;
757  int pdg_1 = hadMothers.at(idx).pdgId();
758 
759  if ( (int)hadMothersIndices.size() <= idx ) {
760  if ( verbose ) {
761  printf ( " Stopping checking particle %d. No mothers are stored.\n",idx );
762  }
763  return -1; // Skipping if no mothers are stored for this particle
764  }
765 
766  if(std::abs(hadMothers.at( idx ).pdgId()) > 10 && std::abs(hadMothers.at( idx ).pdgId()) < 19) printf("Lepton: %d\n", hadMothers.at( idx ).pdgId());
767 
768  std::vector<int> mothers = hadMothersIndices.at(idx);
769  unsigned int nMothers = mothers.size();
770  bool isCorrect=false; // Whether current particle is what is being searched
771  if ( verbose ) {
772  if ( std::abs( hadMothers.at(idx).pdgId() ) ==2212 ) {
773  printf ( "Chk: %d\tpdg: %d\tstatus: %d",idx, hadMothers.at(idx).pdgId(), hadMothers.at(idx).status() );
774  } else {
775  printf ( " Chk: %d(%d mothers)\tpdg: %d\tstatus: %d\tPt: %.3f\tEta: %.3f",idx, nMothers, hadMothers.at(idx).pdgId(), hadMothers.at(idx).status(), hadMothers.at(idx).pt(),hadMothers.at(idx).eta() );
776  }
777  }
778  bool hasCorrectMothers = true;
779  if(nMothers<1) hasCorrectMothers=false; else if(mothers.at(0)<0) hasCorrectMothers=false;
780  if(!hasCorrectMothers) {
781  if(verbose) printf(" NO CORRECT MOTHER\n");
782  return -1;
783  }
784  // Stopping if the particular particle has been found
785  if(stopId>=0 && idx == stopId) return idx;
786 
787  // Stopping if the hadron of particular flavour has been found
788  if(pdgId%111111==0 && pdgId!=0) {
789  if(isHadronPdgId(pdgId/111111, hadMothers.at(idx).pdgId())) {
790  return idx;
791  }
792  }
793 
794  // Checking whether current mother satisfies selection criteria
795  if ( ( ( hadMothers.at(idx).pdgId() == pdgId && pdgAbs==false )
796  || ( std::abs( hadMothers.at(idx).pdgId() ) == std::abs( pdgId ) && pdgAbs==true ) )
797  && ( hadMothers.at(idx).status() == status || status==0 )
798  && hasCorrectMothers ) {
799  isCorrect=true;
800  // Adding to the list of candidates if not there and if mother of this quark has correct flavour sign
801  const bool inList = std::find(mothChains.begin(), mothChains.end(), idx) != mothChains.end();
802  if ( !inList && mothers.at(0) >= 0 && ( hadMothers.at(idx).pdgId()*pdgId > 0 || !pdgAbs ) ) {
803  if ( firstLast==0 || firstLast==1 ) {
804  mothChains.push_back(idx);
805  }
806  if ( verbose ) {
807  printf ( " *" );
808  }
809  }
810  if ( verbose ) {
811  printf ( " +++" );
812  }
813  }
814  if ( verbose ) {
815  printf ( "\n" );
816  }
817 
818  if ( isCorrect && firstLast==1 ) {
819  return -1; // Stopping if only the first particle in the chain is looked for
820  }
821 
822  // Checking next level mothers
823  unsigned int nDifferingMothers = 0;
824  for ( unsigned int i = 0; i < nMothers; i++ ) {
825  int idx2 = mothers.at(i);
826  if ( idx2 < 0 ) {
827  if(verbose) printf("^^^ Has no mother\n");
828  continue; // Skipping if mother's id is -1 (no mother), that means current particle is a proton
829  }
830  if ( idx2 == idx ) {
831  if(verbose) printf("^^^ Stored as its own mother\n");
832  continue; // Skipping if particle is stored as its own mother
833  }
834  int pdg_2 = hadMothers[idx2].pdgId();
835  // Inverting the flavour if bb oscillation detected
836  if ( isHadronPdgId(pdgId, pdg_1) && isHadronPdgId(pdgId, pdg_2) && pdg_1*pdg_2 < 0 ) {
837  pdgId *= -1;
838  if(verbose) printf("######### Inverting flavour of the hadron\n");
839  }
840  // Counting how many mothers are different from this particle
841  if ( ( std::abs( pdg_2 ) != std::abs( pdgId ) && pdgAbs==true ) ||
842  ( pdg_2 != pdgId && pdgAbs == false ) ) {
843  nDifferingMothers++;
844  }
845 
846  // Checking next level mother
847  if ( verbose ) {
848  printf ( "Checking mother %d out of %d mothers (%d -> %d), looking for pdgId: %d\n",i,nMothers,idx, idx2, pdgId );
849  }
850  if(firstLast==2 && pdg_1 != pdg_2) continue; // Requiring the same flavour when finding the last quark
851  foundStopId = findInMothers ( idx2, mothChains, hadMothersIndices, hadMothers, status, pdgId, pdgAbs, stopId, firstLast, verbose );
852  }
853  // Storing this particle if all its mothers are of another type and the last of its kind should be stored
854  if(firstLast==2 && isCorrect && nDifferingMothers >= nMothers) {
855  if ( verbose ) {
856  printf ( "Checking particle %d once more to store it as the last quark\n",idx);
857  }
858  foundStopId = findInMothers ( idx, mothChains, hadMothersIndices, hadMothers, 0, pdgId, pdgAbs, stopId, 1, verbose );
859  }
860 
861  return foundStopId;
862 }
bool isHadronPdgId(const int flavour, const int pdgId)
Check the pdgId if it represents a hadron of particular flavour.
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
int findInMothers(int idx, std::vector< int > &mothChains, const std::vector< std::vector< int > > &hadMothersIndices, const std::vector< reco::GenParticle > &hadMothers, int status, int pdgId, bool pdgAbs, int stopId, int firstLast, bool verbose)
helper function to find indices of particles with particular pdgId and status from the list of mother...
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool GenHFHadronMatcher::fixExtraSameFlavours ( const unsigned int  hadId,
const std::vector< int > &  hadIndices,
const std::vector< reco::GenParticle > &  hadMothers,
const std::vector< std::vector< int > > &  hadMothersIndices,
const std::vector< int > &  isFromTopWeakDecay,
const std::vector< std::vector< int > > &  LastQuarkIds,
const std::vector< std::vector< int > > &  LastQuarkMotherIds,
std::vector< int > &  lastQuarkIndices,
std::vector< int > &  hadronFlavour,
std::set< int > &  checkedHadronIds,
const int  lastQuarkIndex 
)
private

Finds hadrons that have the same flavour and origin and resolve this ambiguity fixExtraSameFlavours

Parameters
hadIdIndex of the hadron being checked
hadIndexVector of indices of each hadron pointing to the hadMothers
hadMothersVector of gen particles (contain hadrons and all its ancestors)
hadMothersIndicesVector of indices for each particle from hadMothers
isFromTopWeakDecayVector of values showing whether particle comes from the top weak decay
LastQuarkIdsVector of indices of last quarks for each hadron
LastQuarkMotherIdsVector of indices of mothers for each last quark from LastQuarkIds
lastQuakIndicesVector of indices pointing to particular index from the LastQuarkIds and LastQuarkMotherIds to be used for each hadron
lastQuarkIndexIndex from the LastQuarkIds and LastQuarkMotherIds for this particular hadron with index hadId
Returns
Whether other mother with unique association has been found for the hadrons

Definition at line 897 of file GenHFHadronMatcher.cc.

References funct::abs(), DEFINE_FWK_MODULE, and isNeutralPdg().

Referenced by findHadronJets().

903 {
904  if(checkedHadronIds.count(hadId) != 0) return false; // Hadron already checked previously and should be skipped
905  checkedHadronIds.insert(hadId); // Putting hadron to the list of checked ones in this run
906 
907  if(lastQuarkIndex<0) return false;
908  if((int)LastQuarkIds.at(hadId).size()<lastQuarkIndex+1) return false;
909  int LastQuarkId = LastQuarkIds.at(hadId).at(lastQuarkIndex);
910  int LastQuarkMotherId = LastQuarkMotherIds.at( hadId ).at( lastQuarkIndex );
911  int qmFlav = hadMothers.at(LastQuarkId).pdgId() < 0 ? -1 : 1;
912  int hadFlavour = qmFlav*std::abs( hadMothers.at( LastQuarkMotherId ).pdgId() );
913  bool ambiguityResolved = true;
914  // If last quark has inconsistent flavour with its mother, setting the hadron flavour to gluon
915  if( (hadMothers.at(LastQuarkId).pdgId()*hadMothers.at(LastQuarkMotherId).pdgId() < 0 && !isNeutralPdg(hadMothers.at(LastQuarkMotherId).pdgId())) ||
916  // If particle not coming from the Top weak decay has Top flavour
917  (std::abs(hadronFlavour.at(hadId))==6 && isFromTopWeakDecay.at(hadId)==0) ) {
918  if((int)LastQuarkIds.at(hadId).size()>lastQuarkIndex+1) fixExtraSameFlavours(hadId, hadIndices, hadMothers, hadMothersIndices, isFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadronFlavour, checkedHadronIds, lastQuarkIndex+1);
919  else hadronFlavour.at(hadId) = qmFlav*21;
920  return true;
921  }
922 
923  int nSameFlavourHadrons = 0;
924  // Looping over all previous hadrons
925  for(unsigned int iHad = 0; iHad<hadronFlavour.size(); iHad++) {
926  if(iHad==hadId) continue;
927  int theLastQuarkIndex = lastQuarkIndices.at(iHad);
928  if(theLastQuarkIndex<0) continue;
929  if((int)LastQuarkMotherIds.at( iHad ).size() <= theLastQuarkIndex) continue;
930  int theLastQuarkMotherId = LastQuarkMotherIds.at( iHad ).at( theLastQuarkIndex );
931  int theHadFlavour = hadronFlavour.at(iHad);
932  // Skipping hadrons with different flavour
933  if(theHadFlavour==0 || std::abs(theHadFlavour)==21) continue;
934  if(theHadFlavour != hadFlavour || theLastQuarkMotherId != LastQuarkMotherId) continue;
935  ambiguityResolved = false;
936  nSameFlavourHadrons++;
937 
938  // Checking other b-quark mother candidates of this hadron
939  if((int)LastQuarkIds.at(hadId).size() > lastQuarkIndex+1) {
940  if(fixExtraSameFlavours(hadId, hadIndices, hadMothers, hadMothersIndices, isFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadronFlavour, checkedHadronIds, lastQuarkIndex+1) ) {
941  ambiguityResolved = true;
942  break;
943  }
944  } else
945  // Checking other b-quark mother candidates of the particular previous hadron
946  if((int)LastQuarkIds.at(iHad).size() > theLastQuarkIndex+1) {
947  if(fixExtraSameFlavours(iHad, hadIndices, hadMothers, hadMothersIndices, isFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadronFlavour, checkedHadronIds, theLastQuarkIndex+1) ) {
948  ambiguityResolved = true;
949  break;
950  };
951  }
952 
953  } // End of loop over all previous hadrons
954 
955  checkedHadronIds.erase(hadId); // Removing the hadron from the list of checked hadrons
956  if(nSameFlavourHadrons>0 && !ambiguityResolved) {
957  hadronFlavour.at(hadId) = qmFlav*21;
958  return true;
959  }
960  lastQuarkIndices.at(hadId) = lastQuarkIndex;
961  hadronFlavour.at(hadId) = hadFlavour;
962  return true;
963 }
bool isNeutralPdg(int pdgId)
Check whether a given pdgId represents neutral particle.
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool fixExtraSameFlavours(const unsigned int hadId, const std::vector< int > &hadIndices, const std::vector< reco::GenParticle > &hadMothers, const std::vector< std::vector< int > > &hadMothersIndices, const std::vector< int > &isFromTopWeakDecay, const std::vector< std::vector< int > > &LastQuarkIds, const std::vector< std::vector< int > > &LastQuarkMotherIds, std::vector< int > &lastQuarkIndices, std::vector< int > &hadronFlavour, std::set< int > &checkedHadronIds, const int lastQuarkIndex)
int GenHFHadronMatcher::flavourSign ( const int  pdgId)
private

Sign of the flavour (matter/antimatter)

Parameters
[in]pdgIdpdgId to be checked
Returns
+1/-1/0 matter/antimatter/undefined

Definition at line 552 of file GenHFHadronMatcher.cc.

References funct::abs(), and isMesonPdgId().

Referenced by findHadronJets().

553 {
554  int flavourSign = pdgId / std::abs(pdgId);
555  // B mesons have opposite sign
556  if( isMesonPdgId(5, pdgId) ) flavourSign *= -1;
557  // Returning 0 for bb/cc resonances
558  if(pdgId % 1000 / 10 / 11 > 0) flavourSign = 0;
559 
560  return flavourSign;
561 }
int flavourSign(const int pdgId)
Sign of the flavour (matter/antimatter)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool isMesonPdgId(const int flavour, const int pdgId)
Check the pdgId if it represents a meson of particular flavour.
std::string GenHFHadronMatcher::getParticleName ( int  id) const
private
bool GenHFHadronMatcher::hasHadronDaughter ( const int  flavour,
const reco::Candidate thisParticle 
)
private

Check if the particle has bHadron among daughters.

Parameters
[in]flavourflavour of a hadron that is being searched (5-B, 4-C)
[in]thisParticlea particle that is to be analysed
Returns
whether the particle has a hadron among its daughters

Definition at line 572 of file GenHFHadronMatcher.cc.

References reco::Candidate::daughter(), createfilelist::int, isHadron(), gen::k, and reco::Candidate::numberOfDaughters().

Referenced by findHadronJets().

573 {
574  // Looping through daughters of the particle
575  bool hasDaughter = false;
576  for ( int k=0; k< (int)thisParticle->numberOfDaughters(); k++ ) {
577  if ( !isHadron( flavour, thisParticle->daughter(k) ) ) {
578  continue;
579  }
580  hasDaughter = true;
581  break;
582  }
583 
584  return hasDaughter;
585 }
virtual const Candidate * daughter(size_type i) const =0
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
int k[5][pyjets_maxn]
bool isHadron(const int flavour, const reco::Candidate *thisParticle)
Check the pdgId of a given particle if it is a hadron.
virtual size_type numberOfDaughters() const =0
number of daughters
int GenHFHadronMatcher::idInList ( std::vector< const reco::Candidate * >  particleList,
const reco::Candidate particle 
)
private

Check if the cpecified particle is already in the list of particles.

Parameters
[in]particleListlist of particles to be checked
[in]particleparticle that should be checked
Returns
the index of the particle in the list [-1 if particle not found]

Definition at line 454 of file GenHFHadronMatcher.cc.

References spr::find(), and position.

Referenced by analyzeMothers().

455 {
456  const unsigned int position = std::find(particleList.begin(), particleList.end(), particle) - particleList.begin();
457  if( position >= particleList.size() ) return -1;
458 
459  return position;
460 }
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
static int position[264][3]
Definition: ReadPGInfo.cc:509
int GenHFHadronMatcher::idInList ( std::vector< int >  list,
const int  value 
)
private

Definition at line 462 of file GenHFHadronMatcher.cc.

References spr::find(), position, and relativeConstraints::value.

463 {
464  const unsigned int position = std::find(list.begin(), list.end(), value) - list.begin();
465  if( position >= list.size() ) return -1;
466 
467  return position;
468 }
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
static int position[264][3]
Definition: ReadPGInfo.cc:509
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger list("!*","!HLTx*"if it matches 2 triggers or more) will accept the event if all the matching triggers are FAIL.It will reject the event if any of the triggers are PASS or EXCEPTION(this matches the behavior of"!*"before the partial wildcard feature was incorporated).Triggers which are in the READY state are completely ignored.(READY should never be returned since the trigger paths have been run
bool GenHFHadronMatcher::isBaryonPdgId ( const int  flavour,
const int  pdgId 
)
private

Check the pdgId if it represents a baryon of particular flavour.

Parameters
[in]flavourflavour of a hadron that is being searched (5-B, 4-C)
[in]pdgIdpdgId to be checked
Returns
true if the pdgId represents a baryon of specified flavour

Definition at line 533 of file GenHFHadronMatcher.cc.

References funct::abs().

Referenced by isHadronPdgId().

534 {
535  const int flavour_abs = std::abs(flavour);
536  if(flavour_abs != 5 && flavour_abs != 4) return false;
537  const int pdgId_abs = std::abs(pdgId);
538 
539  if ( pdgId_abs/1000 != flavour_abs) return false;
540 
541  return true;
542 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool GenHFHadronMatcher::isHadron ( const int  flavour,
const reco::Candidate thisParticle 
)
private

Check the pdgId of a given particle if it is a hadron.

Parameters
[in]flavourflavour of a hadron that is being searched (5-B, 4-C)
[in]thisParticlea particle that is to be analysed
Returns
whether the particle is a hadron of specified flavour

Definition at line 479 of file GenHFHadronMatcher.cc.

References isHadronPdgId(), and reco::Candidate::pdgId().

Referenced by findHadronJets(), and hasHadronDaughter().

480 {
481  return isHadronPdgId(flavour, thisParticle->pdgId());
482 }
bool isHadronPdgId(const int flavour, const int pdgId)
Check the pdgId if it represents a hadron of particular flavour.
virtual int pdgId() const =0
PDG identifier.
bool GenHFHadronMatcher::isHadronPdgId ( const int  flavour,
const int  pdgId 
)
private

Check the pdgId if it represents a hadron of particular flavour.

Parameters
[in]flavourflavour of a hadron that is being searched (5-B, 4-C)
[in]pdgIdpdgId to be checked
Returns
true if the pdgId represents a hadron of specified flavour

Definition at line 493 of file GenHFHadronMatcher.cc.

References isBaryonPdgId(), and isMesonPdgId().

Referenced by findInMothers(), and isHadron().

494 {
495  if( isBaryonPdgId(flavour, pdgId) || isMesonPdgId(flavour, pdgId) ) return true;
496 
497  return false;
498 }
bool isMesonPdgId(const int flavour, const int pdgId)
Check the pdgId if it represents a meson of particular flavour.
bool isBaryonPdgId(const int flavour, const int pdgId)
Check the pdgId if it represents a baryon of particular flavour.
bool GenHFHadronMatcher::isMesonPdgId ( const int  flavour,
const int  pdgId 
)
private

Check the pdgId if it represents a meson of particular flavour.

Parameters
[in]flavourflavour of a hadron that is being searched (5-B, 4-C)
[in]pdgIdpdgId to be checked
Returns
true if the pdgId represents a meson of specified flavour

Definition at line 509 of file GenHFHadronMatcher.cc.

References funct::abs(), and noBBbarResonances_.

Referenced by flavourSign(), and isHadronPdgId().

510 {
511  const int flavour_abs = std::abs(flavour);
512  if(flavour_abs != 5 && flavour_abs != 4) return false;
513  const int pdgId_abs = std::abs(pdgId);
514 
515  if( pdgId_abs/100%10 != flavour_abs) return false;
516  // Excluding baryons
517  if ( pdgId_abs/1000 == flavour_abs) return false;
518  // Excluding bb/cc resonances if required
519  if ( noBBbarResonances_ && pdgId_abs/10%100 == 11*flavour_abs ) return false;
520 
521  return true;
522 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool GenHFHadronMatcher::isNeutralPdg ( int  pdgId)
private

Check whether a given pdgId represents neutral particle.

Parameters
[in]pdgId
[in]thisParticle- a particle that is to be analysed
Returns
if the particle has a hadron among its daughters

Definition at line 873 of file GenHFHadronMatcher.cc.

References funct::abs(), and spr::find().

Referenced by fixExtraSameFlavours().

874 {
875  const int neutralPdgs_array[] = {9, 21, 22, 23, 25};
876  const std::vector<int> neutralPdgs( neutralPdgs_array, neutralPdgs_array + sizeof(neutralPdgs_array) / sizeof(int) );
877  if( std::find( neutralPdgs.begin(), neutralPdgs.end(), std::abs(pdgId) ) == neutralPdgs.end() ) return false;
878 
879  return true;
880 }
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void GenHFHadronMatcher::produce ( edm::Event evt,
const edm::EventSetup setup 
)
privatevirtual

Definition at line 190 of file GenHFHadronMatcher.cc.

References findHadronJets(), flavourStr_, GenHFHadronMatcher_cfi::genParticles, genParticlesToken_, edm::Event::getByToken(), edm::EventSetup::getData(), GenHFHadronMatcher_cfi::jetFlavourInfos, jetFlavourInfosToken_, eostools::move(), pdt_, edm::Handle< T >::product(), and edm::Event::put().

191 {
192  setup.getData ( pdt_ );
193 
194  using namespace edm;
195 
197  evt.getByToken(genParticlesToken_, genParticles);
198 
200  evt.getByToken(jetFlavourInfosToken_, jetFlavourInfos);
201 
202  // Defining adron matching variables
203  auto hadMothers = std::make_unique<std::vector<reco::GenParticle> >();
204  auto hadMothersIndices = std::make_unique<std::vector<std::vector<int> > >();
205  auto hadIndex = std::make_unique<std::vector<int> >();
206  auto hadFlavour = std::make_unique<std::vector<int> >();
207  auto hadJetIndex = std::make_unique<std::vector<int> >();
208  auto hadLeptonIndex = std::make_unique<std::vector<int> >();
209  auto hadLeptonHadIndex = std::make_unique<std::vector<int> >();
210  auto hadLeptonViaTau = std::make_unique<std::vector<int> >();
211  auto hadFromTopWeakDecay = std::make_unique<std::vector<int> >();
212  auto hadBHadronId = std::make_unique<std::vector<int> >();
213 
214  *hadJetIndex = findHadronJets (genParticles.product(), jetFlavourInfos.product(), *hadIndex, *hadMothers, *hadMothersIndices, *hadLeptonIndex, *hadLeptonHadIndex, *hadLeptonViaTau, *hadFlavour, *hadFromTopWeakDecay, *hadBHadronId );
215 
216  // Putting products to the event
217  evt.put(std::move(hadMothers), "gen"+flavourStr_+"HadPlusMothers" );
218  evt.put(std::move(hadMothersIndices), "gen"+flavourStr_+"HadPlusMothersIndices" );
219  evt.put(std::move(hadIndex), "gen"+flavourStr_+"HadIndex" );
220  evt.put(std::move(hadFlavour), "gen"+flavourStr_+"HadFlavour" );
221  evt.put(std::move(hadJetIndex), "gen"+flavourStr_+"HadJetIndex" );
222  evt.put(std::move(hadLeptonIndex), "gen"+flavourStr_+"HadLeptonIndex" );
223  evt.put(std::move(hadLeptonHadIndex), "gen"+flavourStr_+"HadLeptonHadronIndex" );
224  evt.put(std::move(hadLeptonViaTau), "gen"+flavourStr_+"HadLeptonViaTau" );
225  evt.put(std::move(hadFromTopWeakDecay),"gen"+flavourStr_+"HadFromTopWeakDecay" );
226  evt.put(std::move(hadBHadronId), "gen"+flavourStr_+"HadBHadronId" );
227 }
edm::ESHandle< ParticleDataTable > pdt_
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:122
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:460
void getData(T &iHolder) const
Definition: EventSetup.h:79
edm::EDGetTokenT< reco::JetFlavourInfoMatchingCollection > jetFlavourInfosToken_
T const * product() const
Definition: Handle.h:81
edm::EDGetTokenT< reco::GenParticleCollection > genParticlesToken_
HLT enums.
def move(src, dest)
Definition: eostools.py:510
std::vector< int > findHadronJets(const reco::GenParticleCollection *genParticles, const reco::JetFlavourInfoMatchingCollection *jetFlavourInfos, std::vector< int > &hadIndex, std::vector< reco::GenParticle > &hadMothersGenPart, std::vector< std::vector< int > > &hadMothersIndices, std::vector< int > &hadLeptonIndex, std::vector< int > &hadLeptonHadIndex, std::vector< int > &hadLeptonViaTau, std::vector< int > &hadFlavour, std::vector< int > &hadFromTopWeakDecay, std::vector< int > &hadBHadronId)
identify the jets that contain b-hadrons
bool GenHFHadronMatcher::putMotherIndex ( std::vector< std::vector< int > > &  hadMothersIndices,
int  partIndex,
int  mothIndex 
)
private

puts mother index to the list of mothers of particle, if it isn't there already

Parameters
[in]hadMothersIndicesvector of indices of mothers for each particle
[in]partIndexindex of the particle for which the mother index should be stored
[in]mothIndexindex of mother that should be stored for current particle
Returns
whether the particle index was alreade in the list

Definition at line 699 of file GenHFHadronMatcher.cc.

References createfilelist::int, and gen::k.

Referenced by analyzeMothers().

700 {
701  // Putting vector of mothers indices for the given particle
702  bool inList=false;
703  if ( partIndex<0 ) {
704  return false;
705  }
706 
707  while ( (int)hadMothersIndices.size() <= partIndex ) { // If there is no list of mothers for current particle yet
708  std::vector<int> mothersIndices;
709  hadMothersIndices.push_back ( mothersIndices );
710  }
711 
712  std::vector<int> *hadMotherIndices = &hadMothersIndices.at(partIndex);
713  // Removing other mothers if particle must have no mothers
714  if ( mothIndex == -1 ) {
715  hadMotherIndices->clear();
716  } else {
717  // Checking if current mother is already in the list of theParticle's mothers
718  for ( int k = 0; k < (int)hadMotherIndices->size(); k++ ) {
719  if ( hadMotherIndices->at(k) != mothIndex && hadMotherIndices->at(k) != -1 ) {
720  continue;
721  }
722  inList=true;
723  break;
724  }
725  }
726  // Adding current mother to the list of mothers of this particle
727  if ( !inList ) {
728  hadMotherIndices->push_back(mothIndex);
729  }
730 
731  return inList;
732 }
int k[5][pyjets_maxn]

Member Data Documentation

int GenHFHadronMatcher::flavour_
private

Definition at line 103 of file GenHFHadronMatcher.cc.

Referenced by findHadronJets(), and GenHFHadronMatcher().

std::string GenHFHadronMatcher::flavourStr_
private

Definition at line 107 of file GenHFHadronMatcher.cc.

Referenced by GenHFHadronMatcher(), and produce().

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

Definition at line 101 of file GenHFHadronMatcher.cc.

Referenced by produce().

edm::EDGetTokenT<reco::JetFlavourInfoMatchingCollection> GenHFHadronMatcher::jetFlavourInfosToken_
private

Definition at line 102 of file GenHFHadronMatcher.cc.

Referenced by produce().

bool GenHFHadronMatcher::noBBbarResonances_
private

Definition at line 104 of file GenHFHadronMatcher.cc.

Referenced by GenHFHadronMatcher(), and isMesonPdgId().

bool GenHFHadronMatcher::onlyJetClusteredHadrons_
private

Definition at line 105 of file GenHFHadronMatcher.cc.

Referenced by findHadronJets(), and GenHFHadronMatcher().

edm::ESHandle<ParticleDataTable> GenHFHadronMatcher::pdt_
private

Definition at line 109 of file GenHFHadronMatcher.cc.

Referenced by produce().