CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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
 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 ()
 

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...
 
virtual void beginJob ()
 
virtual void beginLuminosityBlock (edm::LuminosityBlock &, edm::EventSetup const &)
 
virtual void beginRun (edm::Run &, edm::EventSetup const &)
 
bool checkForLoop (std::vector< const reco::Candidate * > &particleChain, const reco::Candidate *particle)
 
virtual void endJob ()
 
virtual void endLuminosityBlock (edm::LuminosityBlock &, edm::EventSetup const &)
 
virtual void endRun (edm::Run &, edm::EventSetup const &)
 
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::GenParticleCollection
genParticlesToken_
 
edm::EDGetTokenT
< reco::JetFlavourInfoMatchingCollection
jetFlavourInfosToken_
 
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
 
- 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

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 139 of file GenHFHadronMatcher.cc.

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

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

170 {
171 }

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 644 of file GenHFHadronMatcher.cc.

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

Referenced by findHadronJets().

645 {
646  // Getting the index of the particle which is a hadron in the first call
647  int hadronIndex=-1; // Index of the hadron that is returned by this function
648  int index = idInList( hadMothers, thisParticle );
649  if ( index<0 ) { // If hadron is not in the list of mothers yet
650  hadMothers.push_back ( thisParticle );
651  hadronIndex=hadMothers.size()-1;
652  } else { // If hadron is in the list of mothers already
653  hadronIndex=index;
654  }
655 
656  int partIndex = -1; // Index of particle being checked in the list of mothers
657  partIndex = idInList( hadMothers, thisParticle );
658 
659  // Checking whether this particle is already in the chain of analyzed particles in order to identify a loop
660  bool isLoop = false;
661  if ( !analyzedParticles ) {
662  analyzedParticles = new std::set<const reco::Candidate*>;
663  }
664  for ( unsigned int i=0; i<analyzedParticles->size(); i++ ) {
665  if ( analyzedParticles->count ( thisParticle ) <=0 ) {
666  continue;
667  }
668  isLoop = true;
669  break;
670  }
671 
672  // If a loop has been detected
673  if ( isLoop ) {
674  if ( prevPartIndex>=0 ) {
675  putMotherIndex ( hadMothersIndices, prevPartIndex, -1 ); // Setting mother index of previous particle to -1
676  }
677  return hadronIndex; // Stopping further processing of the current chain
678  }
679  analyzedParticles->insert ( thisParticle );
680 
681  // Putting the mothers to the list of mothers
682  for ( size_t iMother = 0; iMother < thisParticle->numberOfMothers(); ++iMother ) {
683  const reco::Candidate* mother = thisParticle->mother ( iMother );
684  int mothIndex = idInList( hadMothers, mother );
685  if ( mothIndex == partIndex && partIndex>=0 ) {
686  continue; // Skipping the mother that is its own daughter
687  }
688 
689  // If this mother isn't yet in the list and hadron or lepton is in the list
690  if ( mothIndex<0 ) {
691  hadMothers.push_back ( mother );
692  mothIndex=hadMothers.size()-1;
693  }
694  // If hadron has already been found in current chain and the mother isn't a duplicate of the particle being checked
695  if ( mothIndex!=partIndex && partIndex>=0 ) {
696  putMotherIndex ( hadMothersIndices, partIndex, mothIndex ); // Putting the index of mother for current particle
697  }
698  analyzeMothers ( mother, topDaughterQId, topBarDaughterQId, hadMothers, hadMothersIndices, analyzedParticles, partIndex );
699  // Setting the id of the particle which is a quark from the top decay
700  if(std::abs(mother->pdgId())==6) {
701  int& bId = mother->pdgId() < 0 ? topBarDaughterQId : topDaughterQId;
702  int thisFlav = std::abs(thisParticle->pdgId());
703  if( bId<0){
704  if(thisFlav <= 5) bId = partIndex;
705  } else {
706  int bIdFlav = std::abs(hadMothers.at(bId)->pdgId());
707  if( bIdFlav != 5 && thisFlav == 5) bId = partIndex;
708  else if( thisFlav == 5 && thisParticle->pt() > hadMothers.at(bId)->pt() ) bId = partIndex;
709  } // If daughter quark of the top not found yet
710  } // If the mother is a top quark and hadron has been found
711  } // End of loop over mothers
712 
713  analyzedParticles->erase ( thisParticle ); // Removing current particle from the current chain that is being analyzed
714 
715  if ( partIndex<0 ) {
716  return hadronIndex; // Safety check
717  }
718 
719  // Adding -1 to the list of mother indices for current particle if it has no mothers (for consistency between numbering of indices and mothers)
720  if ( (int)thisParticle->numberOfMothers() <=0) {
721  putMotherIndex ( hadMothersIndices, partIndex, -1 );
722  }
723 
724  return hadronIndex;
725 
726 }
int i
Definition: DBlmapReader.cc:9
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
virtual size_type numberOfMothers() const =0
number of mothers (zero or one in most of but not all the cases)
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 float pt() const =0
transverse momentum
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 int pdgId() const =0
PDG identifier.
void GenHFHadronMatcher::beginJob ( void  )
privatevirtual

Reimplemented from edm::EDProducer.

Definition at line 237 of file GenHFHadronMatcher.cc.

238 {
239 }
void GenHFHadronMatcher::beginLuminosityBlock ( edm::LuminosityBlock ,
edm::EventSetup const &   
)
privatevirtual

Definition at line 258 of file GenHFHadronMatcher.cc.

259 {
260 }
void GenHFHadronMatcher::beginRun ( edm::Run ,
edm::EventSetup const &   
)
privatevirtual

Definition at line 247 of file GenHFHadronMatcher.cc.

248 {
249 }
bool GenHFHadronMatcher::checkForLoop ( std::vector< const reco::Candidate * > &  particleChain,
const reco::Candidate particle 
)
private
void GenHFHadronMatcher::endJob ( void  )
privatevirtual

Reimplemented from edm::EDProducer.

Definition at line 242 of file GenHFHadronMatcher.cc.

243 {
244 }
void GenHFHadronMatcher::endLuminosityBlock ( edm::LuminosityBlock ,
edm::EventSetup const &   
)
privatevirtual

Definition at line 263 of file GenHFHadronMatcher.cc.

264 {
265 }
void GenHFHadronMatcher::endRun ( edm::Run ,
edm::EventSetup const &   
)
privatevirtual

Definition at line 253 of file GenHFHadronMatcher.cc.

254 {
255 }
void GenHFHadronMatcher::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

description of the run-time parameters

Definition at line 179 of file GenHFHadronMatcher.cc.

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

180 {
182  desc.add<edm::InputTag>("genParticles")->setComment( "Collection of GenParticle objects which contains all particles produced in the event" );
183  desc.add<edm::InputTag>("jetFlavourInfos")->setComment( "Output from the JetFlavour tool. Contains information about partons/hadrons/leptons associated to jets" );
184  desc.add<bool> ( "noBBbarResonances",true )->setComment ( "Whether resonances should not be treated as hadrons" );
185  desc.add<bool> ( "onlyJetClusteredHadrons",false )->setComment ( "Whether only hadrons that are matched to jets should be analysed. Runs x1000 faster in Sherpa" );
186  desc.add<int> ( "flavour",5 )->setComment ( "Flavour of weakly decaying hadron that should be analysed (4-c, 5-b)" );
187  descriptions.add ( "matchGenHFHadron",desc );
188 }
ParameterDescriptionBase * add(U const &iLabel, T const &value)
void add(std::string const &label, ParameterSetDescription const &psetDescription)
volatile std::atomic< bool > shutdown_flag false
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 291 of file GenHFHadronMatcher.cc.

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

Referenced by produce().

298 {
299  std::vector<int> hadJetIndex;
300  std::vector<const reco::Candidate*> hadMothersCand;
301 
302  int topDaughterQId = -1;
303  int topBarDaughterQId = -1;
304 
305  // Looping over all jets to get hadrons associated to them
306  for(reco::JetFlavourInfoMatchingCollection::const_iterator i_info = jetFlavourInfos->begin(); i_info != jetFlavourInfos->end(); ++i_info){
307  reco::JetFlavourInfo jetInfo = i_info->second;
308  const int jetIndex = i_info - jetFlavourInfos->begin();
309  // Looping over each hadron associated with the jet and finding its origin
310  const reco::GenParticleRefVector& hadronsInJet = flavour_==5 ? jetInfo.getbHadrons() : flavour_==4 ? jetInfo.getcHadrons() : reco::GenParticleRefVector();
311  for(reco::GenParticleRefVector::const_iterator hadron = hadronsInJet.begin(); hadron != hadronsInJet.end(); ++hadron) {
312  // Check that the hadron satisfies criteria configured in the module
313  if(!isHadron ( flavour_, (&**hadron) )) continue;
314  if(hasHadronDaughter ( flavour_, (reco::Candidate*)(&**hadron) )) continue;
315  // Scanning the chain starting from the hadron
316  int hadronIndex = analyzeMothers ( (reco::Candidate*)(&**hadron), topDaughterQId, topBarDaughterQId, hadMothersCand, hadMothersIndices, 0, -1 );
317  // Storing the index of the hadron to the list
318  hadIndex.push_back ( hadronIndex );
319  hadJetIndex.push_back ( jetIndex ); // Putting jet index to the result list
320  }
321  } // End of loop over jets
322 
323  // Access all hadrons which are not associated with jets, if requested
325  for(reco::GenParticleCollection::const_iterator i_particle = genParticles->begin(); i_particle != genParticles->end(); ++i_particle){
326  const reco::GenParticle* thisParticle = &*i_particle;
327  if(!isHadron(flavour_, thisParticle)) continue;
328  // Skipping the hadron if it was already found directly from jets
329  if(std::find(hadMothersCand.begin(), hadMothersCand.end(), thisParticle) != hadMothersCand.end()) continue;
330 
331  // Scanning the chain starting from the hadron
332  int hadronIndex = analyzeMothers ( thisParticle, topDaughterQId, topBarDaughterQId, hadMothersCand, hadMothersIndices, 0, -1 );
333  // Storing the index of the hadron to the list
334  hadIndex.push_back ( hadronIndex );
335  hadJetIndex.push_back ( -1 ); // Jet index undefined
336  }
337  }
338 
339  // Transfering Candidates to the list of processed particles for further analysis
340  for ( int i=0; i< (int)hadMothersCand.size(); i++ ) {
341  const reco::GenParticle* particle = dynamic_cast<const reco::GenParticle*>( hadMothersCand.at(i) );
342  hadMothers.push_back(*particle);
343  }
344 
345  // Adding leptons from hadron decays
346  for(reco::GenParticleCollection::const_iterator i_particle = genParticles->begin(); i_particle != genParticles->end(); ++i_particle){
347  const reco::GenParticle lepton = *i_particle;
348  const int pdg_abs = lepton.pdgId();
349  // Skipping if not a lepton: e/mu
350  if(pdg_abs != 11 && pdg_abs != 13) continue;
351  bool leptonViaTau = false;
352  const reco::Candidate* leptonMother = lepton.mother();
353  if(!leptonMother) continue;
354  // Taking next mother if direct mother is a tau
355  if(std::abs(leptonMother->pdgId()) == 15) {
356  leptonViaTau = 1;
357  leptonMother = leptonMother->mother();
358  }
359  // Skipping this lepton if its mother is not a proper hadron
360  if(!isHadron(flavour_, leptonMother)) continue;
361  // Finding the index of this hadron in the list of analysed particles
362  size_t leptonHadronParticleIndex = std::find(hadMothersCand.begin(), hadMothersCand.end(), leptonMother) - hadMothersCand.begin();
363  if(leptonHadronParticleIndex >= hadMothersCand.size()) continue;
364  // Finding the actual hadron index among those that were found
365  size_t leptonHadronIndex = std::find(hadIndex.begin(), hadIndex.end(), leptonHadronParticleIndex) - hadIndex.begin();
366  if(leptonHadronIndex >= hadIndex.size()) continue;
367  // Putting the lepton, its index and hadron index to the corresponding lists
368  hadMothers.push_back(lepton);
369  const int leptonIndex = hadMothersCand.size()-1;
370  hadLeptonIndex.push_back(leptonIndex);
371  hadLeptonViaTau.push_back(leptonViaTau);
372  hadLeptonHadIndex.push_back(leptonHadronIndex);
373  }
374 
375  // Checking mothers of hadrons in order to assign flags (where the hadron comes from)
376  unsigned int nHad = hadIndex.size();
377 
378  std::vector<std::vector<int> > LastQuarkMotherIds;
379  std::vector<std::vector<int> > LastQuarkIds;
380  std::vector<int> lastQuarkIndices(nHad, -1);
381 
382  // Looping over all hadrons
383  for ( unsigned int hadNum=0; hadNum<nHad; hadNum++ ) {
384  unsigned int hadIdx = hadIndex.at(hadNum); // Index of hadron in the hadMothers
385 
386  std::vector <int> FirstQuarkId;
387  std::vector <int> LastQuarkId;
388  std::vector <int> LastQuarkMotherId;
389 
390  const int hadronFlavourSign = flavourSign( hadMothers.at(hadIdx).pdgId() );
391 
392  // Searching only first quark in the chain with the same flavour as hadron
393  if(hadronFlavourSign != 0) {
394  findInMothers( hadIdx, FirstQuarkId, hadMothersIndices, hadMothers, 0, hadronFlavourSign*flavour_, false, -1, 1, false );
395  }
396  // Searching for quarks with both flavours since it is a bb/cc resonance
397  else {
398  findInMothers( hadIdx, FirstQuarkId, hadMothersIndices, hadMothers, 0, flavour_, false, -1, 1, false );
399  findInMothers( hadIdx, FirstQuarkId, hadMothersIndices, hadMothers, 0, -1*flavour_, false, -1, 1, false );
400  }
401 
402  // Finding last quark for each first quark
403  for ( unsigned int qId=0; qId<FirstQuarkId.size(); qId++ ) {
404  // Identifying the flavour of the first quark to find the last quark of the same flavour
405  const int quarkFlavourSign = flavourSign( hadMothers.at(FirstQuarkId.at(qId)).pdgId() );
406  // Finding last quark of the hadron starting from the first quark
407  findInMothers( FirstQuarkId.at(qId), LastQuarkId, hadMothersIndices, hadMothers, 0, quarkFlavourSign*flavour_, false, -1, 2, false );
408  } // End of loop over all first quarks of the hadron
409 
410 
411  // Setting initial flavour of the hadron
412  int hadronFlavour = 0;
413 
414  // Initialising pairs of last quark index and its distance from the hadron (to sort by distance)
415  std::vector<std::pair<double, int> > lastQuark_dR_id_pairs;
416 
417  // Finding the closest quark in dR
418  for ( unsigned int qId = 0; qId < LastQuarkId.size(); qId++ ) {
419  int qIdx = LastQuarkId.at(qId);
420  // Calculating the dR between hadron and quark
421  float dR = deltaR ( hadMothers.at(hadIdx).eta(),hadMothers.at(hadIdx).phi(),hadMothers.at(qIdx).eta(),hadMothers.at(qIdx).phi() );
422 
423  std::pair<double, int> dR_hadId_pair(dR,qIdx);
424  lastQuark_dR_id_pairs.push_back(dR_hadId_pair);
425  } // End of loop over all last quarks of the hadron
426 
427  std::sort(lastQuark_dR_id_pairs.begin(), lastQuark_dR_id_pairs.end());
428 
429  if(lastQuark_dR_id_pairs.size() > 1) {
430  double dRratio = (lastQuark_dR_id_pairs.at(1).first - lastQuark_dR_id_pairs.at(0).first)/lastQuark_dR_id_pairs.at(1).first;
431  int qIdx_closest = lastQuark_dR_id_pairs.at(0).second;
432  LastQuarkId.clear();
433  if(dRratio>0.5) LastQuarkId.push_back(qIdx_closest);
434  else for(std::pair<double, int> qIdDrPair : lastQuark_dR_id_pairs) LastQuarkId.push_back(qIdDrPair.second);
435  }
436  for(int qIdx : LastQuarkId) {
437  int qmIdx = hadMothersIndices.at(qIdx).at(0);
438  LastQuarkMotherId.push_back(qmIdx);
439  }
440 
441  if((int)LastQuarkId.size()>0) lastQuarkIndices.at(hadNum) = 0; // Setting the first quark in array as a candidate if it exists
442 
443  LastQuarkIds.push_back( LastQuarkId );
444 
445  LastQuarkMotherIds.push_back ( LastQuarkMotherId );
446 
447  if(LastQuarkMotherId.size()<1) {
448  hadronFlavour = 0;
449  } else {
450  int qIdx = LastQuarkId.at( lastQuarkIndices.at(hadNum) );
451  int qFlav = flavourSign( hadMothers.at(qIdx).pdgId() );
452  hadronFlavour = qFlav*std::abs( hadMothers.at( LastQuarkMotherId.at( lastQuarkIndices.at(hadNum) ) ).pdgId() );
453  }
454  hadFlavour.push_back(hadronFlavour); // Adding hadron flavour to the list of flavours
455 
456  // Checking whether hadron comes from the Top weak decay
457  int isFromTopWeakDecay = 1;
458  std::vector <int> checkedParticles;
459  if(hadFlavour.at(hadNum) != 0) {
460  int lastQIndex = LastQuarkId.at(lastQuarkIndices.at(hadNum));
461  bool fromTB = topDaughterQId >= 0 ? findInMothers( lastQIndex, checkedParticles, hadMothersIndices, hadMothers, -1, 0, false, topDaughterQId, 2, false ) >= 0 : false;
462  checkedParticles.clear();
463  bool fromTbarB = topBarDaughterQId >= 0 ? findInMothers( lastQIndex, checkedParticles, hadMothersIndices, hadMothers, -1, 0, false, topBarDaughterQId, 2, false) >= 0 : false;
464  checkedParticles.clear();
465  if(!fromTB && !fromTbarB) {
466  isFromTopWeakDecay = 0;
467  }
468  } else isFromTopWeakDecay = 2;
469  hadFromTopWeakDecay.push_back(isFromTopWeakDecay);
470  int bHadronMotherId = findInMothers( hadIdx, checkedParticles, hadMothersIndices, hadMothers, 0, 555555, true, -1, 1, false );
471  hadBHadronId.push_back(bHadronMotherId);
472 
473 
474  if(LastQuarkMotherId.size()>0) {
475  std::set<int> checkedHadronIds;
476  fixExtraSameFlavours(hadNum, hadIndex, hadMothers, hadMothersIndices, hadFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadFlavour, checkedHadronIds, 0);
477  }
478 
479  } // End of loop over all hadrons
480 
481  return hadJetIndex;
482 }
int i
Definition: DBlmapReader.cc:9
virtual int pdgId() const
PDG identifier.
virtual const Candidate * mother(size_type i=0) const =0
return pointer to mother
transient_vector_type::const_iterator const_iterator
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
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:249
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:244
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...
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
virtual int pdgId() const =0
PDG identifier.
edm::RefVector< GenParticleCollection > GenParticleRefVector
vector of reference to GenParticle in the same collection
volatile std::atomic< bool > shutdown_flag false
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 791 of file GenHFHadronMatcher.cc.

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

Referenced by findHadronJets().

794 {
795  int foundStopId = -1;
796  int pdg_1 = hadMothers.at(idx).pdgId();
797 
798  if ( (int)hadMothersIndices.size() <= idx ) {
799  if ( verbose ) {
800  printf ( " Stopping checking particle %d. No mothers are stored.\n",idx );
801  }
802  return -1; // Skipping if no mothers are stored for this particle
803  }
804 
805  if(std::abs(hadMothers.at( idx ).pdgId()) > 10 && std::abs(hadMothers.at( idx ).pdgId()) < 19) printf("Lepton: %d\n", hadMothers.at( idx ).pdgId());
806 
807  std::vector<int> mothers = hadMothersIndices.at(idx);
808  unsigned int nMothers = mothers.size();
809  bool isCorrect=false; // Whether current particle is what is being searched
810  if ( verbose ) {
811  if ( std::abs( hadMothers.at(idx).pdgId() ) ==2212 ) {
812  printf ( "Chk: %d\tpdg: %d\tstatus: %d",idx, hadMothers.at(idx).pdgId(), hadMothers.at(idx).status() );
813  } else {
814  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() );
815  }
816  }
817  bool hasCorrectMothers = true;
818  if(nMothers<1) hasCorrectMothers=false; else if(mothers.at(0)<0) hasCorrectMothers=false;
819  if(!hasCorrectMothers) {
820  if(verbose) printf(" NO CORRECT MOTHER\n");
821  return -1;
822  }
823  // Stopping if the particular particle has been found
824  if(stopId>=0 && idx == stopId) return idx;
825 
826  // Stopping if the hadron of particular flavour has been found
827  if(pdgId%111111==0 && pdgId!=0) {
828  if(isHadronPdgId(pdgId/111111, hadMothers.at(idx).pdgId())) {
829  return idx;
830  }
831  }
832 
833  // Checking whether current mother satisfies selection criteria
834  if ( ( ( hadMothers.at(idx).pdgId() == pdgId && pdgAbs==false )
835  || ( std::abs( hadMothers.at(idx).pdgId() ) == std::abs( pdgId ) && pdgAbs==true ) )
836  && ( hadMothers.at(idx).status() == status || status==0 )
837  && hasCorrectMothers ) {
838  isCorrect=true;
839  // Adding to the list of candidates if not there and if mother of this quark has correct flavour sign
840  const bool inList = std::find(mothChains.begin(), mothChains.end(), idx) != mothChains.end();
841  if ( !inList && mothers.at(0) >= 0 && ( hadMothers.at(idx).pdgId()*pdgId > 0 || !pdgAbs ) ) {
842  if ( firstLast==0 || firstLast==1 ) {
843  mothChains.push_back(idx);
844  }
845  if ( verbose ) {
846  printf ( " *" );
847  }
848  }
849  if ( verbose ) {
850  printf ( " +++" );
851  }
852  }
853  if ( verbose ) {
854  printf ( "\n" );
855  }
856 
857  if ( isCorrect && firstLast==1 ) {
858  return -1; // Stopping if only the first particle in the chain is looked for
859  }
860 
861  // Checking next level mothers
862  unsigned int nDifferingMothers = 0;
863  for ( unsigned int i = 0; i < nMothers; i++ ) {
864  int idx2 = mothers.at(i);
865  if ( idx2 < 0 ) {
866  if(verbose) printf("^^^ Has no mother\n");
867  continue; // Skipping if mother's id is -1 (no mother), that means current particle is a proton
868  }
869  if ( idx2 == idx ) {
870  if(verbose) printf("^^^ Stored as its own mother\n");
871  continue; // Skipping if particle is stored as its own mother
872  }
873  int pdg_2 = hadMothers[idx2].pdgId();
874  // Inverting the flavour if bb oscillation detected
875  if ( isHadronPdgId(pdgId, pdg_1) && isHadronPdgId(pdgId, pdg_2) && pdg_1*pdg_2 < 0 ) {
876  pdgId *= -1;
877  if(verbose) printf("######### Inverting flavour of the hadron\n");
878  }
879  // Counting how many mothers are different from this particle
880  if ( ( std::abs( pdg_2 ) != std::abs( pdgId ) && pdgAbs==true ) ||
881  ( pdg_2 != pdgId && pdgAbs == false ) ) {
882  nDifferingMothers++;
883  }
884 
885  // Checking next level mother
886  if ( verbose ) {
887  printf ( "Checking mother %d out of %d mothers (%d -> %d), looking for pdgId: %d\n",i,nMothers,idx, idx2, pdgId );
888  }
889  if(firstLast==2 && pdg_1 != pdg_2) continue; // Requiring the same flavour when finding the last quark
890  foundStopId = findInMothers ( idx2, mothChains, hadMothersIndices, hadMothers, status, pdgId, pdgAbs, stopId, firstLast, verbose );
891  }
892  // Storing this particle if all its mothers are of another type and the last of its kind should be stored
893  if(firstLast==2 && isCorrect && nDifferingMothers >= nMothers) {
894  if ( verbose ) {
895  printf ( "Checking particle %d once more to store it as the last quark\n",idx);
896  }
897  foundStopId = findInMothers ( idx, mothChains, hadMothersIndices, hadMothers, 0, pdgId, pdgAbs, stopId, 1, verbose );
898  }
899 
900  return foundStopId;
901 }
int i
Definition: DBlmapReader.cc:9
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:7
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
tuple idx
DEBUGGING if hasattr(process,&quot;trackMonIterativeTracking2012&quot;): print &quot;trackMonIterativeTracking2012 D...
volatile std::atomic< bool > shutdown_flag false
tuple status
Definition: ntuplemaker.py:245
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 936 of file GenHFHadronMatcher.cc.

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

Referenced by findHadronJets().

942 {
943  if(checkedHadronIds.count(hadId) != 0) return false; // Hadron already checked previously and should be skipped
944  checkedHadronIds.insert(hadId); // Putting hadron to the list of checked ones in this run
945 
946  if(lastQuarkIndex<0) return false;
947  if((int)LastQuarkIds.at(hadId).size()<lastQuarkIndex+1) return false;
948  int LastQuarkId = LastQuarkIds.at(hadId).at(lastQuarkIndex);
949  int LastQuarkMotherId = LastQuarkMotherIds.at( hadId ).at( lastQuarkIndex );
950  int qmFlav = hadMothers.at(LastQuarkId).pdgId() < 0 ? -1 : 1;
951  int hadFlavour = qmFlav*std::abs( hadMothers.at( LastQuarkMotherId ).pdgId() );
952  bool ambiguityResolved = true;
953  // If last quark has inconsistent flavour with its mother, setting the hadron flavour to gluon
954  if( (hadMothers.at(LastQuarkId).pdgId()*hadMothers.at(LastQuarkMotherId).pdgId() < 0 && !isNeutralPdg(hadMothers.at(LastQuarkMotherId).pdgId())) ||
955  // If particle not coming from the Top weak decay has Top flavour
956  (std::abs(hadronFlavour.at(hadId))==6 && isFromTopWeakDecay.at(hadId)==0) ) {
957  if((int)LastQuarkIds.at(hadId).size()>lastQuarkIndex+1) fixExtraSameFlavours(hadId, hadIndices, hadMothers, hadMothersIndices, isFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadronFlavour, checkedHadronIds, lastQuarkIndex+1);
958  else hadronFlavour.at(hadId) = qmFlav*21;
959  return true;
960  }
961 
962  int nSameFlavourHadrons = 0;
963  // Looping over all previous hadrons
964  for(unsigned int iHad = 0; iHad<hadronFlavour.size(); iHad++) {
965  if(iHad==hadId) continue;
966  int theLastQuarkIndex = lastQuarkIndices.at(iHad);
967  if(theLastQuarkIndex<0) continue;
968  if((int)LastQuarkMotherIds.at( iHad ).size() <= theLastQuarkIndex) continue;
969  int theLastQuarkMotherId = LastQuarkMotherIds.at( iHad ).at( theLastQuarkIndex );
970  int theHadFlavour = hadronFlavour.at(iHad);
971  // Skipping hadrons with different flavour
972  if(theHadFlavour==0 || std::abs(theHadFlavour)==21) continue;
973  if(theHadFlavour != hadFlavour || theLastQuarkMotherId != LastQuarkMotherId) continue;
974  ambiguityResolved = false;
975  nSameFlavourHadrons++;
976 
977  // Checking other b-quark mother candidates of this hadron
978  if((int)LastQuarkIds.at(hadId).size() > lastQuarkIndex+1) {
979  if(fixExtraSameFlavours(hadId, hadIndices, hadMothers, hadMothersIndices, isFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadronFlavour, checkedHadronIds, lastQuarkIndex+1) ) {
980  ambiguityResolved = true;
981  break;
982  }
983  } else
984  // Checking other b-quark mother candidates of the particular previous hadron
985  if((int)LastQuarkIds.at(iHad).size() > theLastQuarkIndex+1) {
986  if(fixExtraSameFlavours(iHad, hadIndices, hadMothers, hadMothersIndices, isFromTopWeakDecay, LastQuarkIds, LastQuarkMotherIds, lastQuarkIndices, hadronFlavour, checkedHadronIds, theLastQuarkIndex+1) ) {
987  ambiguityResolved = true;
988  break;
989  };
990  }
991 
992  } // End of loop over all previous hadrons
993 
994  checkedHadronIds.erase(hadId); // Removing the hadron from the list of checked hadrons
995  if(nSameFlavourHadrons>0 && !ambiguityResolved) {
996  hadronFlavour.at(hadId) = qmFlav*21;
997  return true;
998  }
999  lastQuarkIndices.at(hadId) = lastQuarkIndex;
1000  hadronFlavour.at(hadId) = hadFlavour;
1001  return true;
1002 }
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 591 of file GenHFHadronMatcher.cc.

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

Referenced by findHadronJets().

592 {
593  int flavourSign = pdgId / std::abs(pdgId);
594  // B mesons have opposite sign
595  if( isMesonPdgId(5, pdgId) ) flavourSign *= -1;
596  // Returning 0 for bb/cc resonances
597  if(pdgId % 1000 / 10 / 11 > 0) flavourSign = 0;
598 
599  return flavourSign;
600 }
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 611 of file GenHFHadronMatcher.cc.

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

Referenced by findHadronJets().

612 {
613  // Looping through daughters of the particle
614  bool hasDaughter = false;
615  for ( int k=0; k< (int)thisParticle->numberOfDaughters(); k++ ) {
616  if ( !isHadron( flavour, thisParticle->daughter(k) ) ) {
617  continue;
618  }
619  hasDaughter = true;
620  break;
621  }
622 
623  return hasDaughter;
624 }
virtual const Candidate * daughter(size_type i) const =0
return daughter at a given position, i = 0, ... numberOfDaughters() - 1 (read only mode) ...
virtual size_type numberOfDaughters() const =0
number of daughters
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.
int flavour(const Candidate &part)
Definition: pdgIdUtils.h:31
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 493 of file GenHFHadronMatcher.cc.

References spr::find(), and position.

Referenced by analyzeMothers().

