Commit 5bfa11ed authored by Chad Lantz's avatar Chad Lantz
Browse files

Added a test to see if optical photons are internally reflected after their first step

parent 2c824904
...@@ -40,6 +40,8 @@ ...@@ -40,6 +40,8 @@
#include "G4ParticleTypes.hh" #include "G4ParticleTypes.hh"
#include "G4ParticleDefinition.hh" #include "G4ParticleDefinition.hh"
#include "G4SystemOfUnits.hh" #include "G4SystemOfUnits.hh"
#include "G4OpBoundaryProcess.hh"
#include "G4OpProcessSubType.hh"
#include <string> #include <string>
#include <iostream> #include <iostream>
...@@ -110,7 +112,7 @@ void FiberSD::Initialize(G4HCofThisEvent* HCE){ ...@@ -110,7 +112,7 @@ void FiberSD::Initialize(G4HCofThisEvent* HCE){
G4bool FiberSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){ G4bool FiberSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){
G4ParticleDefinition *particle = aStep->GetTrack()->GetDefinition(); G4ParticleDefinition *particle = aStep->GetTrack()->GetDefinition();
//Get the number of Cherenkov photons created in this step //Get the number of Cherenkov photons created in this step
int capturedPhotons = 0; int capturedPhotons = 0;
const std::vector<const G4Track*>* secVec = aStep->GetSecondaryInCurrentStep(); const std::vector<const G4Track*>* secVec = aStep->GetSecondaryInCurrentStep();
...@@ -127,31 +129,67 @@ G4bool FiberSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){ ...@@ -127,31 +129,67 @@ G4bool FiberSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){
G4int rodNum = aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(0); G4int rodNum = aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(0);
G4ThreeVector pos = aStep->GetTrack()->GetPosition(); G4ThreeVector pos = aStep->GetTrack()->GetPosition();
G4ThreeVector momentum = aStep->GetPreStepPoint()->GetMomentum();
// G4ParticleDefinition *particle = aStep->GetTrack()->GetDefinition(); // G4ParticleDefinition *particle = aStep->GetTrack()->GetDefinition();
// If OPTICAL is true, determine if the photon has reached the top of the topOfVolume // If OPTICAL is true, check for photons that are totally internally reflected within the fiber
// and add the hit to the collection if it has
if(OPTICAL){ if(OPTICAL){
if( aStep->GetTrack()->GetDefinition() == G4OpticalPhoton::OpticalPhotonDefinition() && pos.y() >= m_topOfVolume - 0.1*mm){ static G4ParticleDefinition* opticalphoton =
if(REDUCED_TREE){ G4OpticalPhoton::OpticalPhotonDefinition();
m_cherenkovVec->at(rodNum)++; if( aStep->GetTrack()->GetDefinition() == opticalphoton && momentum.y() > 0.0 ){
FillTimeVector( rodNum, aStep->GetTrack()->GetGlobalTime() );
//Make sure we're on a boundary before checking for total internal reflection.
m_nHits++; //Otherwise let the simulation continue to a step that is on a boundary
return true; G4StepPoint* endPoint = aStep->GetPostStepPoint();
if(endPoint->GetStepStatus() != fGeomBoundary) return true;
//Get the refractive index and check if this photon will be reflected back down at the top of the fiber
G4MaterialPropertiesTable* MPT = aStep->GetPreStepPoint()->GetMaterial()->GetMaterialPropertiesTable();
//Make sure the PreStepPoint has a material defined, otherwise let the simulation continue
if(!MPT) return true;
G4double Rindex;
size_t index = 0;
//Make sure the material has a refractive index defined, otherwise let the simulation continue
G4MaterialPropertyVector* RindexMPV = MPT->GetProperty(kRINDEX);
if(!RindexMPV) return true;
Rindex = RindexMPV->Value(aStep->GetTrack()->GetDynamicParticle()->GetTotalMomentum(), index);
//Use the refractive index to see if the photon will be transmitted at the top of the fiber
//instead of reflected back downward
// if( 3.1415/2.0 - atan( momentum.y() / ( 1-momentum.y() ) ) >= asin(1.0/Rindex) ){
if( momentum.y() > 0.0 ){
//Now check for total internal reflection
G4ProcessVector* postStepDoItVector = opticalphoton->GetProcessManager()->GetPostStepProcessVector(typeDoIt);
G4int n_proc = postStepDoItVector->entries();
for(G4int i = 0; i < n_proc; ++i){
G4OpBoundaryProcess* opProc = dynamic_cast<G4OpBoundaryProcess*>((*postStepDoItVector)[i]);
if(opProc && opProc->GetStatus() == TotalInternalReflection){
if(REDUCED_TREE){
m_cherenkovVec->at(rodNum)++;
FillTimeVector( rodNum, aStep->GetTrack()->GetGlobalTime() );
m_nHits++;
return true;
}//end if REDUCED_TREE
FiberHit* newHit = new FiberHit();
newHit->setPos ( pos );
newHit->setOrigin ( aStep->GetTrack()->GetVertexPosition() );
newHit->setMomentum ( momentum );
newHit->setEnergy ( aStep->GetPreStepPoint()->GetTotalEnergy() );
newHit->setTime ( aStep->GetTrack()->GetGlobalTime() );
newHit->setRodNb ( rodNum );
fiberCollection->insert ( newHit );
m_nHits++;
}//end if opProc && TotalInternalReflection
}//end for n_proc
} }
FiberHit* newHit = new FiberHit();
newHit->setPos ( pos );
newHit->setOrigin ( aStep->GetTrack()->GetVertexPosition() );
newHit->setMomentum ( aStep->GetPreStepPoint()->GetMomentum() );
newHit->setEnergy ( aStep->GetPreStepPoint()->GetTotalEnergy() );
newHit->setTime ( aStep->GetTrack()->GetGlobalTime() );
newHit->setRodNb ( rodNum );
fiberCollection->insert ( newHit );
m_nHits++;
aStep->GetTrack()->SetTrackStatus( fStopAndKill ); //Kill the track so we only record it once //Kill all photons
aStep->GetTrack()->SetTrackStatus( fStopAndKill );
return true; return true;
} }
}else{ // Otherwise record all hits }else{ // Otherwise record all hits
...@@ -174,7 +212,7 @@ G4bool FiberSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){ ...@@ -174,7 +212,7 @@ G4bool FiberSD::ProcessHits(G4Step* aStep,G4TouchableHistory*){
newHit->setPos ( pos ); newHit->setPos ( pos );
newHit->setParticle ( particle ); newHit->setParticle ( particle );
newHit->setEnergy ( aStep->GetPreStepPoint()->GetTotalEnergy() ); newHit->setEnergy ( aStep->GetPreStepPoint()->GetTotalEnergy() );
newHit->setMomentum ( aStep->GetPreStepPoint()->GetMomentum() ); newHit->setMomentum ( momentum );
newHit->setTime ( aStep->GetTrack()->GetGlobalTime() ); newHit->setTime ( aStep->GetTrack()->GetGlobalTime() );
newHit->setNCherenkovs ( capturedPhotons ); newHit->setNCherenkovs ( capturedPhotons );
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment