6 modified files
slic/include
diff -u -r1.44 -r1.45
--- LcioManager.hh 31 Aug 2005 19:38:56 -0000 1.44
+++ LcioManager.hh 14 Sep 2005 23:10:37 -0000 1.45
@@ -103,10 +103,6 @@
// formerly writeEvent()
void endEvent(const G4Event* anEvent);
- // set the current G4Event
- void setCurrentG4Event(const G4Event* anEvent) { m_currentG4Event = anEvent; }
- const G4Event* getCurrentG4Event() { return m_currentG4Event; }
-
// set the current LCEvent
void setCurrentLCEvent(IMPL::LCEventImpl* anLCEvent) { m_currentLCEvent = anLCEvent; }
IMPL::LCEventImpl* getCurrentLCEvent() { return m_currentLCEvent; }
@@ -211,9 +207,6 @@
// current LCEvent
IMPL::LCEventImpl* m_currentLCEvent;
- // current G4Event
- const G4Event* m_currentG4Event;
-
// messenger
LcioMessenger* m_messenger;
slic/include
diff -u -r1.36 -r1.37
--- LcioMcpManager.hh 14 Sep 2005 18:03:06 -0000 1.36
+++ LcioMcpManager.hh 14 Sep 2005 23:10:37 -0000 1.37
@@ -22,13 +22,12 @@
{
class LcioMcpMessenger;
+ class LcioMcpFactory;
class Trajectory;
/**
* @class LcioMcpManager
- * @brief This singleton manager class is the workhorse for post-event
- * creation of the Lcio MCParticle collection from stored G4Event data.
- * @note Uses input Mcps, trajectories and primaries to create the output Mcp coll.
+ * @brief This singleton manages creation of the LCIO MCParticle collection from G4Event data.
*/
class LcioMcpManager : public Verbose
{
@@ -58,43 +57,50 @@
// action at end of event
void endEvent(const G4Event* anEvent);
- // set the current G4Event
- void setCurrentG4Event(const G4Event* anEvent) { m_currentG4Event = anEvent; }
- const G4Event* getCurrentG4Event() { return m_currentG4Event; }
-
- // set the current LCEvent
- void setCurrentLCEvent(IMPL::LCEventImpl* anLCEvent) { m_currentLCEvent = anLCEvent; }
- IMPL::LCEventImpl* getCurrentLCEvent() { return m_currentLCEvent; }
-
// connect an initial Mcp to a primary
void addMcpToPrimaryLink(EVENT::MCParticle* mcp, G4PrimaryParticle* primary);
// connect track to MCParticle via track ID
void addTrackIDToMcpLink(G4int trkID, IMPL::MCParticleImpl* mcp);
+ // connect initial Mcp to final Mcp
+ void addInitialMcpToFinalMcpLink( EVENT::MCParticle* mcpInit, IMPL::MCParticleImpl* mcpFinal);
+
// find an MCParticle from the trackID
IMPL::MCParticleImpl* findMcpFromTrackID( G4int trkID );
// find matching primary particle given an MCP ptr
G4PrimaryParticle* findPrimaryFromMcp(EVENT::MCParticle* mcp);
- // get the TrackToMcpMap
- inline const TrackToMcpMap& getTrackToMcpMap() const;
+ // loop over input collection to find matching Mcp primary for this trajectory
+ EVENT::MCParticle* findPrimaryInitialMcpFromTrajectory(Trajectory* trj);
+
+ // find a dau Mcp given a G4PrimaryParticle dau
+ EVENT::MCParticle* findDaughterMcpFromPrimary(EVENT::MCParticle* mcpInit,
+ G4PrimaryParticle* primary);
- // get the McpToPrimaryMap
- inline const McpToPrimaryMap& getMcpToPrimaryMap() const;
+ // given an input particle from stdhep, find the proper output particle
+ IMPL::MCParticleImpl* findFinalParticleFromInitial( EVENT::MCParticle* mcpInit);
- // get MCP collection generated from StdHep
+ /*
+ * Create an empty MCParticle collection indexed by collName.
+ * SIDE EFFECT: Automatically registers collection into the map.
+ */
+ EVENT::LCCollection* createMcpCollection(const std::string& collName);
+
+ /**
+ * Get the MCP collection generated from StdHep.
+ * Creates empty collection if doesn't exist.
+ */
inline EVENT::LCCollection* getInitialMcpCollection();
- // get MCP collection generated at EndOfEvent
+ /**
+ * Get the MCP collection generated at EndOfEvent.
+ * Creates empty collection if doesn't exist.
+ */
inline EVENT::LCCollection* getFinalMcpCollection();
- // get Mcp coll map
- LcioMcpManager::McpCollectionMap* getMcpCollectionMap();
-
void setMinimumTrackingDistance(double minDist);
-
double getMinimumTrackingDistance();
// print maps
@@ -103,97 +109,33 @@
void printMcpToPrimaryMap();
void printInitMcpToMcpMap();
- // reset at end of event
+ // reset the manager at the end of event
void reset();
// verbosity
virtual inline void setVerbose(Verbose::LevelType lt);
- // distance computation based on Mcp 3vect positions
- static G4double computeMcpDistance(IMPL::MCParticleImpl*, IMPL::MCParticleImpl*);
- static G4double computeDistance(const double *v1, const double *v2);
-
// set initial mcp coll externally
void setInitialMcpCollection(EVENT::LCCollection* mcpColl);
- // create new, blank Mcp coll
- EVENT::LCCollection* createEmptyMcpCollection(const std::string& collName, bool errorOnExist = true);
+ /* Register an LCCollection object containing MCParticles, by name. */
+ void registerMcpCollection(const std::string& collName, EVENT::LCCollection*);
- // find the Mcp coll
+ /* Find the first MCParticle coll in the LCevent. */
EVENT::LCCollection* findMcpCollection(EVENT::LCEvent* event);
+ /* Find a registered MCParticle coll. */
+ EVENT::LCCollection* findMcpCollection(const std::string& collName);
+
// get the Mcp printer
LcioMcpPrinter* getMcpPrinter();
private:
- // convert to CLHEP vector and return vertex.isNear(parEndpoint)
- G4bool isNear(const double* vertex, const double* parEndpoint);
-
- // find a dau Mcp given a G4PrimaryParticle dau
- EVENT::MCParticle* findDaughterMcpFromPrimary(EVENT::MCParticle* mcpInit,
- G4PrimaryParticle* primary);
-
- // given an input particle from stdhep, find the proper output particle
- IMPL::MCParticleImpl* findFinalParticleFromInitial( EVENT::MCParticle* mcpInit);
-
- // loop over input collection to find matching Mcp primary for this trajectory
- EVENT::MCParticle* findPrimaryInitialMcpFromTrajectory(Trajectory* trj);
-
- // Mcp with no parents
- inline bool isPrimary(IMPL::MCParticleImpl* mcp);
-
- // chooses to create from initial or traj cont
- void createFinalMcpCollection();
-
- // use StdHep-generated mcpVec to create final mcpVec for an event
- void createFinalMcpCollectionFromInitial(EVENT::LCCollection* mcpVecInitial);
-
- // These two functions create the Mcp coll from trajectories only in case of G4 GPS or gun
- void createFinalMcpCollectionFromTrajectoryContainer(G4TrajectoryContainer* trjCont);
-
- /* Create Mcp daughters of a trajectory with the given track ID. */
- void addMcpDaughtersFromTrajectoryContainer(IMPL::MCParticleImpl* parMcp, G4int parTrkID);
-
- // create an Mcp from a trajectory only
- IMPL::MCParticleImpl* createMcpFromTrajectory( Trajectory* );
-
- // create Mcps recursively given an initial Mcp
- IMPL::MCParticleImpl* createMcpFromInitialRecurse(EVENT::MCParticle* mcp);
-
- // create an Mcp from another, associated initial Mcp only
- IMPL::MCParticleImpl* createMcpFromInitialOnly(EVENT::MCParticle* mcpInit);
-
- // create an Mcp from a G4PrimaryParticle
- IMPL::MCParticleImpl* createMcpFromPrimary(G4PrimaryParticle*, EVENT::MCParticle* );
-
- // create an Mcp from primary and input Mcp
- IMPL::MCParticleImpl* createMcpFromInitialAndPrimary(G4PrimaryParticle* primary, EVENT::MCParticle* mcpInit);
-
- // create an Mcp from a trajectory and an input Mcp
- IMPL::MCParticleImpl* createMcpFromInitialAndTrajectory( Trajectory*, EVENT::MCParticle* mcpInit );
-
- // create an Mcp from a G4PrimaryParticle, shallow copy only
- IMPL::MCParticleImpl* createMcpFromPrimaryShallowCopy(G4PrimaryParticle*);
-
- // create new Mcp and shallow copy from input Mcp
- IMPL::MCParticleImpl* createMcpShallowCopy(EVENT::MCParticle* mcp);
-
- // create and add MCP daughters based on G4PrimaryParticle and associated MCP from StdHep
- void createDaughtersFromPrimary(G4PrimaryParticle* primary,
- EVENT::MCParticle* mcpInit,
- IMPL::MCParticleImpl* mcpPar);
-
- // connect initial Mcp to final Mcp
- void addInitialMcpToFinalMcpLink( EVENT::MCParticle* mcpInit, IMPL::MCParticleImpl* mcpFinal);
-
- // add daughters to a MCP based on associated intial MCP from StdHep
- void addMcpDaughtersFromInitial( IMPL::MCParticleImpl* mcpNew, EVENT::MCParticle* mcpInit );
-
// clear the maps
void clearMaps();
- // delete the initial collection as it isn't attached to an LCEvent
+ // delete the initial collection
void deleteInitialMcpCollection();
private:
@@ -204,9 +146,6 @@
// current LCEvent
IMPL::LCEventImpl* m_currentLCEvent;
- // current G4Event
- const G4Event* m_currentG4Event;
-
// trackID -> MCP
TrackToMcpMap m_trackToMcp;
@@ -228,6 +167,9 @@
// map of Mcp colls being tracked by the manager
McpCollectionMap m_mcpColls;
+
+ /* Factory for creating MCParticles from Geant4 runtime data. */
+ LcioMcpFactory* m_factory;
};
};
slic/src
diff -u -r1.59 -r1.60
--- LcioManager.cc 7 Sep 2005 22:23:33 -0000 1.59
+++ LcioManager.cc 14 Sep 2005 23:10:37 -0000 1.60
@@ -1,4 +1,4 @@
-// $Header: /cvs/lcd/slic/src/LcioManager.cc,v 1.59 2005/09/07 22:23:33 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/LcioManager.cc,v 1.60 2005/09/14 23:10:37 jeremy Exp $
#include "LcioManager.hh"
// SLIC
@@ -28,11 +28,12 @@
#include "IMPL/LCCollectionVec.h"
#include "UTIL/LCTOOLS.h"
-// G4
+// geant4
#include "G4Run.hh"
#include "G4VHitsCollection.hh"
#include "G4SDManager.hh"
#include "G4TrajectoryContainer.hh"
+#include "G4EventManager.hh"
//
#include <ctime>
@@ -339,14 +340,12 @@
LCEventImpl* LcioManager::createLCEvent()
{
- return createLCEvent( getCurrentG4Event() );
+ return createLCEvent( G4EventManager::GetEventManager()
+ ->GetNonconstCurrentEvent() );
}
- void LcioManager::endEvent(const G4Event* anEvent)
+ void LcioManager::endEvent(const G4Event*)
{
- // set current G4Event
- setCurrentG4Event( anEvent );
-
// create LCEvent
createLCEvent();
@@ -383,13 +382,11 @@
void LcioManager::addHCs()
{
- m_HCBuilder->createHCsFromG4Event( m_currentG4Event, m_currentLCEvent );
+ m_HCBuilder->createHCsFromG4Event( G4EventManager::GetEventManager()->GetNonconstCurrentEvent(), m_currentLCEvent );
}
void LcioManager::reset()
{
- m_currentG4Event = 0;
-
// delete transient event container
delete m_currentLCEvent;
m_currentLCEvent = 0;
@@ -401,7 +398,7 @@
void LcioManager::createFinalMcpCollection()
{
// create Mcp coll in LcioMcpManager
- m_mcpManager->endEvent( getCurrentG4Event() );
+ m_mcpManager->endEvent( G4EventManager::GetEventManager()->GetNonconstCurrentEvent() );
// add Mcp coll to current event
getCurrentLCEvent()->addCollection( m_mcpManager->getFinalMcpCollection(), "MCParticle" );
slic/src
diff -u -r1.49 -r1.50
--- LcioMcpManager.cc 14 Sep 2005 18:03:06 -0000 1.49
+++ LcioMcpManager.cc 14 Sep 2005 23:10:37 -0000 1.50
@@ -1,10 +1,10 @@
-// $Header: /cvs/lcd/slic/src/LcioMcpManager.cc,v 1.49 2005/09/14 18:03:06 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/LcioMcpManager.cc,v 1.50 2005/09/14 23:10:37 jeremy Exp $
#include "LcioMcpManager.hh"
// SLIC
#include "LcioMcpMessenger.hh"
#include "LcioManager.hh"
-#include "LcioMcpStatusSetter.hh"
+#include "LcioMcpFactory.hh"
#include "TrajectoryManager.hh"
#include "PrimaryGeneratorAction.hh"
#include "EventSourceManager.hh"
@@ -51,15 +51,18 @@
// setup mcp printer
m_mcpPrinter = new LcioMcpPrinter();
+
+ /* Setup the MCParticle factory. */
+ m_factory = new LcioMcpFactory(this);
}
- void LcioMcpManager::endEvent(const G4Event* anEvent)
+ void LcioMcpManager::endEvent(const G4Event*)
{
- // set current G4Event
- setCurrentG4Event( anEvent );
-
- // create the Mcp collection
- createFinalMcpCollection();
+ /*
+ * Create the Mcp collection. Factory will
+ * retrieve current event from G4EventManager.
+ */
+ m_factory->createFinalMcpCollection();
#ifdef SLIC_VERBOSE
if ( verbose() > 2 ) {
@@ -100,151 +103,57 @@
return mcpColl;
}
- void LcioMcpManager::createFinalMcpCollection()
+ EVENT::LCCollection* LcioMcpManager::findMcpCollection(const std::string& collName)
{
- // new coll for final Mcps
- createEmptyMcpCollection("final");
-
- // using StdHep or LCIO e.g. a file-based source?
- if ( EventSourceManager::instance()->isFileSource() ) {
-
- // create from the initial coll of Mcp
- createFinalMcpCollectionFromInitial( getInitialMcpCollection() );
- }
- // G4 source using traj only
- else {
-
- // create from trj cont only
- createFinalMcpCollectionFromTrajectoryContainer( getCurrentG4Event()->GetTrajectoryContainer() );
+ EVENT::LCCollection* coll = 0;
+ if ( m_mcpColls[collName] ) {
+ coll = m_mcpColls[collName];
}
+ return coll;
}
- void LcioMcpManager::createFinalMcpCollectionFromInitial(EVENT::LCCollection* mcpVecInitial)
+ void LcioMcpManager::setInitialMcpCollection(EVENT::LCCollection* mcpColl)
{
- // coll must exist
- if ( mcpVecInitial ) {
-
- // get num Mcp in coll
- int numInitMcp = mcpVecInitial->getNumberOfElements();
-
- // coll must have members
- if ( numInitMcp > 0 ) {
-
- // loop over input Mcp coll
- for ( int i=0; i < numInitMcp; i++ ) {
-
- // get Mcp
- MCParticleImpl* mcp = static_cast<MCParticleImpl*> ( mcpVecInitial->getElementAt( i ) );
-
- // create only from Mcp primaries, e.g. no parents
- if ( isPrimary( mcp ) ) {
- createMcpFromInitialRecurse( mcp );
- }
- }
- }
- else {
- G4Exception("Initial McpVec has no members.");
- }
+ if ( mcpColl ) {
+ m_mcpColls["initial"] = mcpColl;
}
else {
- G4Exception("Initial McpVec ptr is null.");
+ G4Exception("FATAL ERROR: Initial collection cannot be set to null.");
}
}
- EVENT::LCCollection* LcioMcpManager::createEmptyMcpCollection(const std::string& collName,
- bool errorOnExist)
+ EVENT::LCCollection* LcioMcpManager::createMcpCollection(const std::string& collName)
{
#ifdef SLIC_VERBOSE
if ( verbose() ) {
- std::cout << "creating empty mcp coll <" << collName << ">" << std::endl;
+ std::cout << "creating mcp coll <" << collName << ">" << std::endl;
}
#endif
- EVENT::LCCollection* mcpVec = 0;
-
- // found existing
- if ( m_mcpColls[collName] ) {
- if ( errorOnExist ) {
- G4Exception("Attempting to create duplicate Mcp collection.");
- }
- else {
- mcpVec = m_mcpColls[collName];
- }
+ LCCollection* coll = 0;
+ if ( findMcpCollection(collName) ) {
+ G4Exception("Attempting to create duplicate Mcp collection.");
}
- // new mcp coll
else {
- mcpVec = new LCCollectionVec( LCIO::MCPARTICLE );
- m_mcpColls[collName] = mcpVec;
+ coll = new LCCollectionVec(LCIO::MCPARTICLE);
+ registerMcpCollection( collName, coll );
}
- return mcpVec;
+ return coll;
}
- void LcioMcpManager::setInitialMcpCollection(EVENT::LCCollection* mcpColl)
+ void LcioMcpManager::registerMcpCollection(const std::string& collName, EVENT::LCCollection* mcpColl)
{
- if ( mcpColl ) {
- m_mcpColls["initial"] = mcpColl;
+ if (mcpColl->getTypeName() != LCIO::MCPARTICLE ) {
+ G4Exception("ERROR: The LCCollection does not contain MCParticle objects.");
}
- else {
- G4Exception("FATAL ERROR: Initial collection cannot be set to null.");
- }
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromTrajectory( Trajectory* trj)
- {
- // new Mcp
- MCParticleImpl* mcp = new MCParticleImpl();
-
- // PDG
- mcp->setPDG( trj->GetPDGEncoding() );
-
- // mom
- G4ThreeVector mom = trj->GetInitialMomentum() / GeV;
- float fmom[3] = { mom.x(), mom.y(), mom.z() };
- mcp->setMomentum( fmom );
-
- // mass
- mcp->setMass( trj->GetParticleDefinition()->GetPDGMass() / GeV );
-
- // charge
- mcp->setCharge( trj->GetCharge() );
-
- // begin and end points
- G4int npts = trj->GetPointEntries();
-
- // begin pnt
- G4VTrajectoryPoint* beginTrjPnt = trj->GetPoint( 0 );
- G4ThreeVector beginPos = beginTrjPnt->GetPosition() / mm;
- double fBeginPos[3] = { beginPos.x(), beginPos.y(), beginPos.z() };
- mcp->setVertex( fBeginPos );
-
- // end pnt
- G4VTrajectoryPoint* endTrjPnt = trj->GetPoint( npts - 1 );
- G4ThreeVector endPos = endTrjPnt->GetPosition() / mm;
- double fEndPos[3] = { endPos.x(), endPos.y(), endPos.z() };
- mcp->setEndpoint( fEndPos );
-
- // add to map of trkID -> mcp
- addTrackIDToMcpLink( trj->GetTrackID(), mcp );
-
- // set status codes from trajectory info (with no dep on any associated Mcps)
- LcioMcpStatusSetter::setMcpStatusCodesFromTrajectory( trj, mcp );
- return mcp;
- }
-
- G4bool LcioMcpManager::isNear(const double* vertex, const double* parEndpoint)
- {
- G4ThreeVector vtx( vertex[0], vertex[1], vertex[2] );
- G4ThreeVector endpnt( parEndpoint[0], parEndpoint[1], parEndpoint[2] );
-
- G4bool r = true;
-
- if ( !vtx.isNear(endpnt) ) {
- r = false;
+ if ( m_mcpColls[collName] ) {
+ std::cerr << "WARNING: Mcp collection " + collName + " is already registered" << std::endl;
+ }
+ else {
+ m_mcpColls[collName] = mcpColl;
}
-
- return r;
}
void LcioMcpManager::addTrackIDToMcpLink(G4int trkID,
@@ -270,192 +179,6 @@
return mcp;
}
- G4double LcioMcpManager::computeMcpDistance(IMPL::MCParticleImpl* mcp1,
- IMPL::MCParticleImpl* mcp2 )
- {
- return computeDistance( mcp1->getVertex(), mcp2->getVertex() );
- }
-
- G4double LcioMcpManager::computeDistance(const double *v1, const double *v2)
- {
- // distance = sqrt ( x^2 + y^2 + z^2 )
-
- // x^2
- G4double x = (v1[0] - v2[0]);
- G4double x2 = pow( x, 2 );
-
- // y^2
- G4double y = (v1[1] - v2[1]);
- G4double y2 = pow( y, 2 );
-
- // z^2
- G4double z = (v1[2] - v2[2]);
- G4double z2 = pow( z, 2 );
-
- return sqrt( x2 + y2 + z2 );
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromInitialRecurse(EVENT::MCParticle* mcpInit)
- {
- MCParticleImpl* mcp = 0;
-
- // find associated primary
- G4PrimaryParticle* g4primary = findPrimaryFromMcp( mcpInit );
-
- // no associated primary with this initial mcp?
- if( !g4primary ) {
-
- // create from initial Mcp only
- mcp = createMcpFromInitialOnly( mcpInit );
- }
- // create from primary
- else {
-
- // create it from primary
- mcp = createMcpFromPrimary( g4primary, mcpInit );
- }
-
- // return the output MCP
- return mcp;
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromInitialOnly(EVENT::MCParticle* mcpInit)
- {
- // copy shallow the initial
- MCParticleImpl* mcp = createMcpShallowCopy( mcpInit );
-
- // should add to Mcp map here?
- addInitialMcpToFinalMcpLink( mcpInit, mcp);
-
- // given initial, create and add daughters
- addMcpDaughtersFromInitial( mcp, mcpInit );
-
- // add to final coll
- getFinalMcpCollection()->addElement( mcp );
-
- return mcp;
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromPrimary(G4PrimaryParticle* primary,
- EVENT::MCParticle* mcpInit)
- {
- G4int trkID = primary->GetTrackID();
- MCParticleImpl* mcp = 0;
- Trajectory* trj = TrajectoryManager::instance()->findTrajectory( trkID );
-
- /*
- * No trajectory?
- *
- * This is a trap for G4Primaries with unrecognized PDGid,
- * aka the dreaded "unknown particle". In this case, no
- * trajectory was created, so the final Mcp must be composed
- * from trajectory and initial Mcp information.
- *
- */
- if ( !trj ) {
-
- /* Recursively create Mcps from primary and initial Mcp. */
- mcp = createMcpFromInitialAndPrimary( primary, mcpInit );
-
- }
- /* Found a trajectory. */
- else {
-
- /* Recursively create Mcps from trajectory and initial Mcp. */
- mcp = createMcpFromInitialAndTrajectory( trj, mcpInit );
-
- }
-
- /*
- * In either case, the generator status is set from the initial particle.
- */
- mcp->setGeneratorStatus( mcpInit->getGeneratorStatus() );
-
- /*
- * The particle did not originate in the simulator, because it
- * has a matching primary, which means that it had a corresponding
- * record from the event generator.
- */
- mcp->setCreatedInSimulation( false );
-
- /*
- * Link the input MCP to an output MCP in order to check later if it was
- * already created.
- */
- addInitialMcpToFinalMcpLink( mcpInit, mcp);
-
- return mcp;
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromInitialAndPrimary(G4PrimaryParticle* primary,
- EVENT::MCParticle* mcpInit)
- {
- /*
- * Create new Mcp from primary particle,
- * shallow copy only (no recursion).
- */
- MCParticleImpl* mcp = createMcpFromPrimaryShallowCopy( primary );
-
- /* Get vertex from initial particle, because primary does not have it -- was immediately decayed */
- double vtx[3] = { mcpInit->getVertex()[0], mcpInit->getVertex()[1], mcpInit->getVertex()[2] };
- mcp->setVertex( vtx );
-
- /* setup primary's daughters */
- createDaughtersFromPrimary( primary, /* this primary */
- mcpInit, /* associated initial Mcp particle */
- mcp /* parent Mcp particle */
- );
-
-#ifdef SLIC_VERBOSE
- if ( mcp->getDaughters().size() == 0 ) {
- if ( verbose() ) {
- std::cout << "No Mcp daughters added." << std::endl;
- }
- }
-#endif
-
- // done with particle; add to coll
- getFinalMcpCollection()->addElement( mcp );
-
- return mcp;
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromInitialAndTrajectory(Trajectory* trj,
- EVENT::MCParticle* mcpInit)
- {
- // shallow copy create with basic info and no daughters
- IMPL::MCParticleImpl* mcp = createMcpFromTrajectory( trj );
-
- // set generator status and created in sim flag
- LcioMcpStatusSetter::setGeneratorStatus( mcpInit, mcp );
-
- // loop over all trajectories to find trajectory daughter
- G4TrajectoryContainer* trjCont = m_currentG4Event->GetTrajectoryContainer();
- G4int numTrj = trjCont->entries();
- for( int j=0; j < numTrj; j++ ) {
-
- Trajectory* thisTrj = static_cast<Trajectory*> ((*trjCont)[j]);
-
- // found daughter
- if ( thisTrj->GetParentID() == trj->GetTrackID() ) {
-
- // find matching primary Mcp input particle
- MCParticle* mcpPrim = findPrimaryInitialMcpFromTrajectory( thisTrj );
-
- // create MCP output particle from this input daughter
- MCParticleImpl* mcpDauNew = createMcpFromInitialAndTrajectory( thisTrj, mcpPrim );
-
- // set parent link of new daughter
- mcpDauNew->addParent( mcp );
- }
- }
-
- // add to collection
- getFinalMcpCollection()->addElement( mcp );
-
- return mcp;
- }
-
EVENT::MCParticle* LcioMcpManager::findPrimaryInitialMcpFromTrajectory(Trajectory* trj)
{
EVENT::MCParticle* mcpPrim = 0;
@@ -488,124 +211,6 @@
return mcpPrim;
}
- // Mcp with no parents
- inline bool LcioMcpManager::isPrimary(IMPL::MCParticleImpl* mcp)
- {
- return ( mcp->getParents().size() == 0 );
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpFromPrimaryShallowCopy(G4PrimaryParticle* primary)
- {
- // new MCP
- MCParticleImpl* mcp = new MCParticleImpl();
-
- // PDG
- mcp->setPDG( primary->GetPDGcode() );
-
- // momentum
- G4ThreeVector pVec = primary->GetMomentum() / GeV;
- float p[3] = { pVec.x(), pVec.y(), pVec.z() };
- mcp->setMomentum( p );
-
- // mass
- mcp->setMass( primary->GetMass() / GeV );
-
- // charge
- mcp->setCharge( primary->GetCharge() );
-
- return mcp;
- }
-
-
- void LcioMcpManager::addMcpDaughtersFromInitial(IMPL::MCParticleImpl* mcpNew,
- EVENT::MCParticle* mcpInit)
- {
- int numDau = mcpInit->getDaughters().size();
-
- for ( int i=0; i < numDau; i++ ) {
-
- // child of input particle
- MCParticle* mcpChildInit = mcpInit->getDaughters()[i];
-
- // use a search function to avoid insert when null!
- MCParticleImpl* mcpChildFinal = findFinalParticleFromInitial( mcpChildInit );
-
- // if not added already
- if ( mcpChildFinal == 0 ) {
-
- // create daughter
- mcpChildFinal = createMcpFromInitialRecurse( mcpChildInit );
- }
-
- // add parent regardless of above; maybe new MCP or already created
- mcpChildFinal->addParent( mcpNew );
- }
- }
-
- IMPL::MCParticleImpl* LcioMcpManager::createMcpShallowCopy(EVENT::MCParticle* mcp)
- {
-
- IMPL::MCParticleImpl* mcpNew = new IMPL::MCParticleImpl();
-
- mcpNew->setPDG( mcp->getPDG() );
-
- float p[3] = { mcp->getMomentum()[0], mcp->getMomentum()[1], mcp->getMomentum()[2] };
- mcpNew->setMomentum( p );
-
- mcpNew->setMass( mcp->getMass() );
-
- double vtx[3] = { mcp->getVertex()[0], mcp->getVertex()[1], mcp->getVertex()[2] };
- mcpNew->setVertex( vtx );
-
- mcpNew->setGeneratorStatus( mcp->getGeneratorStatus() );
-
- mcpNew->setSimulatorStatus( 0 );
-
- mcpNew->setCharge( mcp->getCharge() );
-
- return mcpNew;
- }
-
- // create LCIO MCParticles from daughters, setting proper linkage
- void LcioMcpManager::createDaughtersFromPrimary(G4PrimaryParticle* primary,
- EVENT::MCParticle* mcpInit,
- IMPL::MCParticleImpl* mcpPar)
- {
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 1 ) {
- std::cout << "LcioMcpManager::createDaughtersFromPrimary()" << std::endl;
- std::cout << "primary <" << primary << ">" << std::endl;
- std::cout << "trkID <" << primary->GetTrackID() << ">" << std::endl;
- std::cout << "mcpInit <" << mcpInit << ">" << std::endl;
- std::cout << "mcpPar <" << mcpPar << ">" << std::endl;
- }
-#endif
-
- // loop over primary daughters
- G4PrimaryParticle* primDau = primary->GetDaughter();
- while ( primDau ) {
-
- // find input Mcp for the primary dau
- MCParticle* mcpDau = findDaughterMcpFromPrimary( mcpInit,
- primDau );
-
- // must find it or die
- if ( 0 == mcpDau ) {
- G4Exception( "Mcp daughter was not found." );
- }
-
- // create dau Mcp
- MCParticleImpl* mcpDauNew = createMcpFromPrimary( primDau, mcpDau );
-
- // make mcpPar the parent
- mcpDauNew->addParent( mcpPar );
-
- // get next primary dau
- primDau = primDau->GetNext();
- }
- }
-
EVENT::MCParticle* LcioMcpManager::findDaughterMcpFromPrimary(EVENT::MCParticle* mcpInit,
G4PrimaryParticle* primDau)
{
@@ -627,75 +232,6 @@
return mcpDau;
}
- void LcioMcpManager::createFinalMcpCollectionFromTrajectoryContainer(G4TrajectoryContainer* trjCont)
- {
- if ( trjCont ) {
-
- // loop over trj cont
- G4int n_trj = trjCont->entries();
- for ( G4int i = 0;
- i < n_trj;
- i++ ) {
-
- Trajectory* trj = static_cast<Trajectory*> ((*trjCont)[i]);
-
- // only add primary trajectories; daughters are picked up within
- if ( trj->GetParentID() == 0 ) {
-
- // create mcp output particle from trj
- MCParticleImpl* mcp = createMcpFromTrajectory( trj );
-
- // add this to coll
- getFinalMcpCollection()->addElement( mcp );
-
- // add daughters
- addMcpDaughtersFromTrajectoryContainer( mcp, trj->GetTrackID() );
- }
- }
- }
- else {
- G4Exception("G4TrajectoryContainer is null!");
- }
- }
-
- void LcioMcpManager::addMcpDaughtersFromTrajectoryContainer( MCParticleImpl* parMcp,
- G4int parTrkID)
- {
- // loop over trj cont
- G4TrajectoryContainer* trjCont = m_currentG4Event->GetTrajectoryContainer();
- G4int n_trj = trjCont->entries();
- for ( G4int i = 0;
- i < n_trj;
- i++ ) {
-
- Trajectory* trj = static_cast<Trajectory*> ((*trjCont)[i]);
-
- if ( trj->GetParentID() == parTrkID ) {
-
- // lkp dau Mcp by track ID
- MCParticleImpl* dauMcp = findMcpFromTrackID( trj->GetTrackID() );
-
- if ( dauMcp == 0 ) {
-
- // create from trj if null
- dauMcp = createMcpFromTrajectory( trj );
- }
-
- // must be created by now
- assert( dauMcp );
-
- // add this to final coll
- getFinalMcpCollection()->addElement( dauMcp );
-
- // add daughters
- addMcpDaughtersFromTrajectoryContainer( dauMcp, trj->GetTrackID() );
-
- // set parent in daughter
- dauMcp->addParent( parMcp );
- }
- }
- }
-
void LcioMcpManager::deleteInitialMcpCollection()
{
LCCollection* mcpVec = getInitialMcpCollection();
@@ -705,13 +241,11 @@
}
}
- // reset at end of event
+ /* Does not need to delete Mcp final collection, as it is deleted by the LCEvent's dtor. */
void LcioMcpManager::reset()
{
deleteInitialMcpCollection();
- // Do not need to delete Mcp final collection, as it is deleted by the LCEvent's dtor.
-
clearMaps();
}
@@ -754,32 +288,23 @@
return pp;
}
- // get the TrackToMcpMap
- inline const LcioMcpManager::TrackToMcpMap& LcioMcpManager::getTrackToMcpMap() const
- {
- return m_trackToMcp;
- }
-
- // get the McpToPrimaryMap
- inline const LcioMcpManager::McpToPrimaryMap& LcioMcpManager::getMcpToPrimaryMap() const
- {
- return m_McpToPrimary;
- }
-
// get MCP collection generated at EndOfEvent
inline EVENT::LCCollection* LcioMcpManager::getFinalMcpCollection()
{
+ if ( !findMcpCollection("final") ) {
+ createMcpCollection("final");
+ }
+
return m_mcpColls["final"];
}
inline EVENT::LCCollection* LcioMcpManager::getInitialMcpCollection()
{
- return m_mcpColls["initial"];
- }
+ if ( !findMcpCollection("initial") ) {
+ createMcpCollection("initial");
+ }
- LcioMcpManager::McpCollectionMap* LcioMcpManager::getMcpCollectionMap()
- {
- return &m_mcpColls;
+ return m_mcpColls["initial"];
}
void LcioMcpManager::setMinimumTrackingDistance(double minDist)
slic/src
diff -u -r1.24 -r1.25
--- LcioPrimaryGenerator.cc 17 May 2005 19:06:34 -0000 1.24
+++ LcioPrimaryGenerator.cc 14 Sep 2005 23:10:38 -0000 1.25
@@ -1,10 +1,13 @@
-// $Header: /cvs/lcd/slic/src/LcioPrimaryGenerator.cc,v 1.24 2005/05/17 19:06:34 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/LcioPrimaryGenerator.cc,v 1.25 2005/09/14 23:10:38 jeremy Exp $
+// slic
#include "LcioPrimaryGenerator.hh"
#include "LcioManager.hh"
#include "LcioMcpManager.hh"
+#include "LcioMcpUtil.hh"
#include "StringUtil.hh"
+// std
#include <sstream>
using IMPL::MCParticleImpl;
@@ -17,20 +20,13 @@
LcioPrimaryGenerator::LcioPrimaryGenerator(LcioManager* mgr)
{
m_mgr = mgr;
-
+
m_mcpManager = LcioMcpManager::instance();
}
LcioPrimaryGenerator::~LcioPrimaryGenerator()
{}
- // from G4
- //void LcioPrimaryGenerator::GeneratePrimaryVertex(G4Event* anEvent)
- //{
- // call to mgr creates MCP Coll within mgr from current StdHep event record (for usage later)
- // generatePrimaryVertexFromMcpCollection( LcioMcpManager::instance()->createInitialMcpCollFromStdHep(), anEvent );
- //}
-
void LcioPrimaryGenerator::generatePrimaryVertexFromMcpCollection( LCCollection* mcpVec, G4Event* anEvent )
{
assert( mcpVec );
@@ -55,7 +51,7 @@
std::cout << "i <" << i << ">" << std::endl;
}
#endif
-
+
MCParticleImpl* mcp = dynamic_cast<MCParticleImpl*>
(mcpVec->getElementAt( i ) );
@@ -67,7 +63,7 @@
}
#endif
- G4int numPar = mcp->getParents().size();
+ G4int numPar = mcp->getParents().size();
// handle 4 cases, following RC in LCS's LCSHEPEvtInterface
// case 1: predecay
@@ -77,7 +73,7 @@
// loop vars
MCParticleImpl* parMcp = 0;
- G4PrimaryParticle* g4parent = 0;
+ G4PrimaryParticle* g4parent = 0;
G4PrimaryParticle* thePrimary = 0;
G4PrimaryVertex* theVertex = 0;
@@ -90,7 +86,7 @@
if ( numPar > 0 ) {
parMcp = dynamic_cast<MCParticleImpl*>(mcp->getParents()[0]);
g4parent = m_mcpManager->findPrimaryFromMcp( parMcp );
-
+
if ( g4parent != 0 ) {
#ifdef SLIC_VERBOSE
@@ -111,13 +107,13 @@
#ifdef SLIC_VERBOSE
if ( verbose() > 2 ) {
- std::cout << "case 1: predecay" << std::endl;
+ std::cout << "case 1: predecay" << std::endl;
}
#endif
createPrimary = true;
// no vertex
- }
+ }
else {
// case 2: final state
if ( mcp->getGeneratorStatus() == 1 ) {
@@ -131,16 +127,16 @@
createPrimary = true;
createVertex = true;
}
- else {
+ else {
G4double dist = 0;
-
- // check case 3, if traveled > minDist
+
+ // check case 3, if traveled > minDist
if ( mcp->getDaughters().size() > 0 ) {
IMPL::MCParticleImpl* firstDau = dynamic_cast<IMPL::MCParticleImpl*>
( mcp->getDaughters()[0] );
- dist = LcioMcpManager::computeMcpDistance( mcp, firstDau );
+ dist = LcioMcpUtil::computeMcpDistance( mcp, firstDau );
// case 3
if ( dist > m_mcpManager->getMinimumTrackingDistance() ) {
@@ -178,16 +174,16 @@
#endif
thePrimary = createPrimaryParticleFromMcp( mcp );
-
+
// set daughters for PreDecay
- if ( isPreDecay ) {
+ if ( isPreDecay ) {
#ifdef SLIC_VERBOSE
if ( verbose() > 2 ) {
std::cout << "predecay setup" << std::endl;
}
#endif
-
+
assert( g4parent );
// computation of proper_time from RC
@@ -202,13 +198,13 @@
std::cout << "proper_time <" << proper_time << ">" << std::endl;
}
#endif
-
+
g4parent->SetDaughter( thePrimary );
g4parent->SetProperTime( proper_time );
#ifdef SLIC_VERBOSE
if ( verbose() > 2 ) {
- std::cout << "mcp decay time <" << mcp->getTime() - parMcp->getTime() << ">" << std::endl;
+ std::cout << "mcp decay time <" << mcp->getTime() - parMcp->getTime() << ">" << std::endl;
}
#endif
}
@@ -222,7 +218,7 @@
std::cout << "creating vertex" << std::endl;
}
#endif
-
+
theVertex = createPrimaryVertexFromMcp( mcp );
theVertex->SetPrimary( thePrimary );
anEvent->AddPrimaryVertex( theVertex );
@@ -257,7 +253,7 @@
mcp->getMomentum()[2] * GeV
);
primary->SetMass( mcp->getMass() * GeV );
-
+
return primary;
}
slic/src
diff -u -r1.2 -r1.3
--- StdHepEventSource.cc 26 May 2005 19:34:08 -0000 1.2
+++ StdHepEventSource.cc 14 Sep 2005 23:10:38 -0000 1.3
@@ -1,4 +1,4 @@
-// $Header: /cvs/lcd/slic/src/StdHepEventSource.cc,v 1.2 2005/05/26 19:34:08 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/StdHepEventSource.cc,v 1.3 2005/09/14 23:10:38 jeremy Exp $
#include "StdHepEventSource.hh"
namespace slic
@@ -81,11 +81,8 @@
// not EOF?
if ( !isEOF() ) {
- // create Mcp's initial coll
- EVENT::LCCollection* initMcpColl = LcioMcpManager::instance()->createEmptyMcpCollection("initial");
-
// set convertor's coll
- m_convertor->setCurrentMcpCollection( initMcpColl );
+ m_convertor->setCurrentMcpCollection( LcioMcpManager::instance()->getInitialMcpCollection() );
// run convertor
m_convertor->fillCurrentMcpCollectionFromStdHep();
CVSspam 0.2.8