494 {
495  const unsigned int position = std::find(particleList.begin(), particleList.end(), particle) - particleList.begin();
496  if( position >= particleList.size() ) return -1;
497 
498  return position;
499 }
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
static int position[264][3]
Definition: ReadPGInfo.cc:509
int GenHFHadronMatcher::idInList ( std::vector< int >  list,
const int  value 
)
private

Definition at line 501 of file GenHFHadronMatcher.cc.

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

502 {
503  const unsigned int position = std::find(list.begin(), list.end(), value) - list.begin();
504  if( position >= list.size() ) return -1;
505 
506  return position;
507 }
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
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 572 of file GenHFHadronMatcher.cc.

References funct::abs().

Referenced by isHadronPdgId().

573 {
574  const int flavour_abs = std::abs(flavour);
575  if(flavour_abs != 5 && flavour_abs != 4) return false;
576  const int pdgId_abs = std::abs(pdgId);
577 
578  if ( pdgId_abs/1000 != flavour_abs) return false;
579 
580  return true;
581 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int flavour(const Candidate &part)
Definition: pdgIdUtils.h:31
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 518 of file GenHFHadronMatcher.cc.

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

Referenced by findHadronJets(), and hasHadronDaughter().

519 {
520  return isHadronPdgId(flavour, thisParticle->pdgId());
521 }
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.
int flavour(const Candidate &part)
Definition: pdgIdUtils.h:31
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 532 of file GenHFHadronMatcher.cc.

References isBaryonPdgId(), and isMesonPdgId().

Referenced by findInMothers(), and isHadron().

533 {
534  if( isBaryonPdgId(flavour, pdgId) || isMesonPdgId(flavour, pdgId) ) return true;
535 
536  return false;
537 }
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.
int flavour(const Candidate &part)
Definition: pdgIdUtils.h:31
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 548 of file GenHFHadronMatcher.cc.

References funct::abs(), and noBBbarResonances_.

Referenced by flavourSign(), and isHadronPdgId().

549 {
550  const int flavour_abs = std::abs(flavour);
551  if(flavour_abs != 5 && flavour_abs != 4) return false;
552  const int pdgId_abs = std::abs(pdgId);
553 
554  if( pdgId_abs/100%10 != flavour_abs) return false;
555  // Excluding baryons
556  if ( pdgId_abs/1000 == flavour_abs) return false;
557  // Excluding bb/cc resonances if required
558  if ( noBBbarResonances_ && pdgId_abs/10%100 == 11*flavour_abs ) return false;
559 
560  return true;
561 }
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int flavour(const Candidate &part)
Definition: pdgIdUtils.h:31
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 912 of file GenHFHadronMatcher.cc.

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

Referenced by fixExtraSameFlavours().

913 {
914  const int neutralPdgs_array[] = {9, 21, 22, 23, 25};
915  const std::vector<int> neutralPdgs( neutralPdgs_array, neutralPdgs_array + sizeof(neutralPdgs_array) / sizeof(int) );
916  if( std::find( neutralPdgs.begin(), neutralPdgs.end(), std::abs(pdgId) ) == neutralPdgs.end() ) return false;
917 
918  return true;
919 }
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void GenHFHadronMatcher::produce ( edm::Event evt,
const edm::EventSetup setup 
)
privatevirtual

Implements edm::EDProducer.

Definition at line 197 of file GenHFHadronMatcher.cc.

References findHadronJets(), flavourStr_, genParticleCandidates2GenParticles_cfi::genParticles, genParticlesToken_, edm::Event::getByToken(), edm::EventSetup::getData(), jetFlavourInfosToken_, pdt_, and edm::Event::put().

198 {
199  setup.getData ( pdt_ );
200 
201  using namespace edm;
202 
205 
207  evt.getByToken(jetFlavourInfosToken_, jetFlavourInfos);
208 
209  // Defining adron matching variables
210  std::auto_ptr<std::vector<reco::GenParticle> > hadMothers ( new std::vector<reco::GenParticle> );
211  std::auto_ptr<std::vector<std::vector<int> > > hadMothersIndices ( new std::vector<std::vector<int> > );
212  std::auto_ptr<std::vector<int> > hadIndex ( new std::vector<int> );
213  std::auto_ptr<std::vector<int> > hadFlavour ( new std::vector<int> );
214  std::auto_ptr<std::vector<int> > hadJetIndex ( new std::vector<int> );
215  std::auto_ptr<std::vector<int> > hadLeptonIndex ( new std::vector<int> );
216  std::auto_ptr<std::vector<int> > hadLeptonHadIndex ( new std::vector<int> );
217  std::auto_ptr<std::vector<int> > hadLeptonViaTau( new std::vector<int> );
218  std::auto_ptr<std::vector<int> > hadFromTopWeakDecay ( new std::vector<int> );
219  std::auto_ptr<std::vector<int> > hadBHadronId ( new std::vector<int> );
220 
221  *hadJetIndex = findHadronJets (genParticles.product(), jetFlavourInfos.product(), *hadIndex, *hadMothers, *hadMothersIndices, *hadLeptonIndex, *hadLeptonHadIndex, *hadLeptonViaTau, *hadFlavour, *hadFromTopWeakDecay, *hadBHadronId );
222 
223  // Putting products to the event
224  evt.put ( hadMothers, "gen"+flavourStr_+"HadPlusMothers" );
225  evt.put ( hadMothersIndices, "gen"+flavourStr_+"HadPlusMothersIndices" );
226  evt.put ( hadIndex, "gen"+flavourStr_+"HadIndex" );
227  evt.put ( hadFlavour, "gen"+flavourStr_+"HadFlavour" );
228  evt.put ( hadJetIndex, "gen"+flavourStr_+"HadJetIndex" );
229  evt.put ( hadLeptonIndex, "gen"+flavourStr_+"HadLeptonIndex" );
230  evt.put ( hadLeptonHadIndex, "gen"+flavourStr_+"HadLeptonHadronIndex" );
231  evt.put ( hadLeptonViaTau, "gen"+flavourStr_+"HadLeptonViaTau" );
232  evt.put ( hadFromTopWeakDecay,"gen"+flavourStr_+"HadFromTopWeakDecay" );
233  evt.put ( hadBHadronId, "gen"+flavourStr_+"HadBHadronId" );
234 }
edm::ESHandle< ParticleDataTable > pdt_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:434
void getData(T &iHolder) const
Definition: EventSetup.h:67
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:116
edm::EDGetTokenT< reco::JetFlavourInfoMatchingCollection > jetFlavourInfosToken_
edm::EDGetTokenT< reco::GenParticleCollection > genParticlesToken_
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 738 of file GenHFHadronMatcher.cc.

References gen::k.

Referenced by analyzeMothers().

739 {
740  // Putting vector of mothers indices for the given particle
741  bool inList=false;
742  if ( partIndex<0 ) {
743  return false;
744  }
745 
746  while ( (int)hadMothersIndices.size() <= partIndex ) { // If there is no list of mothers for current particle yet
747  std::vector<int> mothersIndices;
748  hadMothersIndices.push_back ( mothersIndices );
749  }
750 
751  std::vector<int> *hadMotherIndices = &hadMothersIndices.at(partIndex);
752  // Removing other mothers if particle must have no mothers
753  if ( mothIndex == -1 ) {
754  hadMotherIndices->clear();
755  } else {
756  // Checking if current mother is already in the list of theParticle's mothers
757  for ( int k = 0; k < (int)hadMotherIndices->size(); k++ ) {
758  if ( hadMotherIndices->at(k) != mothIndex && hadMotherIndices->at(k) != -1 ) {
759  continue;
760  }
761  inList=true;
762  break;
763  }
764  }
765  // Adding current mother to the list of mothers of this particle
766  if ( !inList ) {
767  hadMotherIndices->push_back(mothIndex);
768  }
769 
770  return inList;
771 }
int k[5][pyjets_maxn]

Member Data Documentation

int GenHFHadronMatcher::flavour_
private

Definition at line 110 of file GenHFHadronMatcher.cc.

Referenced by findHadronJets(), and GenHFHadronMatcher().

std::string GenHFHadronMatcher::flavourStr_
private

Definition at line 114 of file GenHFHadronMatcher.cc.

Referenced by GenHFHadronMatcher(), and produce().

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

Definition at line 108 of file GenHFHadronMatcher.cc.

Referenced by produce().

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

Definition at line 109 of file GenHFHadronMatcher.cc.

Referenced by produce().

bool GenHFHadronMatcher::noBBbarResonances_
private

Definition at line 111 of file GenHFHadronMatcher.cc.

Referenced by GenHFHadronMatcher(), and isMesonPdgId().

bool GenHFHadronMatcher::onlyJetClusteredHadrons_
private

Definition at line 112 of file GenHFHadronMatcher.cc.

Referenced by findHadronJets(), and GenHFHadronMatcher().

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

Definition at line 116 of file GenHFHadronMatcher.cc.

Referenced by produce().