AthenaHepMCtoG4EventAction.cpp

Go to the documentation of this file.

//------------------------------------------------------------------
//
// ClassName:   AthenaHepMCtoG4EventAction
//  
// Description: PrimaryGeneratorAction to convert HepMC event to G4Event
//
// Author:      Charles Leggett
// 
// Date:        3-8-2001
// 
// $Id: AthenaHepMCtoG4EventAction.cpp,v 1.1.1.1 2013/11/22 02:11:07 juxd Exp $
// Simulation/G4Sim/G4Svc tag $Name: G4Svc-00-04-01 $
//
//------------------------------------------------------------------
 
#include "G4Svc/AthenaHepMCtoG4EventAction.h"
 
#include "HepMC/GenEvent.h"
 
#include "G4Event.hh"
#include "G4UImanager.hh"
#include "globals.hh"
 
#include "G4PrimaryParticle.hh"
#include "G4PrimaryVertex.hh"
#include "G4LorentzVector.hh"
#include "G4Geantino.hh"
 
#include "GaudiKernel/IMessageSvc.h"
 
using namespace HepMC;
 
AthenaHepMCtoG4EventAction::AthenaHepMCtoG4EventAction( int level):
  m_logLevel(level) {
}
 
AthenaHepMCtoG4EventAction::AthenaHepMCtoG4EventAction():
m_logLevel( int(MSG::INFO) ) {
}
 
AthenaHepMCtoG4EventAction::~AthenaHepMCtoG4EventAction() {
 
}
 
//* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
 
//void AthenaHepMCtoG4EventAction::SetHepMCEvent(const HepMC::GenEvent *p) { 
void AthenaHepMCtoG4EventAction::SetHepMCEvent(McEventCol::iterator p) { 
  
  if ( m_logLevel <= int(MSG::DEBUG) ) {
    std::cout << "AthenaHepMCtoG4EventAction: - Setting HepMC Event ptr" 
          << std::endl;
  }
 
  //  p_evt = p;
  p_evtCollItr = p;
 
}
 
//* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
 
//  void AthenaHepMCtoG4EventAction::GeneratePrimaries(G4Event* anEvent)
//  {
//    if ( m_logLevel <= int(MSG::DEBUG) ) {
//      std::cout << "AthenaHepMCtoG4EventAction: - Generating Primaries" 
//            << std::endl;
//    }
 
//    p_evt = (*p_evtCollItr)->pGenEvt();
 
//    if (p_evt == 0) { 
 
//      if ( m_logLevel <= int(MSG::WARNING) ) {
//        std::cout << "AthenaHepMCtoG4EventAction: - HepMC ptr == null" 
//          << std::endl;
//      }
 
//      return; 
 
//    }
 
//    int nv = 1;
//    int mnv;
 
//    std::cout << "how many vertices: ";
//    std::cin >> mnv;
 
//    for ( HepMC::GenEvent::vertex_const_iterator v = p_evt->vertices_begin();
//          v != p_evt->vertices_end(); ++v ) {
//      bool IsGood=false;
//      for (HepMC::GenVertex::particle_iterator it=
//         (*v)->particles_begin(HepMC::children);
//       it!=(*v)->particles_end(HepMC::children);
//       it++) {
      
//        if (!(*it)->end_vertex() && (*it)->status()==1 ) { 
//      IsGood=true;
//      break;
//        }
//      }
 
//      if (IsGood) {
//        std::cout << "found a good vertex, barcode: " << (*v)->barcode() 
//          << std::endl;
      
//        G4LorentzVector lv=(*v)->position();
//        G4PrimaryVertex *vert= new G4PrimaryVertex(lv.x(),lv.y(),lv.z(),lv.t());
 
//        for (HepMC::GenVertex::particle_iterator it=
//           (*v)->particles_begin(HepMC::children);
//         it!=(*v)->particles_end(HepMC::children);
//         it++) {
    
//      //  if ( (*it)->status() != 1) continue;
 
//      int pdgcode=(*it)->pdg_id();
//      GenVertex *dcy = (*it)->end_vertex();
//      std::cout << "found a good particle, barcode: " << (*it)->barcode()
//            << "  pdg: " << pdgcode << "  status: " << (*it)->status()
//            << "  end: ";
//      if (dcy != 0) {
//        std::cout << dcy->barcode();
//      } else {
//        std::cout << "nil";
//      }
//      std::cout << std::endl;
 
