G4HEPEvtInterface Class Reference

#include <G4HEPEvtInterface.hh>

Inheritance diagram for G4HEPEvtInterface:

G4VPrimaryGenerator

Public Member Functions

 G4HEPEvtInterface (char *evfile)
 G4HEPEvtInterface (G4String evfile)
 ~G4HEPEvtInterface ()
void GeneratePrimaryVertex (G4Event *evt)

Detailed Description

Definition at line 76 of file G4HEPEvtInterface.hh.


Constructor & Destructor Documentation

G4HEPEvtInterface::G4HEPEvtInterface ( char *  evfile  ) 

Definition at line 43 of file G4HEPEvtInterface.cc.

References FatalException, G4Exception(), G4VPrimaryGenerator::particle_position, and G4VPrimaryGenerator::particle_time.

00044 {
00045   inputFile.open(evfile);
00046   if (inputFile) {
00047     fileName = evfile;
00048   }
00049   else {
00050     G4Exception("G4HEPEvtInterface::G4HEPEvtInterface","Event0201",FatalException,
00051     "G4HEPEvtInterface:: cannot open file.");
00052   }
00053   G4ThreeVector zero;
00054   particle_position = zero;
00055   particle_time = 0.0;
00056 
00057 }

G4HEPEvtInterface::G4HEPEvtInterface ( G4String  evfile  ) 

Definition at line 59 of file G4HEPEvtInterface.cc.

References G4String::data(), FatalException, G4Exception(), G4VPrimaryGenerator::particle_position, and G4VPrimaryGenerator::particle_time.

00060 {
00061   const char* fn = evfile.data();
00062   inputFile.open((char*)fn);
00063   if (inputFile) {
00064     fileName = evfile;
00065   }
00066   else {
00067     G4Exception("G4HEPEvtInterface::G4HEPEvtInterface","Event0201",FatalException,
00068     "G4HEPEvtInterface:: cannot open file.");
00069   }
00070   G4ThreeVector zero;
00071   particle_position = zero;
00072   particle_time = 0.0;
00073 }

G4HEPEvtInterface::~G4HEPEvtInterface (  ) 

Definition at line 75 of file G4HEPEvtInterface.cc.

00076 {;}


Member Function Documentation

void G4HEPEvtInterface::GeneratePrimaryVertex ( G4Event evt  )  [virtual]

Implements G4VPrimaryGenerator.

Definition at line 78 of file G4HEPEvtInterface.cc.

References G4Event::AddPrimaryVertex(), G4Exception(), JustWarning, G4VPrimaryGenerator::particle_position, G4VPrimaryGenerator::particle_time, G4PrimaryParticle::SetDaughter(), G4PrimaryParticle::SetMass(), G4PrimaryParticle::SetMomentum(), and G4PrimaryVertex::SetPrimary().

00079 {
00080   G4int NHEP;  // number of entries
00081   inputFile >> NHEP;
00082   if( inputFile.eof() ) 
00083   {
00084     G4Exception("G4HEPEvtInterface::GeneratePrimaryVertex","Event0202",
00085     JustWarning,"End-Of-File : HEPEvt input file");
00086     return;
00087   }
00088 
00089   for( G4int IHEP=0; IHEP<NHEP; IHEP++ )
00090   {
00091     G4int ISTHEP;   // status code
00092     G4int IDHEP;    // PDG code
00093     G4int JDAHEP1;  // first daughter
00094     G4int JDAHEP2;  // last daughter
00095     G4double PHEP1; // px in GeV
00096     G4double PHEP2; // py in GeV
00097     G4double PHEP3; // pz in GeV
00098     G4double PHEP5; // mass in GeV
00099 
00100     inputFile >> ISTHEP >> IDHEP >> JDAHEP1 >> JDAHEP2
00101        >> PHEP1 >> PHEP2 >> PHEP3 >> PHEP5;
00102 
00103     // create G4PrimaryParticle object
00104     G4PrimaryParticle* particle 
00105       = new G4PrimaryParticle( IDHEP );
00106     particle->SetMass( PHEP5*GeV );
00107     particle->SetMomentum(PHEP1*GeV, PHEP2*GeV, PHEP3*GeV );
00108 
00109     // create G4HEPEvtParticle object
00110     G4HEPEvtParticle* hepParticle
00111       = new G4HEPEvtParticle( particle, ISTHEP, JDAHEP1, JDAHEP2 );
00112 
00113     // Store
00114     HPlist.push_back( hepParticle );
00115   }
00116 
00117   // check if there is at least one particle
00118   if( HPlist.size() == 0 ) return; 
00119 
00120   // make connection between daughter particles decayed from 
00121   // the same mother
00122   for( size_t i=0; i<HPlist.size(); i++ )
00123   {
00124     if( HPlist[i]->GetJDAHEP1() > 0 ) //  it has daughters
00125     {
00126       G4int jda1 = HPlist[i]->GetJDAHEP1()-1; // FORTRAN index starts from 1
00127       G4int jda2 = HPlist[i]->GetJDAHEP2()-1; // but C++ starts from 0.
00128       G4PrimaryParticle* mother = HPlist[i]->GetTheParticle();
00129       for( G4int j=jda1; j<=jda2; j++ )
00130       {
00131         G4PrimaryParticle* daughter = HPlist[j]->GetTheParticle();
00132         if(HPlist[j]->GetISTHEP()>0)
00133         {
00134           mother->SetDaughter( daughter );
00135           HPlist[j]->Done();
00136         }
00137       }
00138     }
00139   }
00140 
00141   // create G4PrimaryVertex object
00142   G4PrimaryVertex* vertex = new G4PrimaryVertex(particle_position,particle_time);
00143 
00144   // put initial particles to the vertex
00145   for( size_t ii=0; ii<HPlist.size(); ii++ )
00146   {
00147     if( HPlist[ii]->GetISTHEP() > 0 ) // ISTHEP of daughters had been 
00148                                        // set to negative
00149     {
00150       G4PrimaryParticle* initialParticle = HPlist[ii]->GetTheParticle();
00151       vertex->SetPrimary( initialParticle );
00152     }
00153   }
00154 
00155   // clear G4HEPEvtParticles
00156   //HPlist.clearAndDestroy();
00157   for(size_t iii=0;iii<HPlist.size();iii++)
00158   { delete HPlist[iii]; }
00159   HPlist.clear();
00160 
00161   // Put the vertex to G4Event object
00162   evt->AddPrimaryVertex( vertex );
00163 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:12 2013 for Geant4 by  doxygen 1.4.7