Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Protected Member Functions | Protected Attributes
HepMCG4Interface Class Reference

#include <HepMCG4Interface.hh>

Inheritance diagram for HepMCG4Interface:
G4VPrimaryGenerator G4VPrimaryGenerator HepMCG4AsciiReader HepMCG4AsciiReader HepMCG4PythiaInterface HepMCG4PythiaInterface

Public Member Functions

 HepMCG4Interface ()
 
virtual ~HepMCG4Interface ()
 
HepMC::GenEvent * GetHepMCGenEvent () const
 
virtual void GeneratePrimaryVertex (G4Event *anEvent)
 
 HepMCG4Interface ()
 
virtual ~HepMCG4Interface ()
 
HepMC::GenEvent * GetHepMCGenEvent () const
 
virtual void GeneratePrimaryVertex (G4Event *anEvent)
 
- Public Member Functions inherited from G4VPrimaryGenerator
 G4VPrimaryGenerator ()
 
virtual ~G4VPrimaryGenerator ()
 
G4ThreeVector GetParticlePosition ()
 
G4double GetParticleTime ()
 
void SetParticlePosition (G4ThreeVector aPosition)
 
void SetParticleTime (G4double aTime)
 

Protected Member Functions

virtual G4bool CheckVertexInsideWorld (const G4ThreeVector &pos) const
 
void HepMC2G4 (const HepMC::GenEvent *hepmcevt, G4Event *g4event)
 
virtual HepMC::GenEvent * GenerateHepMCEvent ()
 
virtual G4bool CheckVertexInsideWorld (const G4ThreeVector &pos) const
 
void HepMC2G4 (const HepMC::GenEvent *hepmcevt, G4Event *g4event)
 
virtual HepMC::GenEvent * GenerateHepMCEvent ()
 

Protected Attributes

HepMC::GenEvent * hepmcEvent
 
- Protected Attributes inherited from G4VPrimaryGenerator
G4ThreeVector particle_position
 
G4double particle_time
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VPrimaryGenerator
static G4bool CheckVertexInsideWorld (const G4ThreeVector &pos)
 

Detailed Description

A base class for primary generation via HepMC object. This class is derived from G4VPrimaryGenerator.

Definition at line 41 of file HepMCEx01/include/HepMCG4Interface.hh.

Constructor & Destructor Documentation

HepMCG4Interface::HepMCG4Interface ( )

Definition at line 44 of file HepMCEx01/src/HepMCG4Interface.cc.

45  : hepmcEvent(0)
46 {
47 }
HepMCG4Interface::~HepMCG4Interface ( )
virtual

Definition at line 50 of file HepMCEx01/src/HepMCG4Interface.cc.

References hepmcEvent.

51 {
52  delete hepmcEvent;
53 }
HepMCG4Interface::HepMCG4Interface ( )
virtual HepMCG4Interface::~HepMCG4Interface ( )
virtual

Member Function Documentation

G4bool HepMCG4Interface::CheckVertexInsideWorld ( const G4ThreeVector pos) const
protectedvirtual

Definition at line 57 of file HepMCEx01/src/HepMCG4Interface.cc.

References G4TransportationManager::GetTransportationManager(), kInside, and write_gdml::navigator.

Referenced by HepMC2G4().

58 {
60  -> GetNavigatorForTracking();
61 
62  G4VPhysicalVolume* world= navigator-> GetWorldVolume();
63  G4VSolid* solid= world-> GetLogicalVolume()-> GetSolid();
64  EInside qinside= solid-> Inside(pos);
65 
66  if( qinside != kInside) return false;
67  else return true;
68 }
static G4TransportationManager * GetTransportationManager()
EInside
Definition: geomdefs.hh:58
virtual G4bool HepMCG4Interface::CheckVertexInsideWorld ( const G4ThreeVector pos) const
protectedvirtual
virtual HepMC::GenEvent* HepMCG4Interface::GenerateHepMCEvent ( )
protectedvirtual
HepMC::GenEvent * HepMCG4Interface::GenerateHepMCEvent ( )
protectedvirtual

Reimplemented in HepMCG4PythiaInterface, HepMCG4PythiaInterface, HepMCG4AsciiReader, and HepMCG4AsciiReader.

Definition at line 119 of file HepMCEx01/src/HepMCG4Interface.cc.

Referenced by GeneratePrimaryVertex().

120 {
121  HepMC::GenEvent* aevent= new HepMC::GenEvent();
122  return aevent;
123 }
void HepMCG4Interface::GeneratePrimaryVertex ( G4Event anEvent)
virtual

