CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TauGenJetProducer.cc
Go to the documentation of this file.
2 
4 
8 
12 
14 
15 using namespace std;
16 using namespace edm;
17 using namespace reco;
18 
20 {
21  inputTagGenParticles_
22  = iConfig.getParameter<InputTag>("GenParticles");
23  tokenGenParticles_
24  = consumes<GenParticleCollection>(inputTagGenParticles_);
25 
26  includeNeutrinos_
27  = iConfig.getParameter<bool>("includeNeutrinos");
28 
29  verbose_ =
30  iConfig.getUntrackedParameter<bool>("verbose",false);
31 
32  produces<GenJetCollection>();
33 }
34 
36 
38 
40  const EventSetup& iSetup) {
41 
43 
44  bool found = iEvent.getByToken( tokenGenParticles_, genParticles);
45 
46  if ( !found ) {
47  std::ostringstream err;
48  err<<" cannot get collection: "
49  <<inputTagGenParticles_<<std::endl;
50  edm::LogError("TauGenJetProducer")<<err.str();
51  throw cms::Exception( "MissingProduct", err.str());
52  }
53 
54  std::auto_ptr<GenJetCollection>
55  pOutVisTaus(new GenJetCollection());
56 
57  using namespace GenParticlesHelper;
58 
59  GenParticleRefVector allStatus2Taus;
60  findParticles( *genParticles,
61  allStatus2Taus, 15, 2);
62 
63  for ( IGR iTau=allStatus2Taus.begin(); iTau!=allStatus2Taus.end(); ++iTau ) {
64 
65  // look for all status 1 (stable) descendents
66  GenParticleRefVector descendents;
67  findDescendents( *iTau, descendents, 1);
68 
69  // CV: skip status 2 taus that radiate-off a photon
70  // --> have a status 2 tau lepton in the list of descendents
71  GenParticleRefVector status2TauDaughters;
72  findDescendents( *iTau, status2TauDaughters, 2, 15 );
73  if ( status2TauDaughters.size() > 0 ) continue;
74 
75  // loop on descendents, and take all except neutrinos
76  math::XYZTLorentzVector sumVisMom;
78  Jet::Constituents constituents;
79 
80  if(verbose_)
81  cout<<"tau "<<(*iTau)<<endl;
82 
83  for(IGR igr = descendents.begin();
84  igr!= descendents.end(); ++igr ) {
85 
86  int absPdgId = abs((*igr)->pdgId());
87 
88  // neutrinos
89  if(!includeNeutrinos_ ) {
90  if( absPdgId == 12 ||
91  absPdgId == 14 ||
92  absPdgId == 16 )
93  continue;
94  }
95 
96  if(verbose_)
97  cout<<"\t"<<(*igr)<<endl;
98 
99  charge += (*igr)->charge();
100  sumVisMom += (*igr)->p4();
101 
102  // need to convert the vector of reference to the constituents
103  // to a vector of pointers to build the genjet
104  constituents.push_back( refToPtr( *igr) );
105  }
106 
107  math::XYZPoint vertex;
109 
110  GenJet jet( sumVisMom, vertex, specific, constituents);
111 
112  if (charge != (*iTau)->charge() )
113  std::cout<<" charge of Tau: " << (*iTau) << " not equal to charge of sum of charge of all descendents. " << std::cout;
114 
115  jet.setCharge(charge);
116  pOutVisTaus->push_back( jet );
117 
118  }
119  iEvent.put( pOutVisTaus );
120 }
121 
122 
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
dictionary specific
virtual void setCharge(Charge q) GCC11_FINAL
set electric charge
Ptr< typename C::value_type > refToPtr(Ref< C, typename C::value_type, refhelper::FindUsingAdvance< C, typename C::value_type > > const &ref)
Definition: RefToPtr.h:18
void findParticles(const reco::GenParticleCollection &sourceParticles, reco::GenParticleRefVector &particleRefs, int pdgId, int status)
find all particles of a given pdgId and status
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:434
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
std::vector< GenJet > GenJetCollection
collection of GenJet objects
std::vector< Constituent > Constituents
Definition: Jet.h:23
int Charge
electric charge type
Definition: Particle.h:25
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:249
double charge(const std::vector< uint8_t > &Ampls)
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:244
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
int iEvent
Definition: GenABIO.cc:243
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:116
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Jets made from MC generator particles.
Definition: GenJet.h:24
virtual void beginJob()
builds a GenJet from the visible daughters of each status 2 tau in the event.
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
virtual void produce(edm::Event &, const edm::EventSetup &)
size_type size() const
Size of the RefVector.
Definition: RefVector.h:89
TauGenJetProducer(const edm::ParameterSet &)
tuple cout
Definition: gather_cfg.py:121
void findDescendents(const reco::GenParticleRef &base, reco::GenParticleRefVector &descendents, int status, int pdgId=0)
find all descendents of a given status and pdgId (recursive)