TEvtGen/Tauola/TauolaHepMCEvent.cxx

   1 #include "TauolaHepMCEvent.h"
   2 #include "TauolaLog.h"
   3
   4 using namespace std;
   5
   6 namespace Tauolapp
   7 {
   8
   9 TauolaHepMCEvent::TauolaHepMCEvent(HepMC::GenEvent * event){
  10   m_event=event;
  11
  12   // Default units
  13   m_momentum_unit = "GEV";
  14   m_length_unit   = "MM";
  15
  16   if(m_event->momentum_unit() != HepMC::Units::GEV) m_momentum_unit = "MEV";
  17   if(m_event->length_unit()   != HepMC::Units::MM ) m_length_unit   = "CM";
  18
  19   // If needed - change units used by HepMC to GEV and MM
  20   if( m_event->momentum_unit() != HepMC::Units::GEV ||
  21       m_event->length_unit()   != HepMC::Units::MM     )
  22   {
  23     m_event->use_units(HepMC::Units::GEV,HepMC::Units::MM);
  24   }
  25 }
  26
  27 TauolaHepMCEvent::~TauolaHepMCEvent(){
  28
  29   while(m_tau_list.size()!=0){
  30     TauolaParticle * temp = m_tau_list.back();
  31     m_tau_list.pop_back();
  32     delete temp;
  33   }
  34
  35 }
  36
  37 HepMC::GenEvent * TauolaHepMCEvent::getEvent(){
  38   return m_event;
  39 }
  40
  41 std::vector<TauolaParticle*> TauolaHepMCEvent::findParticles(int pdg_id){
  42
  43   if(m_tau_list.size()==0){
  44
  45     HepMC::GenEvent::particle_const_iterator part_itr = m_event->particles_begin();
  46     //loop over all particle in the event looking for taus (or other)
  47     for( ; part_itr!=m_event->particles_end(); part_itr++){
  48       if(abs((*part_itr)->pdg_id())==pdg_id)
  49         m_tau_list.push_back(new TauolaHepMCParticle(*part_itr));
  50     }
  51   }
  52   return m_tau_list;
  53 }
  54
  55 std::vector<TauolaParticle*> TauolaHepMCEvent::findStableParticles(int pdg_id){
  56
  57   /**  HepMC::GenEvent::particle_const_iterator part_itr = m_event->particles_begin();
  58   //loop over all particle in the event looking for taus (or other)
  59   for( ; part_itr!=m_event->particles_end(); part_itr++){
  60     if(fabs((*part_itr)->pdg_id())==pdg_id){
  61       if((*part_itr)->end_vertex()){
  62         cout << "WARNING: Particle with pdg code " << (*part_itr)->pdg_id()
  63              << " has end vertex" <<endl;
  64       }
  65       else
  66         list.push_back(new TauolaHepMCParticle(*part_itr));
  67     }
  68     }**/
  69
  70   std::vector<TauolaParticle*> tau_list = findParticles(pdg_id);
  71   std::vector<TauolaParticle*> stable_tau_list;
  72
  73   for(int i=0; i<(int) tau_list.size(); i++){
  74
  75     if(!tau_list.at(i)->hasDaughters())
  76       stable_tau_list.push_back(tau_list.at(i));
  77     else
  78     {
  79       std::vector<TauolaParticle*> t = tau_list.at(i)->getDaughters();
  80       //Ignore taus that we won't be decaying anyway
  81       if(t.size()==1) continue;
  82       if(t.size()==2 && (abs(t[0]->getPdgID())==15 || abs(t[1]->getPdgID())==15) ) continue;
  83       Log::Warning()<<"Particle with pdg code "<<tau_list.at(i)->getPdgID()
  84                     <<" already has daughters" <<endl;
  85     }
  86   }
  87
  88   return stable_tau_list;
  89
  90 }
  91
  92 void TauolaHepMCEvent::eventEndgame(){
  93
  94   //Set output units for the event
  95   string momentum("GEV"),length("MM");
  96
  97   switch(Tauola::momentumUnit)
  98   {
  99     case Tauola::GEV:
 100       momentum = "GEV";
 101       break;
 102     case Tauola::MEV:
 103       momentum = "MEV";
 104       break;
 105     default:
 106       momentum = m_momentum_unit;
 107   }
 108
 109   switch(Tauola::lengthUnit)
 110   {
 111     case Tauola::MM:
 112       length = "MM";
 113       break;
 114     case Tauola::CM:
 115       length = "CM";
 116       break;
 117     default:
 118       length = m_length_unit;
 119   }
 120
 121   m_event->use_units(momentum,length);
 122 }
 123
 124 } // namespace Tauolapp