Implements G4VPrimaryGenerator.

Definition at line 126 of file HepMCEx01/src/HepMCG4Interface.cc.

References G4cout, G4endl, GenerateHepMCEvent(), G4RunManager::GetRunManager(), HepMC2G4(), and hepmcEvent.

127 {
128  // delete previous event object
129  delete hepmcEvent;
130 
131  // generate next event
133  if(! hepmcEvent) {
134  G4cout << "HepMCInterface: no generated particles. run terminated..."
135  << G4endl;
136  G4RunManager::GetRunManager()-> AbortRun();
137  return;
138  }
139  HepMC2G4(hepmcEvent, anEvent);
140 }
virtual HepMC::GenEvent * GenerateHepMCEvent()
G4GLOB_DLL std::ostream G4cout
void HepMC2G4(const HepMC::GenEvent *hepmcevt, G4Event *g4event)
static G4RunManager * GetRunManager()
Definition: G4RunManager.cc:74
#define G4endl
Definition: G4ios.hh:61
virtual void HepMCG4Interface::GeneratePrimaryVertex ( G4Event anEvent)
virtual

Implements G4VPrimaryGenerator.

HepMC::GenEvent* HepMCG4Interface::GetHepMCGenEvent ( ) const
HepMC::GenEvent * HepMCG4Interface::GetHepMCGenEvent ( ) const
inline

Definition at line 74 of file HepMCEx01/include/HepMCG4Interface.hh.

References hepmcEvent.

75 {
76  return hepmcEvent;
77 }
void HepMCG4Interface::HepMC2G4 ( const HepMC::GenEvent *  hepmcevt,
G4Event g4event 
)
protected

Definition at line 71 of file HepMCEx01/src/HepMCG4Interface.cc.

References python.hepunit::c_light, CheckVertexInsideWorld(), python.hepunit::GeV, python.hepunit::mm, and position.

Referenced by GeneratePrimaryVertex().

73 {
74  for(HepMC::GenEvent::vertex_const_iterator vitr= hepmcevt->vertices_begin();
75  vitr != hepmcevt->vertices_end(); ++vitr ) { // loop for vertex ...
76 
77  // real vertex?
78  G4bool qvtx=false;
79  for (HepMC::GenVertex::particle_iterator
80  pitr= (*vitr)->particles_begin(HepMC::children);
81  pitr != (*vitr)->particles_end(HepMC::children); ++pitr) {
82 
83  if (!(*pitr)->end_vertex() && (*pitr)->status()==1) {
84  qvtx=true;
85  break;
86  }
87  }
88  if (!qvtx) continue;
89 
90  // check world boundary
91  HepMC::FourVector pos= (*vitr)-> position();
92  G4LorentzVector xvtx(pos.x(), pos.y(), pos.z(), pos.t());
93  if (! CheckVertexInsideWorld(xvtx.vect()*mm)) continue;
94 
95  // create G4PrimaryVertex and associated G4PrimaryParticles
96  G4PrimaryVertex* g4vtx=
97  new G4PrimaryVertex(xvtx.x()*mm, xvtx.y()*mm, xvtx.z()*mm,
98  xvtx.t()*mm/c_light);
99 
100  for (HepMC::GenVertex::particle_iterator
101  vpitr= (*vitr)->particles_begin(HepMC::children);
102  vpitr != (*vitr)->particles_end(HepMC::children); ++vpitr) {
103 
104  if( (*vpitr)->status() != 1 ) continue;
105 
106  G4int pdgcode= (*vpitr)-> pdg_id();
107  pos= (*vpitr)-> momentum();
108  G4LorentzVector p(pos.px(), pos.py(), pos.pz(), pos.e());
109  G4PrimaryParticle* g4prim=
110  new G4PrimaryParticle(pdgcode, p.x()*GeV, p.y()*GeV, p.z()*GeV);
111 
112  g4vtx-> SetPrimary(g4prim);
113  }
114  g4event-> AddPrimaryVertex(g4vtx);
115  }
116 }
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
bool G4bool
Definition: G4Types.hh:79
virtual G4bool CheckVertexInsideWorld(const G4ThreeVector &pos) const
int position
Definition: filter.cc:7
float c_light
Definition: hepunit.py:257
void HepMCG4Interface::HepMC2G4 ( const HepMC::GenEvent *  hepmcevt,
G4Event g4event 
)
protected

Field Documentation

HepMC::GenEvent * HepMCG4Interface::hepmcEvent
protected

The documentation for this class was generated from the following files: