| Commit in slic/src on MAIN | |||
| LcioMcpFactory.cc | +38 | -14 | 1.14 -> 1.15 |
| LcioMcpMaps.cc | +217 | -204 | 1.5 -> 1.6 |
| StdHepToLcioConvertor.cc | +15 | -12 | 1.37 -> 1.38 |
| TrajectoryManager.cc | +169 | -230 | 1.49 -> 1.50 |
| +439 | -460 | ||
minor cleanup; modifications to debug output
diff -u -r1.14 -r1.15 --- LcioMcpFactory.cc 15 Oct 2009 19:14:39 -0000 1.14 +++ LcioMcpFactory.cc 16 Dec 2009 00:16:04 -0000 1.15 @@ -1,4 +1,4 @@
-// $Header: /cvs/lcd/slic/src/LcioMcpFactory.cc,v 1.14 2009/10/15 19:14:39 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/LcioMcpFactory.cc,v 1.15 2009/12/16 00:16:04 jeremy Exp $
#include "LcioMcpFactory.hh" // slic
@@ -28,7 +28,6 @@
namespace slic
{
-
LcioMcpFactory::LcioMcpFactory(LcioMcpManager* manager)
: Module( "LcioMcpFactory" ),
m_finalColl(0),
@@ -303,10 +302,15 @@
EVENT::MCParticle* mcpInit)
{
#ifdef SLIC_LOG
- log() << LOG::done << "createMcpFromPrimary: " << m_manager->getMCParticleIndex( m_manager->getInitialMcpCollection(), mcpInit) << LOG::done;
+ log() << LOG::done << "createMcpFromPrimary - MCP idx <" << m_manager->getMCParticleIndex( m_manager->getInitialMcpCollection(), mcpInit) << ">" << LOG::done;
#endif
int trkID = primary->GetTrackID();
+ +#ifdef SLIC_LOG + log() << LOG::done << "primary TID <" << trkID << ">" << LOG::done; +#endif +
MCParticleImpl* mcp = 0;
Trajectory* trj = TrajectoryManager::instance()->findTrajectory( trkID );
@@ -343,6 +347,10 @@
m_manager->getMaps()->addInitialMcpToFinalMcpLink( mcpInit, mcp);
+#ifdef SLIC_LOG + log() << LOG::debug << "createMcpFromPrimary() - adding trackId <" << trj->GetTrackID() << "> to MCP <" << mcp << "> link." << LOG::done; +#endif +
// Link the track ID to the output particle.
m_manager->getMaps()->addTrackIDToMcpLink( trkID, mcp );
@@ -371,7 +379,8 @@
); #ifdef SLIC_DEBUG
- if ( mcp->getDaughters().size() == 0 ) {
+ if ( mcp->getDaughters().size() == 0 )
+ {
#ifdef SLIC_LOG
log() << LOG::debug << "No Mcp daughters added." << LOG::done;
#endif
@@ -387,7 +396,9 @@
IMPL::MCParticleImpl* LcioMcpFactory::createMcpFromInitialAndTrajectory(Trajectory* trj,
EVENT::MCParticle* mcpInit)
{
- log() << LOG::done << "createMcpFromInitialAndTrajectory: " << m_manager->getMCParticleIndex( m_manager->getInitialMcpCollection(), mcpInit) << LOG::done;
+#ifdef SLIC_LOG + log() << LOG::done << "createMcpFromInitialAndTrajectory: " << m_manager->getMCParticleIndex( m_manager->getInitialMcpCollection(), mcpInit ) << LOG::done; +#endif
// shallow copy create with basic info and no daughters
IMPL::MCParticleImpl* mcp = createMcpFromTrajectory( trj );
@@ -398,16 +409,29 @@
// loop over all trajectories to find trajectory daughter
int numTrj = m_currentTrajectoryContainer->entries();
+#ifdef SLIC_LOG
log() << LOG::done << "nTrajectoryDau: " << numTrj << LOG::done;
+#endif
- for( int j=0; j < numTrj; j++ ) {
- Trajectory* thisTrj = static_cast<Trajectory*> ((*m_currentTrajectoryContainer)[j]);
+//#ifdef SLIC_LOG
+// if ( thisTrj->GetParentID() == -1 )
+// {
+// log() << LOG::always << "The trajectory with trackID <" << trj->GetTrackID() << "> has a parent with trackID of -1 !" << LOG::done;
+// }
+//#endif
- // found daughter
- if ( thisTrj->GetParentID() == trj->GetTrackID() ) {
+ for( int j=0; j < numTrj; j++ )
+ {
-// log() << LOG::done << "trajDau: " << j << LOG::done;
+ Trajectory* thisTrj = static_cast<Trajectory*> ( ( *m_currentTrajectoryContainer ) [j] );
+
+ // found daughter
+ if ( thisTrj->GetParentID() == trj->GetTrackID() )
+ {
+#ifdef SLIC_LOG
+ log() << LOG::debug << "Trajectory with trackID <" << trj->GetTrackID() << "> has parent with trackID <" << thisTrj->GetParentID() << ">" << LOG::done;
+#endif
// find matching primary Mcp input particle
MCParticle* mcpPrim = m_manager->getMaps()->findPrimaryInitialMcpFromTrajectory( thisTrj );
@@ -419,10 +443,10 @@
mcpDauNew->addParent( mcp );
}
// DEBUG
- else
- {
- log() << LOG::done << "BAD NEWS! Could not find Trajectory! <" << j << "> with TID <" << trj->GetTrackID() << ">" << LOG::done;
- }
+ //else
+ //{
+ // log() << LOG::done << "BAD NEWS! Could not find Trajectory <" << j << "> with TID <" << trj->GetTrackID() << ">" << LOG::done;
+ //}
// END DEBUG
}
diff -u -r1.5 -r1.6 --- LcioMcpMaps.cc 27 Apr 2007 01:54:35 -0000 1.5 +++ LcioMcpMaps.cc 16 Dec 2009 00:16:04 -0000 1.6 @@ -1,4 +1,4 @@
-// $Header: /cvs/lcd/slic/src/LcioMcpMaps.cc,v 1.5 2007/04/27 01:54:35 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/LcioMcpMaps.cc,v 1.6 2009/12/16 00:16:04 jeremy Exp $
#include "LcioMcpMaps.hh" // slic
@@ -12,208 +12,221 @@
namespace slic
{
- LcioMcpMaps::LcioMcpMaps(LcioMcpManager* manager)
- : Module( "LcioMcpMaps" )
- {
- m_manager = manager;
- }
-
- LcioMcpMaps::~LcioMcpMaps()
- {}
-
- void LcioMcpMaps::addTrackIDToMcpLink(G4int trkID,
- MCParticleImpl* mcp)
- {
- assert( mcp );
- assert(trkID > 0);
-
- m_trackToMcp[ trkID ] = mcp;
- }
-
- MCParticleImpl* LcioMcpMaps::findMcpFromTrackID( G4int trkID ) const
- {
- MCParticleImpl* mcp = 0;
- for ( TrackToMcpMap::const_iterator iter = m_trackToMcp.begin();
- iter != m_trackToMcp.end();
- iter++ ) {
- if ( iter->first == trkID ) {
- mcp = iter->second;
- break;
- }
- }
- return mcp;
- }
-
- /* Find the trackID for the final MCParticle. */
- G4int LcioMcpMaps::findTrackIDFromFinalMcp( MCParticle* mcp) const
- {
- G4int trkID = -1;
- for ( TrackToMcpMap::const_iterator iter = m_trackToMcp.begin();
- iter != m_trackToMcp.end();
- iter++ ) {
- if ( iter->second == mcp ) {
- trkID = iter->first;
- break;
- }
- }
- return trkID;
- }
-
- MCParticle* LcioMcpMaps::findPrimaryInitialMcpFromTrajectory(Trajectory* trj) const
- {
- MCParticle* mcpPrim = 0;
-
- G4int trjTrkID = trj->GetTrackID();
-
- /* Get initial collection. */
- EVENT::LCCollection* initMcpColl = m_manager->getInitialMcpCollection();
-
- // loop over input collection
- int numPrePart = initMcpColl->getNumberOfElements();
- for( int k=0; k < numPrePart; k++ ) {
-
- // current MCP
- MCParticle* mcpPre = static_cast<MCParticle*> ( initMcpColl->getElementAt( k ) );
-
- // do we have matching primary for this Mcp?
- G4PrimaryParticle* primaryMatch = findPrimaryFromMcp( mcpPre );
- if ( primaryMatch ) {
-
- // found matching primary to this input trajectory
- if ( primaryMatch->GetTrackID() == trjTrkID ) {
- mcpPrim = mcpPre;
- break;
- }
-
- }
- }
-
- return mcpPrim;
- }
-
- MCParticle* LcioMcpMaps::findDaughterMcpFromPrimary(MCParticle* mcpInit,
- G4PrimaryParticle* primDau) const
- {
- // loop over Mcp input particle's daughters
- MCParticle* mcpDau = 0;
- int numMcpDau = mcpInit->getDaughters().size();
- for ( int i=0; i < numMcpDau; i++ ) {
-
- // does ith Mcp dau match this primary?
- if ( findPrimaryFromMcp( mcpInit->getDaughters()[i] ) == primDau ) {
-
- // found it
- mcpDau = mcpInit->getDaughters()[i];
- break;
- }
- }
-
- // return dau
- return mcpDau;
- }
-
- MCParticleImpl* LcioMcpMaps::findFinalParticleFromInitial( MCParticle* mcpInit) const
- {
- MCParticleImpl* mcpFinal = 0;
- for ( InitMcpToMcpMap::const_iterator iter = m_initMcpToMcp.begin();
- iter != m_initMcpToMcp.end();
- iter++ ) {
- if ( iter->first == mcpInit ) {
- mcpFinal = iter->second;
- break;
- }
- }
-
- return mcpFinal;
- }
-
- G4PrimaryParticle* LcioMcpMaps::findPrimaryFromMcp(MCParticle* mcp) const
- {
- G4PrimaryParticle *pp = 0;
- for ( McpToPrimaryMap::const_iterator iter = m_McpToPrimary.begin();
- iter != m_McpToPrimary.end();
- iter++ ){
- if ( iter->first == mcp ) {
- pp = iter->second;
- break;
- }
- }
-
- return pp;
- }
-
- void LcioMcpMaps::addMcpToPrimaryLink(MCParticle* mcp,
- G4PrimaryParticle* primary)
- {
- m_McpToPrimary[ mcp ] = primary;
- }
-
- void LcioMcpMaps::addInitialMcpToFinalMcpLink(MCParticle* mcpInit,
- MCParticleImpl* mcpFinal)
- {
- m_initMcpToMcp[ mcpInit ] = mcpFinal;
- }
-
- void LcioMcpMaps::clear()
- {
- m_trackToMcp.clear();
- m_McpToPrimary.clear();
- m_initMcpToMcp.clear();
- }
-
- void LcioMcpMaps::printMaps()
- {
- log() << LOG::debug << "Printing maps..." << LOG::done;
-
- printTrackToMcpMap();
- printMcpToPrimaryMap();
- printInitMcpToMcpMap();
- // TrajectoryManager::instance()->printTrackIDToTrajectoryMap();
- }
-
- void LcioMcpMaps::printTrackToMcpMap()
- {
- log() << LOG::debug << "Printing track to Mcp map..." << LOG::done;
-
- log() << LOG::debug << "*** TrackToMcpMap ***" << LOG::done;
- log() << LOG::debug << "trkID | McpPtr" << LOG::done;
- for ( TrackToMcpMap::iterator iter = m_trackToMcp.begin();
- iter != m_trackToMcp.end();
- iter++ ) {
- log() << LOG::debug << iter->first << " " << iter->second << LOG::done;
- }
- log() << LOG::debug << LOG::done;
- }
-
- void LcioMcpMaps::printMcpToPrimaryMap()
- {
- log() << LOG::debug << "Printing Mcp to primary map..." << LOG::done;
-
- log() << LOG::debug << LOG::done;
- log() << LOG::debug << "*** McpToPrimaryMap ***" << LOG::done;
- log() << LOG::debug << "McpPtr | PrimaryPtr" << LOG::done;
-
- for ( McpToPrimaryMap::iterator iter = m_McpToPrimary.begin();
- iter != m_McpToPrimary.end();
- iter++ ) {
- log() << LOG::debug << iter->first << " " << iter->second << LOG::done;
- }
- log() << LOG::debug << LOG::done;
- }
-
- void LcioMcpMaps::printInitMcpToMcpMap()
- {
- log() << LOG::debug << "Printing InitMcp to Mcp map..." << LOG::done;
-
- log() << LOG::debug << LOG::done;
- log() << LOG::debug << "*** McpToMcpMap ***" << LOG::done;
- log() << LOG::debug << "McpPtrInit | McpPtrFinal" << LOG::done;
-
- for ( InitMcpToMcpMap::iterator iter = m_initMcpToMcp.begin();
- iter != m_initMcpToMcp.end();
- iter++ ) {
- log() << LOG::debug << iter->first << " " << iter->second << LOG::done;
+ LcioMcpMaps::LcioMcpMaps(LcioMcpManager* manager)
+ : Module( "LcioMcpMaps" )
+ {
+ m_manager = manager;
+ }
+
+ LcioMcpMaps::~LcioMcpMaps()
+ {}
+
+ void LcioMcpMaps::addTrackIDToMcpLink(G4int trkID,
+ MCParticleImpl* mcp)
+ {
+ assert( mcp );
+ assert(trkID > 0);
+
+ m_trackToMcp[ trkID ] = mcp;
+ }
+
+ MCParticleImpl* LcioMcpMaps::findMcpFromTrackID( G4int trkID ) const
+ {
+ MCParticleImpl* mcp = 0;
+ for ( TrackToMcpMap::const_iterator iter = m_trackToMcp.begin();
+ iter != m_trackToMcp.end();
+ iter++ )
+ {
+ if ( iter->first == trkID ) {
+ mcp = iter->second;
+ break;
+ }
+ }
+ return mcp;
+ }
+
+ /* Find the trackID for the final MCParticle. */
+ G4int LcioMcpMaps::findTrackIDFromFinalMcp( MCParticle* mcp) const
+ {
+ G4int trkID = -1;
+ for ( TrackToMcpMap::const_iterator iter = m_trackToMcp.begin();
+ iter != m_trackToMcp.end();
+ iter++ )
+ {
+ if ( iter->second == mcp )
+ {
+ trkID = iter->first;
+ break;
+ }
+ }
+ return trkID;
+ }
+
+ MCParticle* LcioMcpMaps::findPrimaryInitialMcpFromTrajectory(Trajectory* trj) const
+ {
+ MCParticle* mcpPrim = 0;
+
+ G4int trjTrkID = trj->GetTrackID();
+
+ /* Get initial collection. */
+ EVENT::LCCollection* initMcpColl = m_manager->getInitialMcpCollection();
+
+ // loop over input collection
+ int numPrePart = initMcpColl->getNumberOfElements();
+ for( int k=0; k < numPrePart; k++ )
+ {
+
+ // current MCP
+ MCParticle* mcpPre = static_cast<MCParticle*> ( initMcpColl->getElementAt( k ) );
+
+ // do we have matching primary for this Mcp?
+ G4PrimaryParticle* primaryMatch = findPrimaryFromMcp( mcpPre );
+ if ( primaryMatch )
+ {
+
+ // found matching primary to this input trajectory
+ if ( primaryMatch->GetTrackID() == trjTrkID )
+ {
+ mcpPrim = mcpPre;
+ break;
+ }
+
+ }
+ }
+
+ return mcpPrim;
+ }
+
+ MCParticle* LcioMcpMaps::findDaughterMcpFromPrimary(MCParticle* mcpInit,
+ G4PrimaryParticle* primDau) const
+ {
+ // loop over Mcp input particle's daughters
+ MCParticle* mcpDau = 0;
+ int numMcpDau = mcpInit->getDaughters().size();
+ for ( int i=0; i < numMcpDau; i++ )
+ {
+ // does ith Mcp dau match this primary?
+ if ( findPrimaryFromMcp( mcpInit->getDaughters()[i] ) == primDau )
+ {
+ // found it
+ mcpDau = mcpInit->getDaughters()[i];
+ break;
+ }
+ }
+
+ // return dau
+ return mcpDau;
+ }
+
+ MCParticleImpl* LcioMcpMaps::findFinalParticleFromInitial( MCParticle* mcpInit) const
+ {
+ MCParticleImpl* mcpFinal = 0;
+ for ( InitMcpToMcpMap::const_iterator iter = m_initMcpToMcp.begin();
+ iter != m_initMcpToMcp.end();
+ iter++ )
+ {
+ if ( iter->first == mcpInit )
+ {
+ mcpFinal = iter->second;
+ break;
+ }
+ }
+
+ return mcpFinal;
+ }
+
+ G4PrimaryParticle* LcioMcpMaps::findPrimaryFromMcp(MCParticle* mcp) const
+ {
+ G4PrimaryParticle *pp = 0;
+ for ( McpToPrimaryMap::const_iterator iter = m_McpToPrimary.begin();
+ iter != m_McpToPrimary.end();
+ iter++ )
+ {
+ if ( iter->first == mcp )
+ {
+ pp = iter->second;
+ break;
+ }
+ }
+
+ return pp;
+ }
+
+ void LcioMcpMaps::addMcpToPrimaryLink(MCParticle* mcp,
+ G4PrimaryParticle* primary)
+ {
+ m_McpToPrimary[ mcp ] = primary;
+ }
+
+ void LcioMcpMaps::addInitialMcpToFinalMcpLink(MCParticle* mcpInit,
+ MCParticleImpl* mcpFinal)
+ {
+ m_initMcpToMcp[ mcpInit ] = mcpFinal;
+ }
+
+ void LcioMcpMaps::clear()
+ {
+ m_trackToMcp.clear();
+ m_McpToPrimary.clear();
+ m_initMcpToMcp.clear();
+ }
+
+ void LcioMcpMaps::printMaps()
+ {
+ log() << LOG::debug << "Printing maps..." << LOG::done;
+
+ printTrackToMcpMap();
+ printMcpToPrimaryMap();
+ printInitMcpToMcpMap();
+ // TrajectoryManager::instance()->printTrackIDToTrajectoryMap();
+ }
+
+ void LcioMcpMaps::printTrackToMcpMap()
+ {
+ log() << LOG::debug << "Printing track to Mcp map..." << LOG::done;
+
+ log() << LOG::debug << "*** TrackToMcpMap ***" << LOG::done;
+ log() << LOG::debug << "trkID | McpPtr" << LOG::done;
+ for ( TrackToMcpMap::iterator iter = m_trackToMcp.begin();
+ iter != m_trackToMcp.end();
+ iter++ )
+ {
+ log() << LOG::debug << iter->first << " " << iter->second << LOG::done;
+ }
+ log() << LOG::debug << LOG::done;
+ }
+
+ void LcioMcpMaps::printMcpToPrimaryMap()
+ {
+ log() << LOG::debug << "Printing Mcp to primary map..." << LOG::done;
+
+ log() << LOG::debug << LOG::done;
+ log() << LOG::debug << "*** McpToPrimaryMap ***" << LOG::done;
+ log() << LOG::debug << "McpPtr | PrimaryPtr" << LOG::done;
+
+ for ( McpToPrimaryMap::iterator iter = m_McpToPrimary.begin();
+ iter != m_McpToPrimary.end();
+ iter++ )
+ {
+ log() << LOG::debug << iter->first << " " << iter->second << LOG::done;
+ }
+ log() << LOG::debug << LOG::done;
+ }
+
+ void LcioMcpMaps::printInitMcpToMcpMap()
+ {
+ log() << LOG::debug << "Printing InitMcp to Mcp map..." << LOG::done;
+
+ log() << LOG::debug << LOG::done;
+ log() << LOG::debug << "*** McpToMcpMap ***" << LOG::done;
+ log() << LOG::debug << "McpPtrInit | McpPtrFinal" << LOG::done;
+
+ for ( InitMcpToMcpMap::iterator iter = m_initMcpToMcp.begin();
+ iter != m_initMcpToMcp.end();
+ iter++ )
+ {
+ log() << LOG::debug << iter->first << " " << iter->second << LOG::done;
+ }
+ log() << LOG::debug << LOG::done;
}
- log() << LOG::debug << LOG::done; - }
}
diff -u -r1.37 -r1.38 --- StdHepToLcioConvertor.cc 26 Aug 2008 23:15:31 -0000 1.37 +++ StdHepToLcioConvertor.cc 16 Dec 2009 00:16:05 -0000 1.38 @@ -1,4 +1,4 @@
-// $Header: /cvs/lcd/slic/src/StdHepToLcioConvertor.cc,v 1.37 2008/08/26 23:15:31 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/StdHepToLcioConvertor.cc,v 1.38 2009/12/16 00:16:05 jeremy Exp $
#include "StdHepToLcioConvertor.hh" // slic
@@ -232,20 +232,20 @@
if ( parVec.size() > 0 ) {
for ( std::vector<int>::iterator it = parVec.begin();
it != parVec.end();
- it++ ) {
- if ( *it != ihep ) {
+ it++ )
+ {
+ if ( *it != ihep )
+ {
addMcpParent( *it, mcp );
}
-#ifdef SLIC_LOG
- else {
- log() << LOG::warning << "Ignoring parent pointing to self for ihep <" << ihep << ">" << LOG::done;
- }
-#endif
- }
- }
//#ifdef SLIC_LOG
-// log() << LOG::debug << "particle has no parents" << LOG::done;
+// else
+// {
+// log() << LOG::debug << "Ignoring parent pointing to self for ihep <" << ihep << ">" << LOG::done;
+// }
//#endif
+ } + }
}
void StdHepToLcioConvertor::setupDaughters( int ihep, IMPL::MCParticleImpl* mcp )
@@ -297,9 +297,12 @@
if ( *it != ihep ) {
addMcpDaughter( *it, mcp );
}
- else {
+#ifdef SLIC_LOG
+ else
+ {
log() << LOG::error << "Ignoring daughter pointing to self!" << LOG::done;
}
+#endif
}
}
#ifdef SLIC_LOG
diff -u -r1.49 -r1.50 --- TrajectoryManager.cc 15 Oct 2009 19:14:17 -0000 1.49 +++ TrajectoryManager.cc 16 Dec 2009 00:16:05 -0000 1.50 @@ -1,4 +1,4 @@
-// $Header: /cvs/lcd/slic/src/TrajectoryManager.cc,v 1.49 2009/10/15 19:14:17 jeremy Exp $
+// $Header: /cvs/lcd/slic/src/TrajectoryManager.cc,v 1.50 2009/12/16 00:16:05 jeremy Exp $
// slic #include "TrajectoryManager.hh"
@@ -29,25 +29,22 @@
#ifdef SLIC_LOG
log() << LOG::verbose
- << "-------------------------------" << LOG::endl
<< "TrajectoryManager::preTracking()" << LOG::endl
- << " trackID <" << aTrack->GetTrackID() << ">" << LOG::endl - << " particle <" << aTrack->GetDefinition()->GetParticleName() << ">" << LOG::endl - << " volume <" << aTrack->GetVolume()->GetName() << ">"
+ << "trackID <" << aTrack->GetTrackID() << ">" << LOG::endl + << "particle <" << aTrack->GetDefinition()->GetParticleName() << ">" << LOG::endl + << "volume <" << aTrack->GetVolume()->GetName() << ">"
<< LOG::done;
- //printStatusFlags(aTrack);
+ printStatusFlags(aTrack);
#endif
// default trajectory storage to OFF
m_trackingManager->SetStoreTrajectory( false );
// handle primary or secondary track
- if ( isPrimaryTrack(aTrack) )
- {
+ if ( isPrimaryTrack(aTrack) ) {
handlePrimaryTrack( aTrack );
}
- else
- {
+ else {
handleSecondaryTrack( aTrack );
}
}
@@ -60,15 +57,14 @@
void TrajectoryManager::printStatusFlags(TrackInformation* trkInfo)
{
- if ( trkInfo )
- {
+ if ( trkInfo ) {
log() << LOG::verbose
- << " orig tkr status <" << TrackInformation::getTrackingStatusString( trkInfo->getOriginalTrackingStatus() ) << ">" << LOG::endl - << " curr tkr status <" << TrackInformation::getTrackingStatusString( trkInfo->getTrackingStatus() ) << ">" << LOG::endl - << " backscatter <" << trkInfo->getBackscatter() << ">" << LOG::endl - << " belowThreshold <" << trkInfo->getBelowThreshold() << ">" << LOG::endl - << " has tracker hit <" << trkInfo->hasTrackerHit() << ">" << LOG::endl - << " vtxnotendp <" << trkInfo->getVertexIsNotEndpointOfParent() << ">" << LOG::endl
+ << "tracking status <" << TrackInformation::getTrackingStatusString( trkInfo->getTrackingStatus() ) << ">" << LOG::endl + << "orig tracking status <" << TrackInformation::getTrackingStatusString( trkInfo->getOriginalTrackingStatus() ) << ">" << LOG::endl + << "backscatter flag <" << trkInfo->getBackscatter() << ">" << LOG::endl + << "belowThreshold flag <" << trkInfo->getBelowThreshold() << ">" << LOG::endl + << "tracker hit flag <" << trkInfo->hasTrackerHit() << ">" << LOG::endl + << "vtxnotendp flag <" << trkInfo->getVertexIsNotEndpointOfParent() << ">"
<< LOG::done;
}
#ifdef SLIC_LOG
@@ -81,15 +77,15 @@
void TrajectoryManager::printTrajectory(const Trajectory* trj)
{
log() << LOG::debug
- << " trajectory info" << LOG::endl - << " trackID <" << trj->GetTrackID() << ">" << LOG::endl - << " final status: " << trj->getFinalStatusString() << LOG::endl - << " endp dau: " << trj->getHasEndpointDaughters() << LOG::endl - << " backscatter: " << trj->getBackscatter() << LOG::endl - << " vtxnotendp: " << trj->getVertexIsNotEndpointOfParent() << LOG::endl - << " created in sim: " << trj->getCreatedInSimulation() << LOG::endl - << " vtx: " << trj->GetPoint(0)->GetPosition() << LOG::endl - << " endp: " << trj->GetPoint(trj->GetPointEntries() - 1)->GetPosition() << LOG::endl
+ << "trajectory info" << LOG::endl + << "trackID <" << trj->GetTrackID() << ">" << LOG::endl + << "final status: " << trj->getFinalStatusString() << LOG::endl + << "endp dau: " << trj->getHasEndpointDaughters() << LOG::endl + << "backscatter: " << trj->getBackscatter() << LOG::endl + << "vtxnotendp: " << trj->getVertexIsNotEndpointOfParent() << LOG::endl + << "created in sim: " << trj->getCreatedInSimulation() << LOG::endl + << "vtx: " << trj->GetPoint(0)->GetPosition() << LOG::endl + << "endp: " << trj->GetPoint(trj->GetPointEntries() - 1)->GetPosition() << LOG::endl
<< LOG::done;
}
@@ -97,58 +93,58 @@
{
#ifdef SLIC_LOG
log() << LOG::verbose << "TrajectoryManager::postTracking()" << LOG::endl;
- log() << LOG::verbose << " store trj <" << m_trackingManager->GetStoreTrajectory() << ">" << LOG::done;
#endif
- // Check if needs manual setup. This can happen if the track is flagged as
- // non-tracking status but then makes some tracker hits.
- if ( needsManualTrajectorySetup( aTrack ) )
- {
+ // check if needs manual setup
+ if ( needsManualTrajectorySetup( aTrack ) ) {
#ifdef SLIC_LOG
- log() << LOG::verbose << " manual trajectory setup activated" << LOG::done;
+ log() << LOG::verbose << "needs manual trajectory setup" << LOG::done;
#endif
- // Setup trajectory in default way.
+ // setup trajectory in default way
Trajectory* trj = createTrajectoryFromTrack( aTrack );
- // Do additional, manual setup using this track.
+ // do additional, manual setup using this track
trj->setupManuallyFromTrack( aTrack );
}
- // Copy track info to secondaries.
+ // trk info to secondaries
copySecondaryTrackInformationFromParent( aTrack );
- // Set the final status codes in the current trajectory.
+ // set final status codes in current trajectory
setTrajectoryFinalStatus( aTrack );
+ // print if store traj on/off
#ifdef SLIC_LOG
- printStatusFlags(aTrack);
+ log() << LOG::verbose << "store trajectory <" << m_trackingManager->GetStoreTrajectory() << ">" << LOG::done;
#endif
+ // print status flags
#ifdef SLIC_LOG
- log() << LOG::verbose << " store trajectory <" << m_trackingManager->GetStoreTrajectory() << ">" << LOG::done;
+ printStatusFlags(aTrack);
#endif
- // Store the endpoint energy.
- Trajectory* trj = static_cast<Trajectory*>( m_trackingManager->GimmeTrajectory() );
- if ( 0 != trj )
- {
+ Trajectory* trj = static_cast<Trajectory*>(m_trackingManager->GimmeTrajectory());
+ if ( 0 != trj ) {
trj->setEndpointEnergy( aTrack->GetTotalEnergy() );
}
}
void TrajectoryManager::handlePrimaryTrack(const G4Track* aTrack)
{
- // Set the track info.
+ // set track info
setPrimaryTrackInformation(aTrack);
- // Create a trajectory for this track.
+ // create the trajectory
createTrajectoryFromTrack( aTrack );
}
Trajectory* TrajectoryManager::createTrajectoryFromTrack(const G4Track* aTrack )
{
+#ifdef SLIC_LOG + log() << LOG::verbose << "createTrajectoryFromTrack()" << LOG::done; +#endif
Trajectory* trj = new Trajectory( aTrack );
m_trackingManager->SetStoreTrajectory( true );
m_trackingManager->SetTrajectory( trj );
@@ -170,12 +166,10 @@
void TrajectoryManager::setTrackingFlagFromRegionInformation(TrackInformation* trk_info, G4UserRegionInformation* reg_info)
{
- if ( reg_info->getStoreSecondaries() == true )
- {
+ if ( reg_info->getStoreSecondaries() == true ) {
trk_info->setTrackingStatus(TrackInformation::eInTrackingRegion);
}
- else
- {
+ else {
trk_info->setTrackingStatus(TrackInformation::eInNontrackingRegion);
}
}
@@ -187,10 +181,6 @@
*/
void TrajectoryManager::handleSecondaryTrack(const G4Track* aTrack)
{
-#ifdef SLIC_LOG - log() << LOG::verbose << "handleSecondaryTrack <" << aTrack->GetTrackID() << ">" << LOG::endl; -#endif -
// get track information
TrackInformation* trkInfo = TrackInformation::getTrackInformation( aTrack );
@@ -198,78 +188,75 @@
G4UserRegionInformation* regInfo = G4UserRegionInformation::getRegionInformation ( aTrack );
// must have track and region info to proceed
- if ( regInfo && trkInfo )
- {
+ if ( regInfo && trkInfo ) {
- // Store trajectories if in tracking region.
- if ( regInfo->getStoreSecondaries() )
- {
+ // check if in tracking region
+ if ( regInfo->getStoreSecondaries() ) {
// check if backscatter
- if ( isBackscatter( aTrack ) )
- {
-#ifdef SLIC_LOG
- log() << LOG::verbose << " backscatter in tracking" << LOG::endl;
-#endif
+ if ( isBackscatter( aTrack ) ) {
trkInfo->setBackscatter( true );
}
- // Check if below threshold and continuous process creation, i.e. not an endp dau.
+ /* + Check if below threshold and continuous process creation, + i.e. not an endp dau. + */
else if ( isBelowThreshold( aTrack ) &&
- trkInfo->getVertexIsNotEndpointOfParent() )
- {
-#ifdef SLIC_LOG
- log() << LOG::verbose << " backscatter + vtxnotendp in tracking" << LOG::endl;
-#endif
+ trkInfo->getVertexIsNotEndpointOfParent() ) {
trkInfo->setBelowThreshold( true );
}
- // Track requires a new trajectory and trkInfo.
- // FIXME: Using threshold causes missing mcp contribs due to breaking parentage chain.
- //if ( !isBelowThreshold( aTrack ) )
- else
- {
-#ifdef SLIC_LOG
- log() << LOG::verbose << " setup up trajectory in tracking" << LOG::endl;
-#endif
- // Setup the secondary's track information.
+ /*
+ Process of elimination: above threshold and not backscatter,
+ so requires a new trajectory and trkInfo.
+ */
+ else {
+
+ // setup secondary track info
setupSecondaryTrackInformation( aTrack );
- // Create the trajectory for this secondary.
+ // create the secondary trajectory
createSecondaryTrajectory( aTrack );
}
-// else
-// {
-// if ( isBelowThreshold( aTrack ) )
-// {
-// #ifdef SLIC_LOG
-// log() << LOG::debug << " WARNING: Track is below threshold!" << LOG::endl;
-// #endif
-// }
-// }
- // Set flag for type of region. Need to get track info again, because it may have been recreated above. - TrackInformation::getTrackInformation( aTrack )->setTrackingStatus( TrackInformation::eInTrackingRegion );
+
+ /*
+ // DEBUG
+ if ( trkInfo != TrackInformation::getTrackInformation( aTrack ) ) {
+ log() << "WARNING: track info was reset" << LOG::endl;
+
+ log() << LOG::endl;
+ log() << "old flags" << LOG::endl;
+ printStatusFlags( trkInfo );
+ log() << LOG::endl;
+
+ log() << "new flags" << LOG::endl;
+ printStatusFlags( TrackInformation::getTrackInformation( aTrack ) );
+ log() << LOG::endl;
+ }
+ */
+
+ /*
+ Set flag for in tracking region, re-getting track info,
+ because may have been recreated above.
+ */
+ TrackInformation::getTrackInformation( aTrack )
+ ->setTrackingStatus( TrackInformation::eInTrackingRegion );
}
- // Track starts in a non-tracking region and will not be stored unless it makes a tracker hit.
- else
- {
-#ifdef SLIC_LOG
- log() << LOG::verbose << " non-tracking and no trj" << LOG::endl;
-#endif
- // Set status to non-tracking.
+ // in calorimeter, i.e. non-tracking
+ else {
+
+ // set non tracking status
trkInfo->setTrackingStatus( TrackInformation::eInNontrackingRegion );
trkInfo->setOriginalTrackingStatus( TrackInformation::eInNontrackingRegion );
}
}
- // Fatal error: region info or track info was null.
- else
- {
-#ifdef SLIC_LOG
- log() << LOG::fatal << "FATAL ERROR: Region info or track info were null in handleSecondaryTrack()." << LOG::endl
+ // fatal error: regInfo or trkInfo were null
+ else {
+ log() << LOG::fatal << "FATAL ERROR: region info or track info were null in handleSecondaryTrack()." << LOG::endl
<< "region info <" << regInfo << ">" << LOG::endl
<< "track info <" << trkInfo << ">" << LOG::done;
-#endif
G4Exception("Track info or region info was null.");
}
}
@@ -329,11 +316,9 @@
G4TrackVector* secondaries = m_trackingManager->GimmeSecondaries();
// loop over secondaries to 1) create trk info and 2) set vertexIsNotEndpointOfParent flag
- if ( secondaries )
- {
+ if ( secondaries ) {
size_t nseco = secondaries->size();
- for ( size_t i = 0; i < nseco; i++ )
- {
+ for ( size_t i = 0; i < nseco; i++ ) {
// set current seco track ptr
G4Track* seco = (*secondaries)[i];
@@ -342,30 +327,28 @@
TrackInformation* secoTrkInfo = new TrackInformation( parTrkInfo );
seco->SetUserInformation( secoTrkInfo );
- // Set vtxnotendp flag in track info.
- if ( parTrj )
- {
+ // set vtxnotendp flag in trk info
+ if ( parTrj ) {
- // Get parent's endpoint.
+ // get parent's endpoint
G4ThreeVector parEndpoint = getTrajectoryEndpoint( parTrj );
- // Check if direct parent (flag always on if not) OR endp not near.
+ // check if direct parent (flag always on if not) OR endp not near
if ( secoTrkInfo->getOriginalTrackID() != aTrack->GetTrackID() ||
- vertexIsNotEndpointOfParent( parEndpoint, seco ) )
- {
- // Set vertexIsNotEndpointOfParent flag in trk info.
+ vertexIsNotEndpointOfParent( parEndpoint, seco ) ) {
+
+ // set vertexIsNotEndpointOfParent flag in trk info
secoTrkInfo->setVertexIsNotEndpointOfParent( true );
}
}
- // In this case, vtxnotendp is always set, because parent was not stored.
- else
- {
+ // vtxnotendp always set, because parent is not stored
+ else {
secoTrkInfo->setVertexIsNotEndpointOfParent( true );
}
}
#ifdef SLIC_LOG
- log() << LOG::verbose << " track info passed to <" << nseco << "> secondaries" << LOG::endl;
+ log() << LOG::verbose << "track info passed to <" << nseco << "> secondaries" << LOG::endl;
#endif
}
}
@@ -373,7 +356,8 @@
bool TrajectoryManager::needsManualTrajectorySetup(const G4Track* aTrack)
{
// store != ON and has a corresponding tracker hit
- return ( !m_trackingManager->GetStoreTrajectory() ) && ( TrackInformation::getTrackInformation( aTrack )->hasTrackerHit() );
+ return ( !m_trackingManager->GetStoreTrajectory() ) && + ( TrackInformation::getTrackInformation( aTrack )->hasTrackerHit() );
}
void TrajectoryManager::stepping(const G4Step* aStep)
@@ -381,8 +365,7 @@
const G4Track* theTrack = aStep->GetTrack();
// step has a live track
- if ( isAlive( aStep ) )
- {
+ if ( isAlive( aStep ) ) {
TrackInformation* trkInfo = TrackInformation::getTrackInformation( theTrack );
@@ -390,22 +373,18 @@
G4UserRegionInformation* postReg = G4UserRegionInformation::getRegionInformation( aStep->GetPostStepPoint() );
// set tracking status based on region info
- if ( preReg && postReg )
- {
+ if ( preReg && postReg ) {
// tracking to non-tracking
- if ( preReg->getStoreSecondaries() && ( !postReg->getStoreSecondaries() ) )
- {
+ if ( preReg->getStoreSecondaries() && ( !postReg->getStoreSecondaries() ) ) {
trkInfo->setTrackingStatus( TrackInformation::eInNontrackingRegion );
}
// non-tracking to tracking
- else if ( !preReg->getStoreSecondaries() && postReg->getStoreSecondaries() )
- {
+ else if ( !preReg->getStoreSecondaries() && postReg->getStoreSecondaries() ) {
TrackInformation::ETrackingStatus origStatus = trkInfo->getOriginalTrackingStatus();
// if looped like tracking -> nontracking -> tracking without dying
- if ( ! ( origStatus == TrackInformation::eInTrackingRegion ) )
- {
+ if ( ! ( origStatus == TrackInformation::eInTrackingRegion ) ) {
trkInfo->setBackscatter(true);
}
@@ -422,8 +401,8 @@
Trajectory* trj = getCurrentTrajectory();
// must have trajectory to set status
- if ( trj )
- {
+ if ( trj ) {
+
G4StepPoint* postStep = aTrack->GetStep()->GetPostStepPoint();
// get track info
@@ -432,8 +411,7 @@
// has endpoint daughters
G4int nSecoAtEnd = getNumberOfSecondariesAtEnd( aTrack );
- if ( nSecoAtEnd > 0 )
- {
+ if ( nSecoAtEnd > 0 ) {
trj->setHasEndpointDaughters( true );
}
@@ -444,125 +422,92 @@
trj->setBackscatter( trkInfo->getBackscatter() );
// vertexIsNotEndpointOfParent
- if ( trkInfo->getVertexIsNotEndpointOfParent() )
- {
+ if ( trkInfo->getVertexIsNotEndpointOfParent() ) {
trj->setVertexIsNotEndpointOfParent( true );
}
// no decay -> no secondaries at end
- if ( nSecoAtEnd <= 0 )
- {
+ if ( nSecoAtEnd <= 0 ) {
// left detector
- if ( !postStep->GetMaterial() )
- {
-#ifdef SLIC_LOG
- log() << LOG::debug << " left detector" << LOG::done;
-#endif
+ if ( !postStep->GetMaterial() ) {
finalStatus = Trajectory::eLeftDetector;
}
// stopped
- else if ( aTrack->GetKineticEnergy() <= 0. )
- {
-#ifdef SLIC_LOG
- log() << LOG::debug << " stopped" << LOG::done;
-#endif
+ else if ( aTrack->GetKineticEnergy() <= 0. ) {
finalStatus = Trajectory::eStopped;
}
// unflagged loopers killed by G4
#ifdef SLIC_DEBUG
- else
- {
- log() << LOG::debug << " unknown decay case (probably a looper)" << LOG::done;
+ else {
+ log() << LOG::debug << "Unknown decay case (probably a looper)." << LOG::done;
}
#endif
}
// decayed
- else
- {
+ else {
+
// decayed in tracker
- if ( trkInfo->getTrackingStatus() == TrackInformation::eInTrackingRegion )
- {
-#ifdef SLIC_LOG
- log() << LOG::debug << " decayed in tracker" << LOG::done;
-#endif
+ if ( trkInfo->getTrackingStatus() == TrackInformation::eInTrackingRegion ) {
finalStatus = Trajectory::eDecayedInTracker;
}
// process of elimination: decayed in calorimeter
- else
- {
-#ifdef SLIC_LOG
- log() << LOG::debug << " decayed in calorimeter" << LOG::done;
-#endif
+ else {
finalStatus = Trajectory::eDecayedInCalorimeter;
}
}
// set final status in trajectory
trj->setFinalStatus( finalStatus );
- -#ifdef SLIC_LOG - printTrack( aTrack, true ); -#endif
}
}
void TrajectoryManager::printTrack(const G4Track* trk, bool isSecondary)
{
- log() << LOG::debug << LOG::endl; - log() << LOG::debug - << "track ptr <" << trk << ">" << LOG::endl - << "trackID <" << trk->GetTrackID() << ">" << LOG::endl - << "parent trackID <" << trk->GetParentID() << ">" << LOG::endl - << "PDG <" << trk->GetDefinition()->GetPDGEncoding() << ">" << LOG::endl - << "position <" << trk->GetPosition() << ">" << LOG::endl;
+ log() << LOG::endl; + log() << "track ptr <" << trk << ">" << LOG::endl; + log() << "trackID <" << trk->GetTrackID() << ">" << LOG::endl; + log() << "parent trackID <" << trk->GetParentID() << ">" << LOG::endl; + log() << "PDG <" << trk->GetDefinition()->GetPDGEncoding() << ">" << LOG::endl; + log() << "position <" << trk->GetPosition() << ">" << LOG::endl;
- if ( !isSecondary )
- {
- if ( trk->GetMaterial() )
- {
- log() << LOG::debug << "material <" << trk->GetMaterial()->GetName() << ">" << LOG::endl;
+ if ( !isSecondary ) {
+ if ( trk->GetMaterial() ) {
+ log() << "material <" << trk->GetMaterial()->GetName() << ">" << LOG::endl;
}
- else
- {
- log() << LOG::debug << "NO material" << LOG::endl;
+ else {
+ log() << "NO material" << LOG::endl;
}
- if ( trk->GetVolume() )
- {
- log() << LOG::debug << "volume <" << trk->GetVolume()->GetName() << ">" << LOG::endl;
+ if ( trk->GetVolume() ) {
+ log() << "volume <" << trk->GetVolume()->GetName() << ">" << LOG::endl;
}
- else
- {
- log() << LOG::debug << "NO volume" << LOG::endl;
+ else {
+ log() << "NO volume" << LOG::endl;
}
}
- if ( !isSecondary )
- {
- log() << LOG::debug << "numSecoAtEnd <" << getNumberOfSecondariesAtEnd(trk) << ">" << LOG::endl;
+ if ( !isSecondary ) {
+ log() << "numSecoAtEnd <" << getNumberOfSecondariesAtEnd(trk) << ">" << LOG::endl;
}
- log() << LOG::debug << "kE <" << trk->GetKineticEnergy() << ">" << LOG::endl - << "mom <" << trk->GetMomentum() << ">" << LOG::endl;
+ log() << "kE <" << trk->GetKineticEnergy() << ">" << LOG::endl; + log() << "mom <" << trk->GetMomentum() << ">" << LOG::endl;
printStatusFlags(trk);
- log() << LOG::debug << LOG::endl;
+ log() << LOG::endl;
}
void TrajectoryManager::printPhysVolDbg(G4StepPoint* stepPnt)
{
- if ( stepPnt )
- {
+ if ( stepPnt ) {
G4VPhysicalVolume* pv = stepPnt->GetPhysicalVolume();
- if ( pv )
- {
+ if ( pv ) {
log() << LOG::debug << "PV name <" << pv->GetName() << ">" << LOG::endl;
}
- else
- {
+ else {
log() << LOG::debug << "PV is NULL!" << LOG::endl;
}
}
@@ -575,19 +520,17 @@
G4TrackVector* seco = m_trackingManager->GimmeSecondaries();
G4ThreeVector endPnt = aTrack->GetPosition();
- if ( seco )
- {
+ if ( seco ) {
size_t numSeco = seco->size(); #ifdef SLIC_LOG
- log() << LOG::debug << " numSeco <" << numSeco << ">" << LOG::endl;
+ log() << LOG::debug << "numSeco <" << numSeco << ">" << LOG::endl;
#endif
for( size_t i = 0;
i < numSeco;
- i++ )
- {
- // If seco's endpoint is near track's endpoint, increment counter.
+ i++ ) {
+ // if seco endpoint near track's endpoint
if ( (*seco)[i]->GetPosition().isNear( endPnt ) ) {
++numSecoAtEnd;
}
@@ -595,7 +538,7 @@
} #ifdef SLIC_LOG
- log() << LOG::debug << " numSecoAtEndp <" << numSecoAtEnd << ">" << LOG::endl;
+ log() << LOG::debug << "found numSecoAtEnd <" << numSecoAtEnd << ">" << LOG::endl;
#endif
return numSecoAtEnd;
@@ -604,6 +547,9 @@
void TrajectoryManager::fillTrackIDToTrajectoryMap( const G4Event* anEvent, bool clearIt)
{
+#ifdef SLIC_LOG + log() << LOG::verbose << "fillTrackIDToTrajectoryMap()" << LOG::endl; +#endif
if ( clearIt )
{
clearTrackIDToTrajectoryMap();
@@ -630,12 +576,10 @@
void TrajectoryManager::addTrackIDToTrajectoryLink( Trajectory* trj)
{
- if ( trj )
- {
+ if ( trj ) {
m_trackIDToTrajectory[ trj->GetTrackID() ] = trj;
}
- else
- {
+ else {
log() << LOG::error << "Trajectory is null!" << LOG::endl;
}
}
@@ -650,13 +594,11 @@
G4TrackVector* secondaries = m_trackingManager->GimmeSecondaries();
size_t nseco = secondaries->size();
- for ( size_t i = 0; i < nseco; i++ )
- {
+ for ( size_t i = 0; i < nseco; i++ ) {
+
// set current seco track ptr
G4Track* seco = (*secondaries)[i];
-#ifdef SLIC_LOG
log() << LOG::debug << "secondary track <" << i << ">" << LOG::endl;
-#endif
printTrack( seco, true );
}
}
@@ -664,30 +606,27 @@
// find trajectory by track ID
Trajectory* TrajectoryManager::findTrajectory(G4int trkID)
{
-//#ifdef SLIC_LOG -// log() << LOG::debug << "findTrajectory() - trkID <" << trkID << ">" << LOG::endl; -//#endif
+#ifdef SLIC_LOG + log() << LOG::debug << "findTrajectory() - trkID <" << trkID << ">" << LOG::endl; +#endif
Trajectory* trj = 0;
TrackIDToTrajectoryMap::iterator iter;
iter = m_trackIDToTrajectory.find( trkID );
- if ( iter != m_trackIDToTrajectory.end() )
- {
+ if ( iter != m_trackIDToTrajectory.end() ) {
trj = iter->second;
}
-//#ifdef SLIC_LOG
-// if ( trj )
-// {
-// log() << LOG::debug << "found trajectory" << LOG::endl;
-// }
-// else
-// {
-// log() << LOG::debug << "no trajectory found" << LOG::endl;
-// }
-//#endif
+#ifdef SLIC_LOG
+ if ( trj ) {
+ log() << LOG::debug << "found trajectory" << LOG::endl;
+ }
+ else {
+ log() << LOG::debug << "no trajectory found" << LOG::endl;
+ }
+#endif
return trj;
}