2 added + 3 removed, total 5 files
slic/include
diff -N LcioHitsCollectionBuilder.hh
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ LcioHitsCollectionBuilder.hh 20 Sep 2005 20:15:21 -0000 1.1
@@ -0,0 +1,123 @@
+// $Header: /cvs/lcd/slic/include/LcioHitsCollectionBuilder.hh,v 1.1 2005/09/20 20:15:21 jeremy Exp $
+
+#ifndef slic_LcioHitsCollectionBuilder_hh
+#define slic_LcioHitsCollectionBuilder_hh
+
+// LCDD
+#include "G4SensitiveDetector.hh"
+#include "G4CalorimeterHit.hh"
+#include "G4TrackerHit.hh"
+#include "Verbose.hh"
+
+// LCIO
+#include "EVENT/LCEvent.h"
+#include "IMPL/LCFlagImpl.h"
+#include "IMPL/LCEventImpl.h"
+#include "IMPL/LCCollectionVec.h"
+#include "IMPL/SimCalorimeterHitImpl.h"
+#include "IMPL/SimTrackerHitImpl.h"
+
+// G4
+#include "G4Event.hh"
+
+namespace slic
+{
+
+ class LcioMcpManager;
+
+ /**
+ @class LcioHitsCollectionBuilder
+ @brief Builds Lcio hits collections from G4/SLIC event data.
+ @note Makes use of LcioMcpManager for Mcp hit contributions.
+ */
+ class LcioHitsCollectionBuilder : public Verbose
+ {
+ public:
+ LcioHitsCollectionBuilder();
+ virtual ~LcioHitsCollectionBuilder();
+
+ public:
+
+ // interface function to create an Lcio event from a G4 event
+ EVENT::LCEvent* createHCsFromG4Event(const G4Event* g4evt, EVENT::LCEvent* lcevt);
+
+ // flags for hits
+ void setLongFlag(bool s = true);
+ void setPDGFlag(bool s = true);
+
+ void setStoreMomentum( bool s = true )
+ {
+ m_storeMomentum = s;
+ }
+
+ bool getStoreMomentum()
+ {
+ return m_storeMomentum;
+ }
+
+ private:
+ // create the hit collections in Lcio Event from the current G4Event
+ void createHitCollections();
+
+ // create the LCCollectionVec
+ IMPL::LCCollectionVec* createCollectionVec(G4VHitsCollection* g4HC,
+ G4SensitiveDetector::EType SDtype);
+
+ // create a tracker coll
+ IMPL::LCCollectionVec* createTrackerCollectionVec(G4VHitsCollection* g4HC);
+
+ // create a calorimeter coll
+ IMPL::LCCollectionVec* createCalorimeterCollectionVec(G4VHitsCollection* g4HC);
+
+ // save cal hits
+ void saveHits(G4CalorimeterHitsCollection* calHits,
+ IMPL::LCCollectionVec* lcioColl);
+
+ // save trk hits
+ void saveHits(G4TrackerHitsCollection* trkHits,
+ IMPL::LCCollectionVec* lcioColl);
+
+ // create cal hit from G4 cal hit
+ IMPL::SimCalorimeterHitImpl* createHit(G4CalorimeterHit* calHit);
+
+ // create trk hit from G4 trk hit
+ IMPL::SimTrackerHitImpl* createHit(G4TrackerHit* trkHit);
+
+ // copy MCParticle hit contributions from G4 cal hit to the LCIO cal hit
+ void addMcpContribs( G4CalorimeterHit*, IMPL::SimCalorimeterHitImpl* );
+
+ // setup default flag settings; called in ctor
+ void setCalFlagDefaults();
+
+ // combine Mcp hit contribs using track ID
+ void combineMcpHitContribs(const McpHitContribList& long_contrib,
+ McpHitContribList& combined_contrib);
+
+ // find single contrib by track ID in the combined list
+ McpHitContrib* findMcpHitContribByTrackID(int trk_id,
+ McpHitContribList& contribs);
+
+ // set endcap bit in trk or cal flag depending on SD type
+ void setEndcapFlag(G4SensitiveDetector* g4sd);
+
+ private:
+
+ // current Lcio event
+ EVENT::LCEvent* m_currentLCEvent;
+
+ // current G4 event
+ const G4Event* m_currentG4Event;
+
+ // default flags
+ IMPL::LCFlagImpl m_trkCollFlag;
+ IMPL::LCFlagImpl m_calCollFlag;
+
+ // store momentum
+ bool m_storeMomentum;
+
+ // ptr to Lcio Mcp manager for lkp of Mcp data
+ LcioMcpManager* m_mcpManager;
+ };
+};
+
+#endif
slic/include
diff -N HitsCollectionUtil.hh
--- HitsCollectionUtil.hh 16 Sep 2005 23:22:46 -0000 1.2
+++ /dev/null 1 Jan 1970 00:00:00 -0000
@@ -1,25 +0,0 @@
-// $Header: /cvs/lcd/slic/include/Attic/HitsCollectionUtil.hh,v 1.2 2005/09/16 23:22:46 jeremy Exp $
-#ifndef slic_HitsCollectionUtil_hh
-#define slic_HitsCollectionUtil_hh 1
-
-#include <vector>
-#include <string>
-
-namespace slic
-{
- /**
- * @class HitsCollectionUtil
- * @brief Static utility methods for Hits Collections.
- */
- class HitsCollectionUtil
- {
- private:
- HitsCollectionUtil();
-
- public:
- static std::vector<int> getHCIDs();
- static std::vector<std::string> getHCNames();
- };
-};
-
-#endif
slic/include
diff -N LcioHCBuilder.hh
--- LcioHCBuilder.hh 27 Jun 2005 19:17:37 -0000 1.15
+++ /dev/null 1 Jan 1970 00:00:00 -0000
@@ -1,123 +0,0 @@
-// $Header: /cvs/lcd/slic/include/Attic/LcioHCBuilder.hh,v 1.15 2005/06/27 19:17:37 jeremy Exp $
-
-#ifndef slic_LcioHCBuilder_hh
-#define slic_LcioHCBuilder_hh
-
-// LCDD
-#include "G4SensitiveDetector.hh"
-#include "G4CalorimeterHit.hh"
-#include "G4TrackerHit.hh"
-#include "Verbose.hh"
-
-// LCIO
-#include "EVENT/LCEvent.h"
-#include "IMPL/LCFlagImpl.h"
-#include "IMPL/LCEventImpl.h"
-#include "IMPL/LCCollectionVec.h"
-#include "IMPL/SimCalorimeterHitImpl.h"
-#include "IMPL/SimTrackerHitImpl.h"
-
-// G4
-#include "G4Event.hh"
-
-namespace slic
-{
-
- class LcioMcpManager;
-
- /**
- @class LcioHCBuilder
- @brief Builds Lcio hits collections from G4/SLIC event data.
- @note Makes use of LcioMcpManager for Mcp hit contributions.
- */
- class LcioHCBuilder : public Verbose
- {
- public:
- LcioHCBuilder();
- virtual ~LcioHCBuilder();
-
- public:
-
- // interface function to create an Lcio event from a G4 event
- EVENT::LCEvent* createHCsFromG4Event(const G4Event* g4evt, EVENT::LCEvent* lcevt);
-
- // flags for hits
- void setLongFlag(bool s = true);
- void setPDGFlag(bool s = true);
-
- void setStoreMomentum( bool s = true )
- {
- m_storeMomentum = s;
- }
-
- bool getStoreMomentum()
- {
- return m_storeMomentum;
- }
-
- private:
- // create the hit collections in Lcio Event from the current G4Event
- void createHitCollections();
-
- // create the LCCollectionVec
- IMPL::LCCollectionVec* createCollectionVec(G4VHitsCollection* g4HC,
- G4SensitiveDetector::EType SDtype);
-
- // create a tracker coll
- IMPL::LCCollectionVec* createTrackerCollectionVec(G4VHitsCollection* g4HC);
-
- // create a calorimeter coll
- IMPL::LCCollectionVec* createCalorimeterCollectionVec(G4VHitsCollection* g4HC);
-
- // save cal hits
- void saveHits(G4CalorimeterHitsCollection* calHits,
- IMPL::LCCollectionVec* lcioColl);
-
- // save trk hits
- void saveHits(G4TrackerHitsCollection* trkHits,
- IMPL::LCCollectionVec* lcioColl);
-
- // create cal hit from G4 cal hit
- IMPL::SimCalorimeterHitImpl* createHit(G4CalorimeterHit* calHit);
-
- // create trk hit from G4 trk hit
- IMPL::SimTrackerHitImpl* createHit(G4TrackerHit* trkHit);
-
- // copy MCParticle hit contributions from G4 cal hit to the LCIO cal hit
- void addMcpContribs( G4CalorimeterHit*, IMPL::SimCalorimeterHitImpl* );
-
- // setup default flag settings; called in ctor
- void setCalFlagDefaults();
-
- // combine Mcp hit contribs using track ID
- void combineMcpHitContribs(const McpHitContribList& long_contrib,
- McpHitContribList& combined_contrib);
-
- // find single contrib by track ID in the combined list
- McpHitContrib* findMcpHitContribByTrackID(int trk_id,
- McpHitContribList& contribs);
-
- // set endcap bit in trk or cal flag depending on SD type
- void setEndcapFlag(G4SensitiveDetector* g4sd);
-
- private:
-
- // current Lcio event
- EVENT::LCEvent* m_currentLCEvent;
-
- // current G4 event
- const G4Event* m_currentG4Event;
-
- // default flags
- IMPL::LCFlagImpl m_trkCollFlag;
- IMPL::LCFlagImpl m_calCollFlag;
-
- // store momentum
- bool m_storeMomentum;
-
- // ptr to Lcio Mcp manager for lkp of Mcp data
- LcioMcpManager* m_mcpManager;
- };
-};
-
-#endif
slic/src
diff -N LcioHitsCollectionBuilder.cc
--- /dev/null 1 Jan 1970 00:00:00 -0000
+++ LcioHitsCollectionBuilder.cc 20 Sep 2005 20:15:22 -0000 1.1
@@ -0,0 +1,507 @@
+// $Header: /cvs/lcd/slic/src/LcioHitsCollectionBuilder.cc,v 1.1 2005/09/20 20:15:22 jeremy Exp $
+#include "LcioHitsCollectionBuilder.hh"
+
+// SLIC
+#include "LcioMcpManager.hh"
+#include "HitsCollectionUtil.hh"
+#include "G4Application.hh"
+#include "TrackerHitMomentum.hh"
+
+// LCIO
+#include "IMPL/LCCollectionVec.h"
+#include "IMPL/LCFlagImpl.h"
+#include "IMPL/MCParticleImpl.h"
+
+// G4
+#include "G4SDManager.hh"
+
+// std
+#include <iostream>
+
+using IMPL::SimCalorimeterHitImpl;
+using IMPL::SimTrackerHitImpl;
+using IMPL::MCParticleImpl;
+using IMPL::LCEventImpl;
+using IMPL::LCFlagImpl;
+using IMPL::LCCollectionVec;
+using EVENT::LCIO;
+
+using namespace std;
+
+namespace slic
+{
+
+ LcioHitsCollectionBuilder::LcioHitsCollectionBuilder()
+ : m_storeMomentum( 0 )
+ {
+ // set local LcioMcpManager ptr
+ m_mcpManager = LcioMcpManager::instance();
+
+ // setup default coll flag for cal hits
+ setCalFlagDefaults();
+ }
+
+ LcioHitsCollectionBuilder::~LcioHitsCollectionBuilder()
+ {}
+
+ // create the hit collections
+ void LcioHitsCollectionBuilder::createHitCollections()
+ {
+ // fetch HCIDs
+ std::vector<int> hcids = HitsCollectionUtil::getHCIDs();
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 1 ) {
+ std::cout << "hcids size <" << hcids.size() << ">" << std::endl;
+ }
+#endif
+
+ // fetch hit collection of event
+ G4HCofThisEvent* HCE = m_currentG4Event->GetHCofThisEvent();
+
+ // HC table
+ G4HCtable* HCtbl = G4SDManager::GetSDMpointer()->GetHCtable();
+
+ // HCID
+ G4int hcid;
+
+ LCCollectionVec* collVect = 0;
+
+ for (std::vector<int>::const_iterator iter = hcids.begin();
+ iter != hcids.end();
+ iter++) {
+
+ hcid = *iter;
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 1 ) {
+ std::cout << "HCID: " << hcid << std::endl;
+ }
+#endif
+
+ // retrieve Sensitive Detector ptr
+ G4SensitiveDetector* SD =
+ static_cast<G4SensitiveDetector*>
+ ( G4SDManager::GetSDMpointer()->FindSensitiveDetector( HCtbl->GetSDname( hcid ) ) );
+
+ // get hits collection
+ G4VHitsCollection* HC = HCE->GetHC( hcid );
+
+ // add a LCCollectionVec if got a HC
+ if ( HC ) {
+
+ // set LCIO endcap bit in the flag (i.e. CHBIT_BARREL ) according to SDs setting
+ setEndcapFlag( SD );
+
+ // create collection vector based on type of SD
+ collVect = createCollectionVec( HC, SD->getType() );
+
+ // add collection vector to LCEvent
+ m_currentLCEvent->addCollection( collVect, HC->GetName() );
+ }
+ else {
+ G4Exception( "LcioHitsCollectionBuilder::createHitCollections() - No collection found for Hits Collection ID");
+ }
+ }
+ }
+
+ // create the CollectionVec (decides which overloaded subfunction to call)
+ IMPL::LCCollectionVec* LcioHitsCollectionBuilder::createCollectionVec(G4VHitsCollection* g4HC,
+ G4SensitiveDetector::EType SDtype)
+ {
+ // vec to create
+ LCCollectionVec* collVec = 0;
+
+ // cal hits
+ if ( SDtype == G4SensitiveDetector::eCalorimeter ) {
+ collVec = createCalorimeterCollectionVec(g4HC);
+ }
+ // tracker hits
+ else if ( SDtype == G4SensitiveDetector::eTracker ) {
+ collVec = createTrackerCollectionVec(g4HC);
+ }
+ // unknown type of hit
+ else {
+ G4Exception( "Unknown HC type." );
+ }
+
+ return collVec;
+ }
+
+ LCCollectionVec* LcioHitsCollectionBuilder::createTrackerCollectionVec(G4VHitsCollection* g4HC)
+ {
+ // create Lcio tracker coll
+ LCCollectionVec* collVec = new LCCollectionVec( LCIO::SIMTRACKERHIT );
+
+ // cast to G4 trk HC
+ G4TrackerHitsCollection* trkHits =
+ dynamic_cast<G4TrackerHitsCollection*> (g4HC);
+
+ // call overloaded save function for trk hits
+ saveHits(trkHits,
+ collVec);
+
+ // save the momentum in parallel vec
+ if ( m_storeMomentum ) {
+ TrackerHitMomentum::makeTrackerHitMomenta(trkHits, m_currentLCEvent);
+ }
+
+ // set trk flags
+ collVec->setFlag( m_trkCollFlag.getFlag() );
+
+ return collVec;
+ }
+
+ LCCollectionVec* LcioHitsCollectionBuilder::createCalorimeterCollectionVec(G4VHitsCollection* g4HC)
+ {
+ // create Lcio cal coll
+ LCCollectionVec* collVec = new LCCollectionVec(LCIO::SIMCALORIMETERHIT);
+
+ // cast to G4 cal HC
+ G4CalorimeterHitsCollection* calHits
+ = dynamic_cast<G4CalorimeterHitsCollection*> (g4HC);
+
+ // call overloaded save function for cal hits
+ saveHits(calHits,
+ collVec);
+
+ // set cal flags
+ collVec->setFlag( m_calCollFlag.getFlag() );
+
+ return collVec;
+ }
+
+ void LcioHitsCollectionBuilder::setCalFlagDefaults()
+ {
+ m_calCollFlag.setBit( LCIO::CHBIT_LONG );
+ m_calCollFlag.setBit( LCIO::CHBIT_ID1 );
+ }
+
+ void LcioHitsCollectionBuilder::setEndcapFlag(G4SensitiveDetector* g4sd)
+ {
+ bool ec_flag = g4sd->getEndcapFlag();
+
+ // set for cal
+ if ( g4sd->getType() == G4SensitiveDetector::eCalorimeter ) {
+ if ( ec_flag ) {
+ m_calCollFlag.unsetBit( LCIO::CHBIT_BARREL );
+ }
+ else {
+ m_calCollFlag.setBit( LCIO::CHBIT_BARREL );
+ }
+ }
+ // set for trk
+ else if ( g4sd->getType() == G4SensitiveDetector::eTracker ) {
+ if ( ec_flag ) {
+ m_trkCollFlag.unsetBit( LCIO::THBIT_BARREL );
+ }
+ else {
+ m_trkCollFlag.setBit( LCIO::THBIT_BARREL );
+ }
+ }
+ }
+
+ // save cal hits
+ void LcioHitsCollectionBuilder::saveHits(G4CalorimeterHitsCollection* calHits,
+ IMPL::LCCollectionVec* lcioColl)
+ {
+ size_t s = calHits->GetSize();
+ for ( size_t i = 0;
+ i < s;
+ i++ ) {
+ G4CalorimeterHit* calHit =
+ static_cast<G4CalorimeterHit*> ( calHits->GetHit(i) );
+ lcioColl->push_back( createHit( calHit ) );
+ }
+ }
+
+ // save trk hits
+ void LcioHitsCollectionBuilder::saveHits(G4TrackerHitsCollection* trkHits,
+ IMPL::LCCollectionVec* lcioColl)
+ {
+ size_t s = trkHits->GetSize();
+ for ( size_t i = 0;
+ i < s;
+ i++ ) {
+ G4TrackerHit* trkHit =
+ static_cast<G4TrackerHit*> ( trkHits->GetHit( i ) );
+ lcioColl->push_back( createHit( trkHit ) );
+ }
+ }
+
+ // create cal hit from G4
+ IMPL::SimCalorimeterHitImpl* LcioHitsCollectionBuilder::createHit(G4CalorimeterHit* calHit)
+ {
+ SimCalorimeterHitImpl* simCalHit = new SimCalorimeterHitImpl();
+
+ // set cellid from cal hit's id64
+ const Id64bit& id64 = calHit->getId64bit();
+ simCalHit->setCellID0( id64.getId0() );
+ simCalHit->setCellID1( id64.getId1() );
+
+ // position
+ const Hep3Vector hitPos = calHit->getPos();
+ float pos[3] = { hitPos.x(), hitPos.y(), hitPos.z() };
+ simCalHit->setPosition( pos );
+
+ // copy Mcp contrib info; energy is also incremented by contrib addition
+ addMcpContribs( calHit, simCalHit );
+
+ // compare edep of calHit with simHit when debugging
+#ifdef G4DEBUG
+ const McpHitContribList& contribs = calHit->getMcpHitContribList();
+ double totE = 0;
+ for ( McpHitContribList::const_iterator iter = contribs.begin();
+ iter != contribs.end();
+ iter++ ){
+ totE += (*iter).getEdep();
+ }
+
+ // sanity check so that new and old edeps must match
+ if (abs( totE/GeV - simCalHit->getEnergy() ) > ( 0.001 * totE ) ) {
+ std::cout << "g4 hit E: " << totE << std::endl;
+ std::cout << "sim hit E: " << simCalHit->getEnergy() << std::endl;
+ std::cout << "FATAL ERROR: LCIO simCalHit E != G4 CalHit E (within tolerance)" << std::endl;
+ assert( 0 );
+ }
+#endif
+
+ return simCalHit;
+ }
+
+ // create trk hit from G4
+ IMPL::SimTrackerHitImpl* LcioHitsCollectionBuilder::createHit(G4TrackerHit* trkHit)
+ {
+ SimTrackerHitImpl* simTrkHit = new SimTrackerHitImpl();
+
+ // position
+ const Hep3Vector hitPos = trkHit->getPos();
+ double pos[3] = { hitPos.x(), hitPos.y(), hitPos.z() };
+ simTrkHit->setPosition( pos );
+
+ /* momentum */
+ const G4ThreeVector& momentum = trkHit->getMomentum();
+ simTrkHit->setMomentum( momentum.x(), momentum.y(), momentum.z() );
+
+ // dEdx in GeV (LCIO units)
+ float edep = trkHit->getEdep();
+ simTrkHit->setdEdx( edep / GeV );
+
+ // time
+ float tEdep = trkHit->getTdep();
+ simTrkHit->setTime(tEdep);
+
+ // id
+ simTrkHit->setCellID( trkHit->getId() );
+
+ // MCP using McpManager
+ MCParticleImpl* mcp = m_mcpManager->getMaps()->findMcpFromTrackID( trkHit->getTrackID() );
+
+ if ( !mcp ) {
+ std::cout << "WARNING: no MCP found for trackID <" << trkHit->getTrackID() << "> for trk hit." << std::endl;
+ }
+ else {
+ simTrkHit->setMCParticle( mcp );
+ }
+
+ return simTrkHit;
+ }
+
+ // add an MCParticle hit contribution from G4 to LCIO
+ void LcioHitsCollectionBuilder::addMcpContribs( G4CalorimeterHit* g4CalHit, IMPL::SimCalorimeterHitImpl* simCalHit)
+ {
+ // create empty hit contrib list
+ McpHitContribList contribs;
+
+ // aggregation of contribs by track ID if CHBIT_PDG is not set
+ if ( !m_calCollFlag.bitSet(LCIO::CHBIT_PDG) ) {
+
+ // pass ref to contrib list, which will get filled
+ combineMcpHitContribs( g4CalHit->getMcpHitContribList(),
+ contribs );
+
+ }
+ // otherwise use the complete list from cal hit
+ else {
+
+ // set contrib list to cal contribs
+ contribs = g4CalHit->getMcpHitContribList();
+ }
+
+ // add contribs to LCIO Mcp
+ size_t ncontrib = 0;
+ for ( McpHitContribList::const_iterator iter = contribs.begin();
+ iter != contribs.end();
+ iter++ ) {
+
+ // get this contrib
+ const McpHitContrib& contrib = ( *iter );
+
+ // fetch associcated MCP pointer by trackID
+ MCParticleImpl* contribMcp = m_mcpManager->getMaps()->findMcpFromTrackID( contrib.getTrackID() );
+
+ if ( contribMcp != 0 ) {
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 1 ) {
+ std::cout << "found contribMcp <" << contribMcp
+ << "> for trackID <" << contrib.getTrackID()
+ << ">" << std::endl;
+ }
+#endif
+
+ // add the Mcp contrib to the Lcio cal hit
+ simCalHit->addMCParticleContribution(contribMcp,
+ contrib.getEdep() / GeV,
+ contrib.getGlobalTime(),
+ contrib.getPDGID()
+ );
+ ++ncontrib;
+ }
+#ifdef G4DEBUG
+ else {
+ std::cerr << "WARNING: No associated MCP found for trackID <" << contrib.getTrackID() << "> from cal hit." << std::endl;
+ }
+#endif
+ }
+
+#ifdef G4DEBUG
+ if ( ncontrib == 0 ) {
+ std::cerr << "WARNING: no hit contribs for sim cal hit." << std::endl;
+ }
+#endif
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 2 ) {
+ std::cout << "added num contribs <" << ncontrib << "> to cal hit" << std::endl;
+ }
+#endif
+ }
+
+ void LcioHitsCollectionBuilder::combineMcpHitContribs(const McpHitContribList& long_contrib,
+ McpHitContribList& combined_contrib)
+ {
+ combined_contrib.clear();
+
+ // iterate over long list (one entry for every hit)
+ for ( McpHitContribList::const_iterator iter = long_contrib.begin();
+ iter != long_contrib.end();
+ iter++ ) {
+
+ int trk_id = (*iter).getTrackID();
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 2 ) {
+ std::cout << "combining hits on trk_id: " << trk_id << std::endl;
+ }
+#endif
+
+ // old track id in new combined list?
+ McpHitContrib* trk_contrib = 0;
+ if ( ( trk_contrib = findMcpHitContribByTrackID( (*iter).getTrackID(), combined_contrib ) ) ) {
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 2 ) {
+ std::cout << "incr edep <" << (*iter).getEdep() << ">" << std::endl;
+ }
+#endif
+
+ // += edep
+ trk_contrib->incrEdep( (*iter).getEdep() );
+
+ // set min time
+ trk_contrib->setMinTime( (*iter).getGlobalTime() );
+ }
+ // no existing contrib
+ else {
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 2 ) {
+ std::cout << "new hit edep <" << (*iter).getEdep() << ">" << std::endl;
+ }
+#endif
+
+ // create new contrib
+ combined_contrib.push_back(McpHitContrib(trk_id,
+ (*iter).getEdep(),
+ 0,
+ (*iter).getGlobalTime() ) );
+ }
+ }
+ }
+
+ McpHitContrib* LcioHitsCollectionBuilder::findMcpHitContribByTrackID(int trk_id,
+ McpHitContribList& contribs)
+ {
+ McpHitContrib* c = 0;
+ for ( McpHitContribList::iterator iter = contribs.begin();
+ iter != contribs.end();
+ iter++ ) {
+ if ( (*iter).getTrackID() == trk_id ) {
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 2 ) {
+ std::cout << "found match on trkID <" << trk_id << ">" << std::endl;
+ }
+#endif
+
+ c = &(*iter);
+ break;
+ }
+ }
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() > 2 ) {
+ std::cout << "returning McpHitContrib ptr <" << c << ">" << std::endl;
+ }
+#endif
+
+ return c;
+ }
+
+ EVENT::LCEvent* LcioHitsCollectionBuilder::createHCsFromG4Event(const G4Event* g4evt, EVENT::LCEvent* lcevt)
+ {
+ // set instance vars
+ m_currentG4Event = g4evt;
+ m_currentLCEvent = lcevt;
+
+ // call real HC creation function
+ createHitCollections();
+
+ // return evt pntr, which is same as input
+ return m_currentLCEvent;
+ }
+
+ void LcioHitsCollectionBuilder::setLongFlag(bool setting)
+ {
+ if ( setting ) {
+ m_calCollFlag.setBit( LCIO::CHBIT_LONG );
+ }
+ else {
+ m_calCollFlag.unsetBit( LCIO::CHBIT_LONG );
+ }
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() ) {
+ std::cout << "Set CHBIT_LONG: " << setting << G4endl;
+ }
+#endif
+ }
+
+ void LcioHitsCollectionBuilder::setPDGFlag(bool setting)
+ {
+ if ( setting ) {
+ m_calCollFlag.setBit( LCIO::CHBIT_PDG );
+ }
+ else {
+ m_calCollFlag.setBit( LCIO::CHBIT_PDG );
+ }
+
+#ifdef SLIC_VERBOSE
+ if ( verbose() ) {
+ std::cout << "Set CHBIT_PDG: " << setting << G4endl;
+ }
+#endif
+ }
+
+};
slic/src
diff -N LcioHCBuilder.cc
--- LcioHCBuilder.cc 19 Sep 2005 17:33:49 -0000 1.31
+++ /dev/null 1 Jan 1970 00:00:00 -0000
@@ -1,507 +0,0 @@
-// $Header: /cvs/lcd/slic/src/Attic/LcioHCBuilder.cc,v 1.31 2005/09/19 17:33:49 jeremy Exp $
-#include "LcioHCBuilder.hh"
-
-// SLIC
-#include "LcioMcpManager.hh"
-#include "HitsCollectionUtil.hh"
-#include "G4Application.hh"
-#include "TrackerHitMomentum.hh"
-
-// LCIO
-#include "IMPL/LCCollectionVec.h"
-#include "IMPL/LCFlagImpl.h"
-#include "IMPL/MCParticleImpl.h"
-
-// G4
-#include "G4SDManager.hh"
-
-// std
-#include <iostream>
-
-using IMPL::SimCalorimeterHitImpl;
-using IMPL::SimTrackerHitImpl;
-using IMPL::MCParticleImpl;
-using IMPL::LCEventImpl;
-using IMPL::LCFlagImpl;
-using IMPL::LCCollectionVec;
-using EVENT::LCIO;
-
-using namespace std;
-
-namespace slic
-{
-
- LcioHCBuilder::LcioHCBuilder()
- : m_storeMomentum( 0 )
- {
- // set local LcioMcpManager ptr
- m_mcpManager = LcioMcpManager::instance();
-
- // setup default coll flag for cal hits
- setCalFlagDefaults();
- }
-
- LcioHCBuilder::~LcioHCBuilder()
- {}
-
- // create the hit collections
- void LcioHCBuilder::createHitCollections()
- {
- // fetch HCIDs
- std::vector<int> hcids = HitsCollectionUtil::getHCIDs();
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 1 ) {
- std::cout << "hcids size <" << hcids.size() << ">" << std::endl;
- }
-#endif
-
- // fetch hit collection of event
- G4HCofThisEvent* HCE = m_currentG4Event->GetHCofThisEvent();
-
- // HC table
- G4HCtable* HCtbl = G4SDManager::GetSDMpointer()->GetHCtable();
-
- // HCID
- G4int hcid;
-
- LCCollectionVec* collVect = 0;
-
- for (std::vector<int>::const_iterator iter = hcids.begin();
- iter != hcids.end();
- iter++) {
-
- hcid = *iter;
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 1 ) {
- std::cout << "HCID: " << hcid << std::endl;
- }
-#endif
-
- // retrieve Sensitive Detector ptr
- G4SensitiveDetector* SD =
- static_cast<G4SensitiveDetector*>
- ( G4SDManager::GetSDMpointer()->FindSensitiveDetector( HCtbl->GetSDname( hcid ) ) );
-
- // get hits collection
- G4VHitsCollection* HC = HCE->GetHC( hcid );
-
- // add a LCCollectionVec if got a HC
- if ( HC ) {
-
- // set LCIO endcap bit in the flag (i.e. CHBIT_BARREL ) according to SDs setting
- setEndcapFlag( SD );
-
- // create collection vector based on type of SD
- collVect = createCollectionVec( HC, SD->getType() );
-
- // add collection vector to LCEvent
- m_currentLCEvent->addCollection( collVect, HC->GetName() );
- }
- else {
- G4Exception( "LcioHCBuilder::createHitCollections() - No collection found for Hits Collection ID");
- }
- }
- }
-
- // create the CollectionVec (decides which overloaded subfunction to call)
- IMPL::LCCollectionVec* LcioHCBuilder::createCollectionVec(G4VHitsCollection* g4HC,
- G4SensitiveDetector::EType SDtype)
- {
- // vec to create
- LCCollectionVec* collVec = 0;
-
- // cal hits
- if ( SDtype == G4SensitiveDetector::eCalorimeter ) {
- collVec = createCalorimeterCollectionVec(g4HC);
- }
- // tracker hits
- else if ( SDtype == G4SensitiveDetector::eTracker ) {
- collVec = createTrackerCollectionVec(g4HC);
- }
- // unknown type of hit
- else {
- G4Exception( "Unknown HC type." );
- }
-
- return collVec;
- }
-
- LCCollectionVec* LcioHCBuilder::createTrackerCollectionVec(G4VHitsCollection* g4HC)
- {
- // create Lcio tracker coll
- LCCollectionVec* collVec = new LCCollectionVec( LCIO::SIMTRACKERHIT );
-
- // cast to G4 trk HC
- G4TrackerHitsCollection* trkHits =
- dynamic_cast<G4TrackerHitsCollection*> (g4HC);
-
- // call overloaded save function for trk hits
- saveHits(trkHits,
- collVec);
-
- // save the momentum in parallel vec
- if ( m_storeMomentum ) {
- TrackerHitMomentum::makeTrackerHitMomenta(trkHits, m_currentLCEvent);
- }
-
- // set trk flags
- collVec->setFlag( m_trkCollFlag.getFlag() );
-
- return collVec;
- }
-
- LCCollectionVec* LcioHCBuilder::createCalorimeterCollectionVec(G4VHitsCollection* g4HC)
- {
- // create Lcio cal coll
- LCCollectionVec* collVec = new LCCollectionVec(LCIO::SIMCALORIMETERHIT);
-
- // cast to G4 cal HC
- G4CalorimeterHitsCollection* calHits
- = dynamic_cast<G4CalorimeterHitsCollection*> (g4HC);
-
- // call overloaded save function for cal hits
- saveHits(calHits,
- collVec);
-
- // set cal flags
- collVec->setFlag( m_calCollFlag.getFlag() );
-
- return collVec;
- }
-
- void LcioHCBuilder::setCalFlagDefaults()
- {
- m_calCollFlag.setBit( LCIO::CHBIT_LONG );
- m_calCollFlag.setBit( LCIO::CHBIT_ID1 );
- }
-
- void LcioHCBuilder::setEndcapFlag(G4SensitiveDetector* g4sd)
- {
- bool ec_flag = g4sd->getEndcapFlag();
-
- // set for cal
- if ( g4sd->getType() == G4SensitiveDetector::eCalorimeter ) {
- if ( ec_flag ) {
- m_calCollFlag.unsetBit( LCIO::CHBIT_BARREL );
- }
- else {
- m_calCollFlag.setBit( LCIO::CHBIT_BARREL );
- }
- }
- // set for trk
- else if ( g4sd->getType() == G4SensitiveDetector::eTracker ) {
- if ( ec_flag ) {
- m_trkCollFlag.unsetBit( LCIO::THBIT_BARREL );
- }
- else {
- m_trkCollFlag.setBit( LCIO::THBIT_BARREL );
- }
- }
- }
-
- // save cal hits
- void LcioHCBuilder::saveHits(G4CalorimeterHitsCollection* calHits,
- IMPL::LCCollectionVec* lcioColl)
- {
- size_t s = calHits->GetSize();
- for ( size_t i = 0;
- i < s;
- i++ ) {
- G4CalorimeterHit* calHit =
- static_cast<G4CalorimeterHit*> ( calHits->GetHit(i) );
- lcioColl->push_back( createHit( calHit ) );
- }
- }
-
- // save trk hits
- void LcioHCBuilder::saveHits(G4TrackerHitsCollection* trkHits,
- IMPL::LCCollectionVec* lcioColl)
- {
- size_t s = trkHits->GetSize();
- for ( size_t i = 0;
- i < s;
- i++ ) {
- G4TrackerHit* trkHit =
- static_cast<G4TrackerHit*> ( trkHits->GetHit( i ) );
- lcioColl->push_back( createHit( trkHit ) );
- }
- }
-
- // create cal hit from G4
- IMPL::SimCalorimeterHitImpl* LcioHCBuilder::createHit(G4CalorimeterHit* calHit)
- {
- SimCalorimeterHitImpl* simCalHit = new SimCalorimeterHitImpl();
-
- // set cellid from cal hit's id64
- const Id64bit& id64 = calHit->getId64bit();
- simCalHit->setCellID0( id64.getId0() );
- simCalHit->setCellID1( id64.getId1() );
-
- // position
- const Hep3Vector hitPos = calHit->getPos();
- float pos[3] = { hitPos.x(), hitPos.y(), hitPos.z() };
- simCalHit->setPosition( pos );
-
- // copy Mcp contrib info; energy is also incremented by contrib addition
- addMcpContribs( calHit, simCalHit );
-
- // compare edep of calHit with simHit when debugging
-#ifdef G4DEBUG
- const McpHitContribList& contribs = calHit->getMcpHitContribList();
- double totE = 0;
- for ( McpHitContribList::const_iterator iter = contribs.begin();
- iter != contribs.end();
- iter++ ){
- totE += (*iter).getEdep();
- }
-
- // sanity check so that new and old edeps must match
- if (abs( totE/GeV - simCalHit->getEnergy() ) > ( 0.001 * totE ) ) {
- std::cout << "g4 hit E: " << totE << std::endl;
- std::cout << "sim hit E: " << simCalHit->getEnergy() << std::endl;
- std::cout << "FATAL ERROR: LCIO simCalHit E != G4 CalHit E (within tolerance)" << std::endl;
- assert( 0 );
- }
-#endif
-
- return simCalHit;
- }
-
- // create trk hit from G4
- IMPL::SimTrackerHitImpl* LcioHCBuilder::createHit(G4TrackerHit* trkHit)
- {
- SimTrackerHitImpl* simTrkHit = new SimTrackerHitImpl();
-
- // position
- const Hep3Vector hitPos = trkHit->getPos();
- double pos[3] = { hitPos.x(), hitPos.y(), hitPos.z() };
- simTrkHit->setPosition( pos );
-
- /* momentum */
- const G4ThreeVector& momentum = trkHit->getMomentum();
- simTrkHit->setMomentum( momentum.x(), momentum.y(), momentum.z() );
-
- // dEdx in GeV (LCIO units)
- float edep = trkHit->getEdep();
- simTrkHit->setdEdx( edep / GeV );
-
- // time
- float tEdep = trkHit->getTdep();
- simTrkHit->setTime(tEdep);
-
- // id
- simTrkHit->setCellID( trkHit->getId() );
-
- // MCP using McpManager
- MCParticleImpl* mcp = m_mcpManager->getMaps()->findMcpFromTrackID( trkHit->getTrackID() );
-
- if ( !mcp ) {
- std::cout << "WARNING: no MCP found for trackID <" << trkHit->getTrackID() << "> for trk hit." << std::endl;
- }
- else {
- simTrkHit->setMCParticle( mcp );
- }
-
- return simTrkHit;
- }
-
- // add an MCParticle hit contribution from G4 to LCIO
- void LcioHCBuilder::addMcpContribs( G4CalorimeterHit* g4CalHit, IMPL::SimCalorimeterHitImpl* simCalHit)
- {
- // create empty hit contrib list
- McpHitContribList contribs;
-
- // aggregation of contribs by track ID if CHBIT_PDG is not set
- if ( !m_calCollFlag.bitSet(LCIO::CHBIT_PDG) ) {
-
- // pass ref to contrib list, which will get filled
- combineMcpHitContribs( g4CalHit->getMcpHitContribList(),
- contribs );
-
- }
- // otherwise use the complete list from cal hit
- else {
-
- // set contrib list to cal contribs
- contribs = g4CalHit->getMcpHitContribList();
- }
-
- // add contribs to LCIO Mcp
- size_t ncontrib = 0;
- for ( McpHitContribList::const_iterator iter = contribs.begin();
- iter != contribs.end();
- iter++ ) {
-
- // get this contrib
- const McpHitContrib& contrib = ( *iter );
-
- // fetch associcated MCP pointer by trackID
- MCParticleImpl* contribMcp = m_mcpManager->getMaps()->findMcpFromTrackID( contrib.getTrackID() );
-
- if ( contribMcp != 0 ) {
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 1 ) {
- std::cout << "found contribMcp <" << contribMcp
- << "> for trackID <" << contrib.getTrackID()
- << ">" << std::endl;
- }
-#endif
-
- // add the Mcp contrib to the Lcio cal hit
- simCalHit->addMCParticleContribution(contribMcp,
- contrib.getEdep() / GeV,
- contrib.getGlobalTime(),
- contrib.getPDGID()
- );
- ++ncontrib;
- }
-#ifdef G4DEBUG
- else {
- std::cerr << "WARNING: No associated MCP found for trackID <" << contrib.getTrackID() << "> from cal hit." << std::endl;
- }
-#endif
- }
-
-#ifdef G4DEBUG
- if ( ncontrib == 0 ) {
- std::cerr << "WARNING: no hit contribs for sim cal hit." << std::endl;
- }
-#endif
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 2 ) {
- std::cout << "added num contribs <" << ncontrib << "> to cal hit" << std::endl;
- }
-#endif
- }
-
- void LcioHCBuilder::combineMcpHitContribs(const McpHitContribList& long_contrib,
- McpHitContribList& combined_contrib)
- {
- combined_contrib.clear();
-
- // iterate over long list (one entry for every hit)
- for ( McpHitContribList::const_iterator iter = long_contrib.begin();
- iter != long_contrib.end();
- iter++ ) {
-
- int trk_id = (*iter).getTrackID();
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 2 ) {
- std::cout << "combining hits on trk_id: " << trk_id << std::endl;
- }
-#endif
-
- // old track id in new combined list?
- McpHitContrib* trk_contrib = 0;
- if ( ( trk_contrib = findMcpHitContribByTrackID( (*iter).getTrackID(), combined_contrib ) ) ) {
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 2 ) {
- std::cout << "incr edep <" << (*iter).getEdep() << ">" << std::endl;
- }
-#endif
-
- // += edep
- trk_contrib->incrEdep( (*iter).getEdep() );
-
- // set min time
- trk_contrib->setMinTime( (*iter).getGlobalTime() );
- }
- // no existing contrib
- else {
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 2 ) {
- std::cout << "new hit edep <" << (*iter).getEdep() << ">" << std::endl;
- }
-#endif
-
- // create new contrib
- combined_contrib.push_back(McpHitContrib(trk_id,
- (*iter).getEdep(),
- 0,
- (*iter).getGlobalTime() ) );
- }
- }
- }
-
- McpHitContrib* LcioHCBuilder::findMcpHitContribByTrackID(int trk_id,
- McpHitContribList& contribs)
- {
- McpHitContrib* c = 0;
- for ( McpHitContribList::iterator iter = contribs.begin();
- iter != contribs.end();
- iter++ ) {
- if ( (*iter).getTrackID() == trk_id ) {
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 2 ) {
- std::cout << "found match on trkID <" << trk_id << ">" << std::endl;
- }
-#endif
-
- c = &(*iter);
- break;
- }
- }
-
-#ifdef SLIC_VERBOSE
- if ( verbose() > 2 ) {
- std::cout << "returning McpHitContrib ptr <" << c << ">" << std::endl;
- }
-#endif
-
- return c;
- }
-
- EVENT::LCEvent* LcioHCBuilder::createHCsFromG4Event(const G4Event* g4evt, EVENT::LCEvent* lcevt)
- {
- // set instance vars
- m_currentG4Event = g4evt;
- m_currentLCEvent = lcevt;
-
- // call real HC creation function
- createHitCollections();
-
- // return evt pntr, which is same as input
- return m_currentLCEvent;
- }
-
- void LcioHCBuilder::setLongFlag(bool setting)
- {
- if ( setting ) {
- m_calCollFlag.setBit( LCIO::CHBIT_LONG );
- }
- else {
- m_calCollFlag.unsetBit( LCIO::CHBIT_LONG );
- }
-
-#ifdef SLIC_VERBOSE
- if ( verbose() ) {
- std::cout << "Set CHBIT_LONG: " << setting << G4endl;
- }
-#endif
- }
-
- void LcioHCBuilder::setPDGFlag(bool setting)
- {
- if ( setting ) {
- m_calCollFlag.setBit( LCIO::CHBIT_PDG );
- }
- else {
- m_calCollFlag.setBit( LCIO::CHBIT_PDG );
- }
-
-#ifdef SLIC_VERBOSE
- if ( verbose() ) {
- std::cout << "Set CHBIT_PDG: " << setting << G4endl;
- }
-#endif
- }
-
-};
CVSspam 0.2.8