UnstableFinalState.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 #include "Rivet/Projections/UnstableFinalState.hh"
00003 #include "Rivet/Tools/ParticleIdUtils.hh"
00004 #include "Rivet/Tools/Logging.hh"
00005 #include "Rivet/Cmp.hh"
00006 
00007 #define IS_PARTON_PDGID(id) ( abs(id) <= 100 && abs(id) != 22 && (abs(id) < 11 || abs(id) > 18) )
00008 
00009 namespace Rivet {
00010 
00011 
00012   int UnstableFinalState::compare(const Projection& p) const {
00013     const UnstableFinalState& other = dynamic_cast<const UnstableFinalState&>(p);
00014     return \
00015       cmp(_etamin, other._etamin) ||
00016       cmp(_etamax, other._etamax) ||
00017       cmp(_ptmin, other._ptmin);
00018   }
00019 
00020 
00021   void UnstableFinalState::project(const Event& e) {
00022     _theParticles.clear();
00023 
00024     foreach (GenParticle* p, Rivet::particles(e.genEvent())) {
00025       const int st = p->status();
00026       bool passed = \
00027         ( st == 1 || (st == 2 && abs(p->pdg_id()) != 22) ) &&
00028         !isZero(p->momentum().perp()) && p->momentum().perp() >= _ptmin &&
00029         p->momentum().eta() > _etamin && p->momentum().eta() < _etamax &&
00030         !IS_PARTON_PDGID(p->pdg_id());
00031 
00032       // Avoid double counting by re-marking as unpassed if particle ID == parent ID
00033       const GenVertex* pv = p->production_vertex();
00034       const GenVertex* dv = p->end_vertex();
00035       if (passed && pv) {
00036         for (GenVertex::particles_in_const_iterator pp = pv->particles_in_const_begin() ;
00037              pp != pv->particles_in_const_end() ; ++pp) {
00038           if ( p->pdg_id() == (*pp)->pdg_id() ) {
00039             passed = false;
00040             break;
00041           }
00042         }
00043       }
00044 
00045       // Add to output particles collection
00046       if (passed) {
00047         _theParticles.push_back(Particle(*p));
00048       }
00049 
00050       // Log parents and children
00051       if (getLog().isActive(Log::TRACE)) {
00052         MSG_TRACE("ID = " << p->pdg_id()
00053                   << ", status = " << st
00054                   << ", pT = " << p->momentum().perp()
00055                   << ", eta = " << p->momentum().eta()
00056                   << ": result = " << std::boolalpha << passed);
00057         if (pv) {
00058           for (GenVertex::particles_in_const_iterator pp = pv->particles_in_const_begin() ;
00059                pp != pv->particles_in_const_end() ; ++pp) {
00060             MSG_TRACE("  parent ID = " << (*pp)->pdg_id());
00061           }
00062         }
00063         if (dv) {
00064           for (GenVertex::particles_out_const_iterator pp = dv->particles_out_const_begin() ;
00065                pp != dv->particles_out_const_end() ; ++pp) {
00066             MSG_TRACE("  child ID  = " << (*pp)->pdg_id());
00067           }
00068         }
00069       }
00070     }
00071     MSG_DEBUG("Number of final-state particles = " << _theParticles.size());
00072   }
00073 
00074 
00075 }