From d0533d03d0c85f2f993f1793a6b9ea2af3391207 Mon Sep 17 00:00:00 2001 From: shivesh Date: Tue, 20 Nov 2018 17:26:02 +0000 Subject: Tue 20 Nov 17:26:02 GMT 2018 --- src/#G4TestEMPhysics.cxx# | 0 src/G4BeamTestDetectorConstruction.cxx | 3 +- src/G4BeamTestEMPhysics.cxx | 62 --------------- src/G4BeamTestEventAction.cxx | 91 ++++++++++++++++++++++ src/G4BeamTestPhysicsList.cxx | 32 ++++---- src/G4BeamTestRunManager.cxx | 131 +++++++++++++++++++++++++++++++ src/G4BeamTestSiHit.cxx | 116 ++++++++++++++++++++++++++++ src/G4BeamTestSiSD.cxx | 137 +++++++++++++++++++++++++++++++++ src/G4BeamTestTank.cxx | 60 ++++++++------- src/G4BeamTestUserStackingAction.cxx | 76 ++++++++++++++++++ src/G4BeamTestUserSteppingAction.cxx | 43 ++++++++++- src/G4BeamTestUserTrackingAction.cxx | 8 +- src/G4Interface.cxx | 8 ++ 13 files changed, 656 insertions(+), 111 deletions(-) create mode 100644 src/#G4TestEMPhysics.cxx# create mode 100644 src/G4BeamTestEventAction.cxx create mode 100644 src/G4BeamTestRunManager.cxx create mode 100644 src/G4BeamTestSiHit.cxx create mode 100644 src/G4BeamTestSiSD.cxx create mode 100644 src/G4BeamTestUserStackingAction.cxx (limited to 'src') diff --git a/src/#G4TestEMPhysics.cxx# b/src/#G4TestEMPhysics.cxx# new file mode 100644 index 0000000..e69de29 diff --git a/src/G4BeamTestDetectorConstruction.cxx b/src/G4BeamTestDetectorConstruction.cxx index 95a2071..ebbd9b8 100644 --- a/src/G4BeamTestDetectorConstruction.cxx +++ b/src/G4BeamTestDetectorConstruction.cxx @@ -60,7 +60,8 @@ G4VPhysicalVolume* G4BeamTestDetectorConstruction::Construct() // Do not create photons or electrons below cherenkov threshold // See also corresponding UserSpecialCuts in Physicslist !!!! G4UserLimits* energyLimit = new G4UserLimits(); - energyLimit->SetUserMinEkine(280.0 * CLHEP::keV); // Cherenkov threshold of electrons in ice TODO(shivesh) + /* energyLimit->SetUserMinEkine(280.0 * CLHEP::keV); // Cherenkov threshold of electrons in ice */ + energyLimit->SetUserMinEkine(1.907 * CLHEP::eV); // Lower threshold of PMT - 600nm worldLog->SetUserLimits(energyLimit); /* snowLog->SetUserLimits(energyLimit); */ diff --git a/src/G4BeamTestEMPhysics.cxx b/src/G4BeamTestEMPhysics.cxx index 1aa9a7c..55ecf99 100644 --- a/src/G4BeamTestEMPhysics.cxx +++ b/src/G4BeamTestEMPhysics.cxx @@ -55,66 +55,4 @@ void G4BeamTestEMPhysics::ConstructProcess() pManager->AddProcess(&positronIonisation, -1, 2, 2); pManager->AddProcess(&positronBremsStrahlung, -1, 3, 3); pManager->AddProcess(&annihilation, 0,-1, 4); - - // Cerenkov Physics - ConstructOp(); -} - -void G4BeamTestEMPhysics::ConstructOp() -{ - G4Cerenkov* cerenkovProcess = &cerenkov; - G4int fMaxNumPhotonStep = 20; - - cerenkovProcess->SetMaxNumPhotonsPerStep(fMaxNumPhotonStep); - G4cout << "MaxNumPhoton" << fMaxNumPhotonStep << G4endl; - - cerenkovProcess->SetMaxBetaChangePerStep(10.0); - cerenkovProcess->SetTrackSecondariesFirst(true); - G4Scintillation* scintillationProcess = new G4Scintillation("Scintillation"); - scintillationProcess->SetScintillationYieldFactor(1.); - scintillationProcess->SetTrackSecondariesFirst(true); - G4OpAbsorption* absorptionProcess = new G4OpAbsorption(); - G4OpRayleigh* rayleighScatteringProcess = new G4OpRayleigh(); - G4OpMieHG* mieHGScatteringProcess = new G4OpMieHG(); - G4OpBoundaryProcess* boundaryProcess = new G4OpBoundaryProcess(); - -/* cerenkovProcess->SetVerboseLevel(fVerboseLebel); */ -/* // scintillationProcess->SetVerboseLevel(fVerboseLebel); */ -/* absorptionProcess->SetVerboseLevel(fVerboseLebel); */ -/* rayleighScatteringProcess->SetVerboseLevel(fVerboseLebel); */ -/* mieHGScatteringProcess->SetVerboseLevel(fVerboseLebel); */ -/* boundaryProcess->SetVerboseLevel(fVerboseLebel); */ - - // Use Birks Correction in the Scintillation process - - if(!G4Threading::IsWorkerThread()) - { - G4EmSaturation* emSaturation = - G4LossTableManager::Instance()->EmSaturation(); - scintillationProcess->AddSaturation(emSaturation); - } - - auto theParticleIterator = GetParticleIterator(); - theParticleIterator->reset(); - while( (*theParticleIterator)() ){ - G4ParticleDefinition* particle = theParticleIterator->value(); - G4ProcessManager* pmanager = particle->GetProcessManager(); - G4String particleName = particle->GetParticleName(); - if (cerenkovProcess->IsApplicable(*particle)) { - pmanager->AddProcess(cerenkovProcess); - pmanager->SetProcessOrdering(cerenkovProcess,idxPostStep); - } - if (scintillationProcess->IsApplicable(*particle)) { - pmanager->AddProcess(scintillationProcess); - pmanager->SetProcessOrderingToLast(scintillationProcess, idxAtRest); - pmanager->SetProcessOrderingToLast(scintillationProcess, idxPostStep); - } - if (particleName == "opticalphoton") { - G4cout << " AddDiscreteProcess to OpticalPhoton " << G4endl; - pmanager->AddDiscreteProcess(absorptionProcess); - pmanager->AddDiscreteProcess(rayleighScatteringProcess); - pmanager->AddDiscreteProcess(mieHGScatteringProcess); - pmanager->AddDiscreteProcess(boundaryProcess); - } - } } diff --git a/src/G4BeamTestEventAction.cxx b/src/G4BeamTestEventAction.cxx new file mode 100644 index 0000000..176991d --- /dev/null +++ b/src/G4BeamTestEventAction.cxx @@ -0,0 +1,91 @@ +#include "G4Event.hh" + +#include "G4HCofThisEvent.hh" +#include "G4VHitsCollection.hh" +#include "G4SDManager.hh" +#include "G4ios.hh" +#include "G4VDigitizerModule.hh" +#include "G4DigiManager.hh" +#include "G4Event.hh" +#include "G4RunManager.hh" + +#include "G4BeamTestSiSD.h" +#include "G4BeamTestSiHit.h" +#include "G4BeamTestEventAction.h" + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +G4BeamTestEventAction::G4BeamTestEventAction() + : G4UserEventAction(), + + SiCollID(0) +{ +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestEventAction::~G4BeamTestEventAction() +{} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +void G4BeamTestEventAction::BeginOfEventAction(const G4Event* event ) +{ + + + G4SDManager * SDman = G4SDManager::GetSDMpointer(); + if(SiCollID<0) + { + G4String colNam; + SiCollID = SDman->GetCollectionID(colNam="G4BeamTestSiSDCollection"); + + } + + +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +void G4BeamTestEventAction::EndOfEventAction(const G4Event* event) +{ + + G4cout << ">>> Summary of Event " << event->GetEventID() << G4endl; + + if(SiCollID<0) return; + + G4HCofThisEvent* HCE = event->GetHCofThisEvent(); + G4BeamTestSiHitsCollection* SiHC = 0; + + if(HCE) + { + SiHC = (G4BeamTestSiHitsCollection*)(HCE->GetHC(SiCollID)); + + } + + + + + if(SiHC) + { + int n_hit = SiHC->entries(); + G4cout << G4endl; + G4cout << "Si hits " << + "--------------------------------------------------------------" + << G4endl; + G4cout << n_hit << " hits are stored in G4BeamTestSiHitsCollection." + << G4endl; + /* G4cout << "List of hits in tracker" << G4endl; */ + for(int i=0;iPrint(); */ + (*SiHC)[i]->Dataout(); + } + + // G4cout << "sid + " << SiCollID << G4endl; + + } + +} + + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + diff --git a/src/G4BeamTestPhysicsList.cxx b/src/G4BeamTestPhysicsList.cxx index 0c27b0b..65c575c 100644 --- a/src/G4BeamTestPhysicsList.cxx +++ b/src/G4BeamTestPhysicsList.cxx @@ -2,12 +2,12 @@ #include #include "G4BeamTestPhysicsList.h" #include "G4BeamTestGeneralPhysics.h" -#if G4VERSION_NUMBER < 1000 -#include "G4BeamTestEMPhysics.h" -#include "G4BeamTestMuonPhysics.h" -#include "G4BeamTestHadronPhysics.h" -#include "G4BeamTestIonPhysics.h" -#else +/* #if G4VERSION_NUMBER < 1000 */ +/* #include "G4BeamTestEMPhysics.h" */ +/* #include "G4BeamTestMuonPhysics.h" */ +/* #include "G4BeamTestHadronPhysics.h" */ +/* #include "G4BeamTestIonPhysics.h" */ +/* #else */ #include "G4DecayPhysics.hh" #include "G4EmStandardPhysics.hh" #include "G4EmExtraPhysics.hh" @@ -18,7 +18,9 @@ #include "G4DataQuestionaire.hh" #include "G4HadronPhysicsFTFP_BERT.hh" #include -#endif +/* #endif */ + +#include #include #include @@ -32,15 +34,17 @@ G4BeamTestPhysicsList::G4BeamTestPhysicsList() SetVerboseLevel(1); RegisterPhysics(new G4BeamTestGeneralPhysics); -#if G4VERSION_NUMBER < 1000 - RegisterPhysics(new G4BeamTestEMPhysics); - RegisterPhysics(new G4BeamTestMuonPhysics); - RegisterPhysics(new G4BeamTestHadronPhysics); - RegisterPhysics(new G4BeamTestIonPhysics); -#else +/* #if G4VERSION_NUMBER < 1000 */ +/* RegisterPhysics(new G4BeamTestEMPhysics); */ +/* RegisterPhysics(new G4BeamTestMuonPhysics); */ +/* RegisterPhysics(new G4BeamTestHadronPhysics); */ +/* RegisterPhysics(new G4BeamTestIonPhysics); */ +/* #else */ // The following is basically Geant4's FTFP_BERT physics list G4DataQuestionaire it(photon); // this checks that G4LEVELGAMMADATA is there + RegisterPhysics(new G4OpticalPhysics()); + RegisterPhysics(new G4EmStandardPhysics()); RegisterPhysics(new G4EmExtraPhysics()); RegisterPhysics(new G4DecayPhysics()); @@ -49,7 +53,7 @@ G4BeamTestPhysicsList::G4BeamTestPhysicsList() RegisterPhysics(new G4StoppingPhysics()); RegisterPhysics(new G4IonPhysics()); //RegisterPhysics(new G4NeutronTrackingCut()); -#endif +/* #endif */ } diff --git a/src/G4BeamTestRunManager.cxx b/src/G4BeamTestRunManager.cxx new file mode 100644 index 0000000..5ed2e05 --- /dev/null +++ b/src/G4BeamTestRunManager.cxx @@ -0,0 +1,131 @@ +// On Sun, to prevent conflict with ObjectSpace, G4Timer.hh has to be +// loaded *before* globals.hh... +#include "G4Timer.hh" + +#include "G4BeamTestRunManager.h" +#include +#include + + +G4BeamTestRunManager::G4BeamTestRunManager(): G4RunManager() +{ + +} + + +void G4BeamTestRunManager::BeamOn(G4int n_event,const char* macroFile,G4int n_select) +{ + G4String text = "BeamOn method is not supported in G4BeamTestRunManager!"; + G4Exception("G4BeamTestRunManager::BeamOn()", "G4BeamTestRunManager001", FatalException, text); +} + + +void G4BeamTestRunManager::InitializeRun() +{ + G4bool cond = ConfirmBeamOnCondition(); + if(cond) + { + // Reset the event counter + numberOfEventToBeProcessed = 0; + ConstructScoringWorlds(); + RunInitialization(); + + if(verboseLevel>0) timer->Start(); + } +} + + +void G4BeamTestRunManager::InjectParticle(G4ParticleGun* particleGun) +{ + if(!currentRun) + { + G4String text = "Run needs to be initialized before injecting a particle."; + G4Exception("G4BeamTestRunManager::InjectParticle()", "G4BeamTestRunManager002", FatalException, text); + } + assert(currentRun); // the G4Exception() above calls abort(). This assert() silences the clang static analyzer + + numberOfEventToBeProcessed++; + currentRun->SetNumberOfEventToBeProcessed(numberOfEventToBeProcessed); + + currentEvent = GenerateEvent(numberOfEventToBeProcessed); + particleGun->GeneratePrimaryVertex(currentEvent); + + eventManager->ProcessOneEvent(currentEvent); + AnalyzeEvent(currentEvent); + Update_Scoring(); + StackPreviousEvent(currentEvent); + currentEvent = 0; + + if(runAborted) TerminateRun(); +} + + +G4Event* G4BeamTestRunManager::GenerateEvent(G4int i_event) +{ + G4Event* anEvent = new G4Event(i_event); + + if(storeRandomNumberStatusToG4Event==1 || storeRandomNumberStatusToG4Event==3) + { + std::ostringstream oss; + CLHEP::HepRandom::saveFullState(oss); + randomNumberStatusForThisEvent = oss.str(); + anEvent->SetRandomNumberStatus(randomNumberStatusForThisEvent); + } + + if(storeRandomNumberStatus) + { + G4String fileN = randomNumberStatusDir + "currentEvent.rndm"; + CLHEP::HepRandom::saveEngineStatus(fileN); + } + + return anEvent; +} + + +void G4BeamTestRunManager::TerminateRun() +{ + if(verboseLevel>0) + { + timer->Stop(); + G4cout << "Run terminated." << G4endl; + G4cout << "Run Summary" << G4endl; + if(runAborted) + { + G4cout << " Run Aborted after " << numberOfEventToBeProcessed << " events processed." << G4endl; + } + else + { + G4cout << " Number of events processed : " << numberOfEventToBeProcessed << G4endl; + } + G4cout << " " << *timer << G4endl; + } + + RunTermination(); +} + + +// +// The following method is an exact copy of +// UpdateScoring which is private in the G4RunManager +// + +#include +#include +#include + +void G4BeamTestRunManager::Update_Scoring() +{ + G4ScoringManager* ScM = G4ScoringManager::GetScoringManagerIfExist(); + if(!ScM) return; + G4int nPar = ScM->GetNumberOfMesh(); + if(nPar<1) return; + + G4HCofThisEvent* HCE = currentEvent->GetHCofThisEvent(); + if(!HCE) return; + G4int nColl = HCE->GetCapacity(); + for(G4int i=0;iGetHC(i); + if(HC) ScM->Accumulate(HC); + } +} diff --git a/src/G4BeamTestSiHit.cxx b/src/G4BeamTestSiHit.cxx new file mode 100644 index 0000000..591274c --- /dev/null +++ b/src/G4BeamTestSiHit.cxx @@ -0,0 +1,116 @@ +#include +#include + +#include "G4BeamTestSiHit.h" +#include "G4UnitsTable.hh" +#include "G4VVisManager.hh" +#include "G4Circle.hh" +#include "G4Colour.hh" +#include "G4VisAttributes.hh" + +G4ThreadLocal G4Allocator* G4BeamTestSiHitAllocator=0; +static std::fstream testnew("./testnew.txt", std::ofstream::out); + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestSiHit::G4BeamTestSiHit() + : G4VHit(), + fTrackID(-1), + + fEdep(0.), + fTime(0.), + fPos(G4ThreeVector()) +{ + /* G4cout << "opening file" << G4endl; */ + /* testnew.open(testnew_out); */ +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestSiHit ::~G4BeamTestSiHit() { + /* G4cout << "closing file" << G4endl; */ + /* testnew.close(); */ +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestSiHit::G4BeamTestSiHit(const G4BeamTestSiHit& right) + : G4VHit() +{ + fTrackID = right.fTrackID; + // fChamberNb = right.fChamberNb; + fEdep = right.fEdep; + fPos = right.fPos; + fTime = right.fTime; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +const G4BeamTestSiHit& G4BeamTestSiHit::operator=(const G4BeamTestSiHit& right) +{ + fTrackID = right.fTrackID; + // fChamberNb = right.fChamberNb; + fEdep = right.fEdep; + fPos = right.fPos; + fTime = right.fTime; + + return *this; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4int G4BeamTestSiHit::operator==(const G4BeamTestSiHit& right) const +{ + return ( this == &right ) ? 1 : 0; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +void G4BeamTestSiHit::Draw() +{ + /* G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance(); */ + /* if(pVVisManager) */ + /* { */ + /* G4Circle circle(fPos); */ + /* circle.SetScreenSize(4.); */ + /* circle.SetFillStyle(G4Circle::filled); */ + /* G4Colour colour(1.,0.,0.); */ + /* G4VisAttributes attribs(colour); */ + /* circle.SetVisAttributes(attribs); */ + /* pVVisManager->Draw(circle); */ + /* } */ +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + + +void G4BeamTestSiHit::Print() +{ + G4cout + << " trackID: " << fTrackID + << "Edep: " + << std::setw(7) << G4BestUnit(fEdep,"Energy") + << " Position: " + << std::setw(7) << G4BestUnit(fPos ,"Length") + << "Time: " + << std::setw(7) << G4BestUnit(fTime,"Time") + << G4endl; +} + +void G4BeamTestSiHit::Dataout() +{ +testnew + << " trackID: " << fTrackID + << "Edep: " + << std::setw(7) << G4BestUnit(fEdep,"Energy") + << " Position: " + << std::setw(7) << G4BestUnit( fPos ,"Length") + << "Time: " + << std::setw(7) << G4BestUnit( fTime,"Time") + << G4endl; + + } + + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + diff --git a/src/G4BeamTestSiSD.cxx b/src/G4BeamTestSiSD.cxx new file mode 100644 index 0000000..b67092b --- /dev/null +++ b/src/G4BeamTestSiSD.cxx @@ -0,0 +1,137 @@ + +#include + +#include "G4UnitsTable.hh" + +#include "G4BeamTestSiSD.h" +#include "G4HCofThisEvent.hh" +#include "G4Step.hh" +#include "G4ThreeVector.hh" +#include "G4SDManager.hh" +#include "G4ios.hh" + +std::ofstream testdata("EVENTDATA.txt"); + + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestSiSD::G4BeamTestSiSD(const G4String& name, + const G4String& hitsCollectionName) + : G4VSensitiveDetector(name) + , fHitsCollection(NULL) +{ + collectionName.insert(hitsCollectionName); +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestSiSD::~G4BeamTestSiSD() +{ // testdata << "Arrival time" << " " << "Energy " << " " << "Distance" << std::endl; + +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +void G4BeamTestSiSD::Initialize(G4HCofThisEvent* hce) +{ + // Create hits collection + + testdata << "*" << std::endl; + + fHitsCollection + = new G4BeamTestSiHitsCollection(SensitiveDetectorName, collectionName[0]); + + // Add this collection in hce + + G4int hcID + = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]); + hce->AddHitsCollection( hcID, fHitsCollection ); +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4bool G4BeamTestSiSD::ProcessHits(G4Step* aStep, + G4TouchableHistory*) +{ + // G4double stepl = 10; + // aStep->SetStepLength(stepl); + + G4String name=aStep->GetTrack()->GetDefinition()->GetParticleName(); + + /* G4cout << " Particle_name = " << name << G4endl; */ + + if (name == "opticalphoton" || name == "gamma") { +/* G4cout << " Particle_name = " << name << G4endl; */ + +// energy deposit + G4double edep = aStep->GetTotalEnergyDeposit(); + +// if (edep==0.) return false; + + /* G4cout << " Particle_name_after_edep = " << name << G4endl; */ + + G4BeamTestSiHit* newHit = new G4BeamTestSiHit(); + // G4StepPoint* preStepPoint = aStep->GetPreStepPoint(); + + newHit->SetTrackID (aStep->GetTrack()->GetTrackID()); + newHit->SetEdep(edep); + newHit->SetPos (aStep->GetPostStepPoint()->GetPosition()); + newHit-> SetTime (aStep->GetPreStepPoint()->GetProperTime()); + fHitsCollection->insert( newHit ); + + + G4double zbef =aStep->GetPreStepPoint()->GetPosition().z(); + G4double zaft =aStep->GetPostStepPoint()->GetPosition().z(); + G4double zdelta = zaft-zbef; + //G4double zintial = aStep + + G4double zfix = (zdelta + (zbef- 4.7))*1000; + + // G4double zfix = (zdelta + (zbef- 375)); + + // G4cout << " Z VALUE = " << zfix << G4endl; + G4double time = aStep->GetPreStepPoint()->GetProperTime(); + G4double gtime = aStep->GetPreStepPoint()->GetGlobalTime(); + // G4cout << "Global Time " << gtime << G4endl; + G4double timefix = time - 0.288861775183071; + +/* G4cout << "TIME == " << time << G4endl; */ + + G4double vel = aStep->GetPreStepPoint()->GetVelocity(); + + G4double vel2 = zfix/time; + + // G4cout << G4BestUnit(vel, "Speed")<< G4endl; + +/* + + testdata << timefix << " " << edep << " " << zfix << std::endl; + + h_bragg9->Fill(zfix,edep); + enetime->Fill(timefix, edep); + // newHit->Print(); + sally->Fill(zfix,vel); + + G4cout << "DISTANCE = " << zfix << " Velocity = " << vel2 << " Velocityold = " << vel <1 ) { + G4int nofHits = fHitsCollection->entries(); + G4cout << G4endl + << "-------->Hits Collection: in this event they are " << nofHits + << " hits in the tracker chambers: " << G4endl; + for ( G4int i=0; iPrint(); + + } +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... diff --git a/src/G4BeamTestTank.cxx b/src/G4BeamTestTank.cxx index 8a5193d..56570a4 100644 --- a/src/G4BeamTestTank.cxx +++ b/src/G4BeamTestTank.cxx @@ -14,7 +14,8 @@ /* #include */ #include "G4BeamTestTank.h" -#include "G4TankIceSD.h" +#include "G4BeamTestSiSD.h" +/* #include "G4TankIceSD.h" */ //prevent gcc to make something stupid with pretended unused variables #ifdef __GNUC__ @@ -55,7 +56,7 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const // See also corresponding UserSpecialCuts in Physicslist !!!! // TODO(shivesh): Maybe do all of this as stepping action ?????? G4UserLimits* energyLimit = new G4UserLimits(); - energyLimit->SetUserMinEkine(264.1 * CLHEP::keV); // Cherenkov threshold of electrons in water + energyLimit->SetUserMinEkine(1.907 * CLHEP::eV); // Lower threshold of PMT - 600nm // std::string tankName=boost::lexical_cast(tankKey_); std::string tankName = "BTT"; @@ -139,7 +140,7 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const new G4LogicalVolume(upperglasssphere, glass, ("log_dom_up_" + omName).c_str(), 0, 0, 0); G4LogicalVolume* logLowerGlass = - new G4LogicalVolume(lowerglasssphere, glass, + new G4LogicalVolume(lowerglasssphere, water, ("log_dom_lo_" + omName).c_str(), 0, 0, 0); G4LogicalVolume* logUpperDOM = new G4LogicalVolume(upperdomsphere, effectiveDOM, @@ -167,11 +168,12 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const logLowerDOM->SetUserLimits(energyLimit); // } - // Define sensitive detector TODO(shivesh): make the PMT the SD + // Define sensitive detector G4SDManager* sdManager = G4SDManager::GetSDMpointer(); - iceSD_ = new G4TankIceSD(("ice_SD_" + tankName).c_str()); + iceSD_ = new G4BeamTestSiSD(("ice_SD_" + tankName).c_str(), "HitsCollection"); sdManager->AddNewDetector(iceSD_); - logWater->SetSensitiveDetector(iceSD_); + /* logLowerDOM->SetSensitiveDetector(iceSD_); */ + logLowerGlass->SetSensitiveDetector(iceSD_); // Instantiation of a set of visualization attributes with red colour G4VisAttributes * tankVisAtt = new G4VisAttributes(G4Colour(1,0,0)); @@ -196,29 +198,29 @@ G4VPhysicalVolume* G4BeamTestTank::InstallTank(G4VPhysicalVolume* mother, const } -double G4BeamTestTank::GetNumCherenkov(/* const OMKey& omKey */) -{ - return std::max(iceSD_->GetNumCherenkov(/* omKey */), 0.); -} - - -double G4BeamTestTank::GetNumCherenkovWeight(/* const OMKey& omKey */) -{ - return std::max(iceSD_->GetNumCherenkovWeight(/* omKey */), 0.); -} - - -double G4BeamTestTank::GetEDep_G4(/* const OMKey& omKey */) -{ - return std::max(iceSD_->GetEDep(/* omKey */), 0.); -} - - -double G4BeamTestTank::GetTime_G4(/* const OMKey& omKey */) -{ - return iceSD_->GetTime(/* omKey */); -} - +/* double G4BeamTestTank::GetNumCherenkov(#<{(| const OMKey& omKey |)}>#) */ +/* { */ +/* return std::max(iceSD_->GetNumCherenkov(#<{(| omKey |)}>#), 0.); */ +/* } */ +/* */ +/* */ +/* double G4BeamTestTank::GetNumCherenkovWeight(#<{(| const OMKey& omKey |)}>#) */ +/* { */ +/* return std::max(iceSD_->GetNumCherenkovWeight(#<{(| omKey |)}>#), 0.); */ +/* } */ +/* */ +/* */ +/* double G4BeamTestTank::GetEDep_G4(#<{(| const OMKey& omKey |)}>#) */ +/* { */ +/* return std::max(iceSD_->GetEDep(#<{(| omKey |)}>#), 0.); */ +/* } */ +/* */ +/* */ +/* double G4BeamTestTank::GetTime_G4(#<{(| const OMKey& omKey |)}>#) */ +/* { */ +/* return iceSD_->GetTime(#<{(| omKey |)}>#); */ +/* } */ +/* */ /* double G4BeamTestTank::GetEDep_I3(const OMKey& omKey) */ /* { */ diff --git a/src/G4BeamTestUserStackingAction.cxx b/src/G4BeamTestUserStackingAction.cxx new file mode 100644 index 0000000..73a354c --- /dev/null +++ b/src/G4BeamTestUserStackingAction.cxx @@ -0,0 +1,76 @@ +#include "G4BeamTestUserStackingAction.h" + +#include "G4VProcess.hh" + +#include "G4ParticleDefinition.hh" +#include "G4ParticleTypes.hh" +#include "G4Track.hh" +#include "G4ios.hh" + +#include "G4Event.hh" + +#include "G4HCofThisEvent.hh" +#include "G4VHitsCollection.hh" +#include "G4SDManager.hh" +#include "G4ios.hh" +#include "G4VDigitizerModule.hh" +#include "G4DigiManager.hh" +#include "G4Event.hh" +#include "G4RunManager.hh" +//#include "TH1F.h" +//#include "TH2F.h" + + + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestUserStackingAction::G4BeamTestUserStackingAction() + : G4UserStackingAction(), + fScintillationCounter(0), + fCerenkovCounter(0) +{} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4BeamTestUserStackingAction::~G4BeamTestUserStackingAction() +{} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +G4ClassificationOfNewTrack +G4BeamTestUserStackingAction::ClassifyNewTrack(const G4Track * aTrack) +{ + if(aTrack->GetDefinition() == G4OpticalPhoton::OpticalPhotonDefinition()) + { // particle is optical photon + if(aTrack->GetParentID()>0) + { // particle is secondary + if(aTrack->GetCreatorProcess()->GetProcessName() == "Scintillation") + fScintillationCounter++; + if(aTrack->GetCreatorProcess()->GetProcessName() == "Cerenkov") + fCerenkovCounter++; + } + } + return fUrgent; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +void G4BeamTestUserStackingAction::NewStage() +{ + G4cout << "Number of Scintillation photons produced in this event : " + << fScintillationCounter << G4endl; + G4cout << "Number of Cerenkov photons produced in this event : " + << fCerenkovCounter << G4endl; +} + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + +void G4BeamTestUserStackingAction::PrepareNewEvent() +{ + fScintillationCounter = 0; + fCerenkovCounter = 0; +} + + +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... + diff --git a/src/G4BeamTestUserSteppingAction.cxx b/src/G4BeamTestUserSteppingAction.cxx index e064628..f81ab62 100644 --- a/src/G4BeamTestUserSteppingAction.cxx +++ b/src/G4BeamTestUserSteppingAction.cxx @@ -5,11 +5,26 @@ #include "G4LogicalVolume.hh" #include "G4UserLimits.hh" #include "G4SteppingManager.hh" +#include "G4OpticalPhoton.hh" +#include "G4RunManager.hh" G4BeamTestUserSteppingAction::G4BeamTestUserSteppingAction(){} -void G4BeamTestUserSteppingAction::UserSteppingAction(const G4Step*) +void G4BeamTestUserSteppingAction::UserSteppingAction(const G4Step* step) { + G4int eventNumber = G4RunManager::GetRunManager()-> + GetCurrentEvent()->GetEventID(); + + if (eventNumber != fEventNumber) { + /* G4cout << " Number of Scintillation Photons in previous event: " */ + /* << fScintillationCounter << G4endl; */ + /* G4cout << " Number of Cerenkov Photons in previous event: " */ + /* << fCerenkovCounter << G4endl; */ + fEventNumber = eventNumber; + fScintillationCounter = 0; + fCerenkovCounter = 0; + } + G4Track* track = fpSteppingManager->GetTrack(); if(track) { @@ -18,11 +33,33 @@ void G4BeamTestUserSteppingAction::UserSteppingAction(const G4Step*) if(!limit) G4cout << "----> G4LogicalVolume: " << volume->GetName() << " has no defined G4UserLimit" << G4endl; G4double threshold = limit->GetUserMinEkine(*track); //check if particle is a gamma - if(track->GetDefinition()->GetParticleName() == "gamma") + G4String particle = track->GetDefinition()->GetParticleName(); + /* G4cout << "particle = " << particle << G4endl; */ + if(particle == "gamma") { //check if particle energy is below threshold; if true, kill the particle G4double energy = track->GetTotalEnergy(); - if(energy < threshold) track->SetTrackStatus(fStopAndKill); + if(energy < threshold){ + G4cout << "SteppingAction: killing particle " << particle << " with energy " << energy << " < " << threshold << G4endl; + track->SetTrackStatus(fStopAndKill); + } } } + + const std::vector* secondaries = + step->GetSecondaryInCurrentStep(); + + if (secondaries->size()>0) { + for(unsigned int i=0; isize(); ++i) { + if (secondaries->at(i)->GetParentID()>0) { + if(secondaries->at(i)->GetDynamicParticle()->GetParticleDefinition() + == G4OpticalPhoton::OpticalPhotonDefinition()){ + if (secondaries->at(i)->GetCreatorProcess()->GetProcessName() + == "Scintillation")fScintillationCounter++; + if (secondaries->at(i)->GetCreatorProcess()->GetProcessName() + == "Cerenkov")fCerenkovCounter++; + } + } + } + } } diff --git a/src/G4BeamTestUserTrackingAction.cxx b/src/G4BeamTestUserTrackingAction.cxx index c6317df..f7c053b 100644 --- a/src/G4BeamTestUserTrackingAction.cxx +++ b/src/G4BeamTestUserTrackingAction.cxx @@ -24,11 +24,15 @@ void G4BeamTestUserTrackingAction::PostUserTrackingAction(const G4Track* track) for(size_t i=0;iGetDefinition()->GetParticleName() == "gamma") + G4String particle = (*secondaries)[i]->GetDefinition()->GetParticleName(); + if(particle == "gamma" || particle == "opticalphoton") { //check if particle energy is below threshold; if true, kill the particle G4double energy = (*secondaries)[i]->GetTotalEnergy(); - if(energy < threshold) (*secondaries)[i]->SetTrackStatus(fStopAndKill); + if(energy < threshold){ + G4cout << "TrackingAction: killing particle " << particle << " with energy " << energy << " < " << threshold << G4endl; + (*secondaries)[i]->SetTrackStatus(fStopAndKill); + } } } } diff --git a/src/G4Interface.cxx b/src/G4Interface.cxx index 75d6f78..e552173 100644 --- a/src/G4Interface.cxx +++ b/src/G4Interface.cxx @@ -5,7 +5,9 @@ #include "G4BeamTestPhysicsList.h" #include "G4BeamTestUserTrackingAction.h" #include "G4BeamTestUserSteppingAction.h" +#include "G4BeamTestUserStackingAction.h" #include "G4BeamTestPrimaryGeneratorAction.h" +#include "G4BeamTestEventAction.h" #include "G4BeamTestRunAction.h" /* #include */ @@ -165,6 +167,12 @@ void G4Interface::Initialize() G4cout << "Init UserSteppingAction ..." << G4endl; runManager_.SetUserAction(new G4BeamTestUserSteppingAction()); + G4cout << "Init UserStackingAction ..." << G4endl; + runManager_.SetUserAction(new G4BeamTestUserStackingAction()); + + G4cout << "Init G4BeamTestEventAction ..." << G4endl; + runManager_.SetUserAction(new G4BeamTestEventAction()); + // Initialize G4 kernel /* log_debug("Init run manager ..."); */ G4cout << "Init run manager ..." << G4endl; -- cgit v1.2.3