//      const G4ThreeVector& p=(*it)->momentum();
//      G4PrimaryParticle *part;
//      if (pdgcode!=999) {
//        part=new G4PrimaryParticle (pdgcode,p.x()*GeV,p.y()*GeV,p.z()*GeV);
//      } else {
//        G4ParticleDefinition *pdef = G4Geantino::GeantinoDefinition();
//        part=new G4PrimaryParticle (pdef,p.x()*GeV,p.y()*GeV,p.z()*GeV);
//      }
//      vert->SetPrimary(part);
//        }
//        anEvent->AddPrimaryVertex(vert);
//        if (nv > mnv) { break; }
//        nv ++;
//      }
//    }       
 
//    std::cout << "done generating event" << std::endl;
 
//    p_evtCollItr++;
  
//  }
 
//* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
 
void AthenaHepMCtoG4EventAction::GeneratePrimaries(G4Event* anEvent)
{
  if ( m_logLevel <= int(MSG::DEBUG) ) {
    std::cout << "AthenaHepMCtoG4EventAction: - Generating Primaries" 
          << std::endl;
  }
 
  
  p_evt = (const HepMC::GenEvent *)(*p_evtCollItr);
 
  
  if (p_evt == 0) { 
 
    if ( m_logLevel <= int(MSG::WARNING) ) {
      std::cout << "AthenaHepMCtoG4EventAction: - HepMC ptr == null" 
        << std::endl;
    }
 
    return; 
 
  }
  
  G4PrimaryVertex* vertex = new G4PrimaryVertex(0.,0.,0.,0.);
  
  int n=0;
  for ( HepMC::GenEvent::vertex_const_iterator v = p_evt->vertices_begin();
        v != p_evt->vertices_end(); ++v ) {
    bool IsGood=false;
    for (HepMC::GenVertex::particle_iterator it=
       (*v)->particles_begin(HepMC::children);
     it!=(*v)->particles_end(HepMC::children);
     it++) {
      
      if (!(*it)->end_vertex() && (*it)->status()==1 ) { 
    IsGood=true;
    break;
      }
    }
 
    if (IsGood) {
//        std::cout << "found a good vertex, barcode: " << (*v)->barcode() 
//          << std::endl;
      
      for (HepMC::GenVertex::particle_iterator it=
         (*v)->particles_begin(HepMC::children);
       it!=(*v)->particles_end(HepMC::children);
       it++) {
    
        if ( (*it)->status() != 1) continue;
 
    int pdgcode=(*it)->pdg_id();
//      GenVertex *dcy = (*it)->end_vertex();
//      std::cout << "found a good particle, barcode: " << (*it)->barcode()
//            << "  pdg: " << pdgcode << "  status: " << (*it)->status()
//            << "  end: ";
//      if (dcy != 0) {
//        std::cout << dcy->barcode();
//      } else {
//        std::cout << "nil";
//      }
//      std::cout << std::endl;
 
    const G4ThreeVector& p=(*it)->momentum();
    G4PrimaryParticle *part;
    if (pdgcode!=999) {
      part=new G4PrimaryParticle (pdgcode,p.x()*GeV,p.y()*GeV,p.z()*GeV);
    } else {
      G4ParticleDefinition *pdef = G4Geantino::GeantinoDefinition();
      part=new G4PrimaryParticle (pdef,p.x()*GeV,p.y()*GeV,p.z()*GeV);
    }
 
      vertex->SetPosition((*v)->position().x(),
                          (*v)->position().y(),
                          (*v)->position().z());
      vertex->SetT0((*v)->position().t());
 
    vertex->SetPrimary(part);
 
    n++;
      }
    }
  }       
  anEvent->AddPrimaryVertex(vertex);
 
  if ( m_logLevel <= int(MSG::DEBUG) ) {
    std::cout << "AthenaHepMCtoG4EventAction: - done. G4Event has " 
          << n << "  particles. "
          << std::endl;
  }
 
  p_evtCollItr++;
  
}