diff --git a/include/EventAction.hh b/include/EventAction.hh
index b01faafe01a21fe92e7861140174eb02c70ef0a1..370cd904d02f657e9ffc97787f780755a756c292 100644
--- a/include/EventAction.hh
+++ b/include/EventAction.hh
@@ -29,8 +29,6 @@
#include "globals.hh"
#include "G4UserEventAction.hh"
-#include <vector>
-
class EventAction : public G4UserEventAction{
public:
@@ -40,14 +38,6 @@ class EventAction : public G4UserEventAction{
virtual void BeginOfEventAction( const G4Event* event );
virtual void EndOfEventAction( const G4Event* event );
- inline std::vector< std::vector<double>* > GetVectors( ){return fPtrVec;}
-
- private:
- G4int hitsCollID;
- G4int fEventNo;
-
- std::vector<double> x,z,xHit,zHit,time,theta,energy;
- std::vector< std::vector<double>* > fPtrVec;
};
#endif
diff --git a/include/PMTSD.hh b/include/PMTSD.hh
index 4f07175a0a09fed6229b8b8ae9d77e506ad2f57c..a1499fb5078fdb07e3908c51ab53d1656e0f20d4 100644
--- a/include/PMTSD.hh
+++ b/include/PMTSD.hh
@@ -48,6 +48,9 @@ public:
private:
int HCID;
HitsCollection* hitCollection;
+
+ std::vector< std::vector<double>* > fPtrVec;
+ bool VISUALIZE;
};
#endif
diff --git a/include/RunAction.hh b/include/RunAction.hh
index 91f99ab809d280c96e4a9a3262dea431f5d2d43e..ef1c1e3d4f2e68500d283e1bfa501eac1dcd6b18 100644
--- a/include/RunAction.hh
+++ b/include/RunAction.hh
@@ -37,6 +37,8 @@
#include "globals.hh"
#include "G4UserRunAction.hh"
+#include <vector>
+
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
class G4Timer;
@@ -52,9 +54,13 @@ class RunAction : public G4UserRunAction
virtual void BeginOfRunAction(const G4Run* aRun);
virtual void EndOfRunAction(const G4Run* aRun);
+ std::vector< std::vector<double>* > GetVectors(){ return fPtrVec; }
+
private:
G4Timer* fTimer;
G4String m_fileName;
+
+ std::vector< std::vector<double>* > fPtrVec;
};
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/src/EventAction.cc b/src/EventAction.cc
index b07a4ad1c41f467b2d59d1f9f53722cb6a88b3ff..2db980b9f2ce9333f3f0601356c37ec40ab53f04 100644
--- a/src/EventAction.cc
+++ b/src/EventAction.cc
@@ -41,62 +41,32 @@
*
*/
EventAction::EventAction()
- : G4UserEventAction(),
- hitsCollID(0){
- fPtrVec.push_back(&x);
- fPtrVec.push_back(&z);
- fPtrVec.push_back(&xHit);
- fPtrVec.push_back(&zHit);
- fPtrVec.push_back(&time);
- fPtrVec.push_back(&theta);
- fPtrVec.push_back(&energy);
-
+ : G4UserEventAction()
+{
}
/*
*
*/
-EventAction::~EventAction(){
-
+EventAction::~EventAction()
+{
}
/*
*
*/
void EventAction::BeginOfEventAction(const G4Event* event){
- fEventNo = event->GetEventID();
- if(fEventNo%100000 == 0) G4cout << "Begin Event " << fEventNo << G4endl;
- hitsCollID = 0;
- for(uint i = 0; i < fPtrVec.size(); i++){
- fPtrVec.at(i)->clear();
- }
+ if(event->GetEventID()%100000 == 0) G4cout << "Begin Event " << event->GetEventID() << G4endl;
+
}
/*
*
*/
void EventAction::EndOfEventAction(const G4Event* event){
- hitsCollID = G4SDManager::GetSDMpointer()->GetCollectionID("MyPMT");
- G4HCofThisEvent* HCE = event->GetHCofThisEvent();
- if(HCE){
- HitsCollection* HC = (HitsCollection*)(HCE->GetHC(hitsCollID));
- if( HC->entries() ){
- for(size_t i = 0; i < HC->entries(); i++ ){
- PMTHit* aHit = (*HC)[i];
- // fill vectors //
- x.push_back ( aHit->getPos().x() );
- z.push_back ( aHit->getPos().z() );
- xHit.push_back ( aHit->getHit().x() );
- zHit.push_back ( aHit->getHit().z() );
- time.push_back ( aHit->getTime() );
- energy.push_back( aHit->getEnergy() );
- G4double pTheta = atan(sqrt(pow(aHit->getMomentum().x(),2) + pow(aHit->getMomentum().z(),2) )/fabs(aHit->getMomentum().y() ));
- theta.push_back( pTheta );
- }
- // fill ntuple //
- G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
- analysisManager->AddNtupleRow();
- }
- }
+ // fill ntuple //
+ G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
+ analysisManager->FillNtupleIColumn( 0, event->GetEventID());
+ analysisManager->AddNtupleRow();
}
diff --git a/src/PMTSD.cc b/src/PMTSD.cc
index 42c90ef238d1f3763898593d31a0b2012f9c943d..07186a7e0ae277c94f587316faa53d90657af23f 100644
--- a/src/PMTSD.cc
+++ b/src/PMTSD.cc
@@ -26,11 +26,14 @@
#include "PMTSD.hh"
#include "PMTHit.hh"
+#include "RunAction.hh"
#include "G4HCofThisEvent.hh"
#include "G4Step.hh"
#include "G4ThreeVector.hh"
#include "G4SDManager.hh"
+#include "G4RunManager.hh"
+#include "G4UImanager.hh"
#include "G4ios.hh"
#include "G4Poisson.hh"
#include "G4ParticleTypes.hh"
@@ -50,6 +53,7 @@ PMTSD::PMTSD(G4String sdName)
:G4VSensitiveDetector(sdName){
collectionName.insert(sdName);
HCID = -1;
+ VISUALIZE = false;
}
@@ -75,6 +79,14 @@ void PMTSD::Initialize(G4HCofThisEvent* HCE){
{ HCID = G4SDManager::GetSDMpointer()->GetCollectionID( name );}
HCE->AddHitsCollection( HCID, hitCollection );
+
+ // If the UI window exists, set visualize to true so we record and draw hits
+ if( G4UImanager::GetUIpointer()->GetG4UIWindow() ) VISUALIZE = true;
+
+ // Get the vectors from run action
+ RunAction* runAction = (RunAction*)G4RunManager::GetRunManager()->GetUserRunAction();
+ fPtrVec = runAction->GetVectors();
+
}
@@ -82,18 +94,29 @@ void PMTSD::Initialize(G4HCofThisEvent* HCE){
*
*/
G4bool PMTSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){
-
G4Track* theTrack = aStep->GetTrack();
- PMTHit* newHit = new PMTHit();
-
- newHit->setTrackID (theTrack->GetTrackID() );
- newHit->setPos (theTrack->GetVertexPosition() );
- newHit->setHit (theTrack->GetPosition() );
- newHit->setEnergy (aStep->GetPreStepPoint()->GetTotalEnergy());
- newHit->setMomentum (aStep->GetPreStepPoint()->GetMomentum());
- newHit->setTime (theTrack->GetGlobalTime());
- hitCollection->insert( newHit );
+ G4ThreeVector origin = theTrack->GetVertexPosition();
+ G4ThreeVector hitPos = theTrack->GetPosition();
+ G4ThreeVector momentum = aStep->GetPreStepPoint()->GetMomentum();
+
+ //Display hits if the UI is on
+ if(VISUALIZE){
+ PMTHit* newHit = new PMTHit();
+ newHit->setPos( theTrack->GetPosition() );
+ hitCollection->insert ( newHit );
+ }
+
+ // Fill ouput vectors
+ // These vectors are created in RunAction
+ fPtrVec[0]->push_back( origin.x() ); //Origin of the photon in x
+ fPtrVec[1]->push_back( origin.z() ); //Origin of the photon in z
+ fPtrVec[2]->push_back( hitPos.x() ); //Where the photon landed on the PMT in x
+ fPtrVec[3]->push_back( hitPos.z() ); //Where the photon landed on the PMT in z
+ fPtrVec[4]->push_back( theTrack->GetGlobalTime() ); //Time of arrival of the photon at the PMT
+ G4double pTheta = atan(sqrt(pow(momentum.x(),2) + pow(momentum.z(),2) )/fabs(momentum.y() ));
+ fPtrVec[5]->push_back( pTheta ); //Incident angle at the PMT about the y-axis
+ fPtrVec[6]->push_back( aStep->GetPreStepPoint()->GetTotalEnergy() ); //Photon energy
return true;
}
diff --git a/src/RunAction.cc b/src/RunAction.cc
index a01ef1634141baf2e173960438ff6a57c5f16f0d..aa87e78cb05b36d3b7bb931fd4d52af5962863b2 100644
--- a/src/RunAction.cc
+++ b/src/RunAction.cc
@@ -52,37 +52,6 @@ RunAction::RunAction(G4String fileName)
G4RunManager::GetRunManager()->SetPrintProgress(1);
- // Create analysis manager. The choice of analysis
- // technology is done via selection of a namespace
- // in B4Analysis.hh
- auto analysisManager = G4AnalysisManager::Instance();
- G4cout << "Using " << analysisManager->GetType() << G4endl;
-
- // Get event action
- const EventAction* constEventAction
- = static_cast<const EventAction*>(G4RunManager::GetRunManager()
- ->GetUserEventAction());
- EventAction* eventAction
- = const_cast<EventAction*>(constEventAction);
-
- std::vector< std::vector<double>* > pVec = eventAction->GetVectors( );
-
- // Create directories
- analysisManager->SetVerboseLevel(1);
- analysisManager->SetNtupleMerging(true);
-
- // Creating ntuple
- analysisManager->CreateNtuple("lightGuide", "pos and momentum");
- analysisManager->CreateNtupleDColumn( "X", *pVec[0] );
- analysisManager->CreateNtupleDColumn( "Z", *pVec[1] );
- analysisManager->CreateNtupleDColumn( "hitX", *pVec[2] );
- analysisManager->CreateNtupleDColumn( "hitZ", *pVec[3] );
- analysisManager->CreateNtupleDColumn( "time", *pVec[4] );
- analysisManager->CreateNtupleDColumn( "theta", *pVec[5] );
- analysisManager->CreateNtupleDColumn( "energy", *pVec[6] );
- analysisManager->CreateNtupleIColumn( "EventNo" );
- analysisManager->FinishNtuple();
-
}
/*
@@ -91,6 +60,7 @@ RunAction::RunAction(G4String fileName)
RunAction::~RunAction(){
delete fTimer;
delete G4AnalysisManager::Instance();
+ for(auto vec : fPtrVec ) delete vec;
}
/*
@@ -100,12 +70,35 @@ void RunAction::BeginOfRunAction(const G4Run* aRun){
G4cout << "### Run " << aRun->GetRunID() << " start." << G4endl;
fTimer->Start();
- // Get analysis manager
+ // Create analysis manager. The choice of analysis
+ // technology is done via selection of a namespace
+ // in B4Analysis.hh
auto analysisManager = G4AnalysisManager::Instance();
+ G4cout << "Using " << analysisManager->GetType() << G4endl;
// Open an output file
if(m_fileName == "") m_fileName = "output";
analysisManager->OpenFile(m_fileName);
+
+
+ // Create directories
+ analysisManager->SetVerboseLevel(1);
+ analysisManager->SetNtupleMerging(true);
+
+ // Creating ntuple
+ // Populate the vector with valid pointers and assign them to branches
+ // These vectors are filled in PMTSD
+ for(int i = 0; i < 7; i ++ ) fPtrVec.push_back( new std::vector< double > );
+ analysisManager->CreateNtuple("lightGuide", "pos and momentum");
+ analysisManager->CreateNtupleIColumn( "EventNo" );
+ analysisManager->CreateNtupleDColumn( "X", *fPtrVec[0] );
+ analysisManager->CreateNtupleDColumn( "Z", *fPtrVec[1] );
+ analysisManager->CreateNtupleDColumn( "hitX", *fPtrVec[2] );
+ analysisManager->CreateNtupleDColumn( "hitZ", *fPtrVec[3] );
+ analysisManager->CreateNtupleDColumn( "time", *fPtrVec[4] );
+ analysisManager->CreateNtupleDColumn( "theta", *fPtrVec[5] );
+ analysisManager->CreateNtupleDColumn( "energy", *fPtrVec[6] );
+ analysisManager->FinishNtuple();
}
/*