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 
13 
15 
16 using namespace std;
17 using namespace edm;
18 using namespace reco;
19 
21 {
22  inputTagGenParticles_
23  = iConfig.getParameter<InputTag>("GenParticles");
24 
25  includeNeutrinos_
26  = iConfig.getParameter<bool>("includeNeutrinos");
27 
28  verbose_ =
29  iConfig.getUntrackedParameter<bool>("verbose",false);
30 
31  produces<GenJetCollection>();
32 }
33 
35 
37 
39  const EventSetup& iSetup) {
40 
42 
43  bool found = iEvent.getByLabel( inputTagGenParticles_, genParticles);
44 
45  if ( !found ) {
46  std::ostringstream err;
47  err<<" cannot get collection: "
48  <<inputTagGenParticles_<<std::endl;
49  edm::LogError("TauGenJetProducer")<<err.str();
50  throw cms::Exception( "MissingProduct", err.str());
51  }
52 
53  std::auto_ptr<GenJetCollection>
54  pOutVisTaus(new GenJetCollection());
55 
56  using namespace GenParticlesHelper;
57 
58  GenParticleRefVector allStatus2Taus;
59  findParticles( *genParticles,
60  allStatus2Taus, 15, 2);
61 
62  for ( IGR iTau=allStatus2Taus.begin(); iTau!=allStatus2Taus.end(); ++iTau ) {
63 
64  // look for all status 1 (stable) descendents
65  GenParticleRefVector descendents;
66  findDescendents( *iTau, descendents, 1);
67 
68  // CV: skip status 2 taus that radiate-off a photon
69  // --> have a status 2 tau lepton in the list of descendents
70  GenParticleRefVector status2TauDaughters;
71  findDescendents( *iTau, status2TauDaughters, 2, 15 );
72  if ( status2TauDaughters.size() > 0 ) continue;
73 
74  // loop on descendents, and take all except neutrinos
75  math::XYZTLorentzVector sumVisMom;
77  Jet::Constituents constituents;
78 
79  if(verbose_)
80  cout<<"tau "<<(*iTau)<<endl;
81 
82  for(IGR igr = descendents.begin();
83  igr!= descendents.end(); ++igr ) {
84 
85  int absPdgId = abs((*igr)->pdgId());
86 
87  // neutrinos
88  if(!includeNeutrinos_ ) {
89  if( absPdgId == 12 ||
90  absPdgId == 14 ||
91  absPdgId == 16 )
92  continue;
93  }
94 
95  if(verbose_)
96  cout<<"\t"<<(*igr)<<endl;
97 
98  charge += (*igr)->charge();
99  sumVisMom += (*igr)->p4();
100 
101  // need to convert the vector of reference to the constituents
102  // to a vector of pointers to build the genjet
103  constituents.push_back( refToPtr( *igr) );
104  }
105 
106  math::XYZPoint vertex;
107  GenJet::Specific specific;
108 
109  GenJet jet( sumVisMom, vertex, specific, constituents);
110 
111  if (charge != (*iTau)->charge() )
112  std::cout<<" charge of Tau: " << (*iTau) << " not equal to charge of sum of charge of all descendents. " << std::cout;
113 
114  jet.setCharge(charge);
115  pOutVisTaus->push_back( jet );
116 
117  }
118  iEvent.put( pOutVisTaus );
119 }
120 
121 
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
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
virtual void setCharge(Charge q)
set electric charge
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
std::vector< GenJet > GenJetCollection
collection of GenJet objects
#define abs(x)
Definition: mlp_lapack.h:159
std::vector< Constituent > Constituents
Definition: Jet.h:24
int Charge
electric charge type
Definition: Particle.h:24
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:243
double charge(const std::vector< uint8_t > &Ampls)
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:238
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:30
int iEvent
Definition: GenABIO.cc:243
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:84
Jets made from MC generator particles.
Definition: GenJet.h:25
virtual void beginJob()
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:355
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:13
virtual void produce(edm::Event &, const edm::EventSetup &)
size_type size() const
Size of the RefVector.
Definition: RefVector.h:84
TauGenJetProducer(const edm::ParameterSet &)
tuple cout
Definition: gather_cfg.py:41
void findDescendents(const reco::GenParticleRef &base, reco::GenParticleRefVector &descendents, int status, int pdgId=0)
find all descendents of a given status and pdgId (recursive)