286 lines
9.4 KiB
C++
286 lines
9.4 KiB
C++
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$//*
|
|
// LOW ENERGY MUON SPIN RELAXATION, ROTATION, RADIATION Geant4 SIMULATION
|
|
// ID : LEMuSRSteppingAction.cc , v 1.0
|
|
// AUTHOR: Taofiq PARAISO
|
|
// DATE : 2004-07-07 11:15
|
|
//
|
|
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$//
|
|
//
|
|
// & &&&&&&&&&& &&&&&&& &&&&&&&&
|
|
// & & && && & &&
|
|
// & & & & & & &&
|
|
// & &&&&&&& & & &&&&&& &&&&&&&&
|
|
// & & & && & & &&
|
|
// & & && & & && && & &
|
|
// &&&&&&&&&& &&&&&&&&&& & &&&&& && &&&&&&& & &&
|
|
// &
|
|
// &
|
|
// &
|
|
// &
|
|
// STEPPING ACTION.CC
|
|
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$//
|
|
|
|
|
|
|
|
|
|
#include "LEMuSRSteppingAction.hh"
|
|
#include "G4SteppingManager.hh"
|
|
#include "G4Transform3D.hh"
|
|
#include "G4DynamicParticle.hh"
|
|
#include "G4ParticleDefinition.hh"
|
|
#include "G4ProcessManager.hh"
|
|
#include "LEMuSRDecay.hh"
|
|
#include "G4VVisManager.hh"
|
|
#include "G4Polyline.hh"
|
|
#include "G4VisAttributes.hh"
|
|
#include "G4Colour.hh"
|
|
|
|
|
|
//! \ct
|
|
LEMuSRSteppingAction::LEMuSRSteppingAction()
|
|
{
|
|
pointer=this ;
|
|
loop=0;
|
|
}
|
|
|
|
//! \dt
|
|
LEMuSRSteppingAction::~LEMuSRSteppingAction()
|
|
{
|
|
;
|
|
}
|
|
|
|
LEMuSRSteppingAction* LEMuSRSteppingAction::pointer=0;
|
|
LEMuSRSteppingAction* LEMuSRSteppingAction::GetInstance()
|
|
{
|
|
return pointer;
|
|
}
|
|
|
|
|
|
//! \mm
|
|
/*!
|
|
* The main role of the stepping action in \lemu simulation is to kill particles that are looping in electromagnetic fields.
|
|
*
|
|
* Is also illustrated here the possiblity of personalizing the color of the trajectories for visualisation.
|
|
*/
|
|
void LEMuSRSteppingAction::UserSteppingAction(const G4Step* aStep)
|
|
|
|
{
|
|
LoopKiller(aStep);
|
|
|
|
G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance();
|
|
G4String p_name;
|
|
p_name = aStep->GetTrack()->GetDefinition()->GetParticleName(); // particle name
|
|
if (pVVisManager) {
|
|
//----- Define a line segment
|
|
G4Polyline polyline;
|
|
G4Colour colour;
|
|
if (p_name == "mu+") colour = G4Colour(1., 0., 0.);
|
|
else if (p_name == "Mu" ) colour = G4Colour(0., 0., 1.);
|
|
else if (p_name == "e+" ) colour = G4Colour(1., 1., 0.);
|
|
else colour = G4Colour(1., 0., 1.);
|
|
G4VisAttributes attribs(colour);
|
|
polyline.SetVisAttributes(attribs);
|
|
polyline.push_back(aStep->GetPreStepPoint()->GetPosition());
|
|
polyline.push_back(aStep->GetPostStepPoint()->GetPosition());
|
|
|
|
//----- Call a drawing method for G4Polyline
|
|
pVVisManager -> Draw(polyline);
|
|
|
|
}
|
|
}
|
|
|
|
|
|
//! Method to kill the looping partiles
|
|
/*!
|
|
* A looping particle can increase the CPU time, and in the worst case even paralize the simulation. To kill the loops, we just limit the steps number of one track to 2500. The loop killer method check the step number and kills the particle when it gets too big. Note that a usual track does not contain more that few hundred steps and experience showed that above 1000 steps, the particles were actually trapped in a loop and rarely recovered a normal trajectory.
|
|
*
|
|
* In order to verify that this method is not abusive, the word "killed" is printed to the screen each time a particle is killed. Usually, this is the case of few particles over hundred thousands.
|
|
*
|
|
* The neutrinos, gammas and electrons tracks are killed directly because they have no influence on the simulation but increasing the calculation time.
|
|
* In particular, electrons are important sources of looping tracks in electromagntic field.
|
|
*
|
|
* Of course, one can select the particles to kill or not.
|
|
*/
|
|
|
|
void LEMuSRSteppingAction::LoopKiller(const G4Step *aStep)
|
|
{
|
|
// loop killa
|
|
if(aStep->GetTrack()->GetCurrentStepNumber()>2500)
|
|
{
|
|
aStep->GetTrack()->SetTrackStatus(fStopAndKill);
|
|
G4cout<<"killed "<<G4endl;
|
|
}
|
|
|
|
/* if(aStep->GetTrack()->GetDefinition()->GetParticleName()=="e-"
|
|
||aStep->GetTrack()->GetDefinition()->GetParticleName()=="gamma"
|
|
||aStep->GetTrack()->GetDefinition()->GetParticleName()=="nu_e"
|
|
||aStep->GetTrack()->GetDefinition()->GetParticleName()=="anti_nu_e"
|
|
||aStep->GetTrack()->GetDefinition()->GetParticleName()=="nu_mu"
|
|
||aStep->GetTrack()->GetDefinition()->GetParticleName()=="anti_nu_mu")
|
|
{
|
|
aStep->GetTrack()->SetTrackStatus(fStopAndKill);
|
|
}
|
|
*/
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
//--------------------------------------------------------------------------------
|
|
//--------------------------------------------------------------------------------
|
|
|
|
|
|
|
|
|
|
|
|
void LEMuSRSteppingAction:: FieldInfo(const G4Step *aStep)
|
|
{
|
|
|
|
G4String p_name;
|
|
p_name = aStep->GetTrack()->GetDefinition()->GetParticleName(); // particle name
|
|
G4String v_name = aStep->GetTrack()->GetVolume()->GetLogicalVolume()->GetName();
|
|
|
|
if( aStep->GetTrack()->GetVolume()->GetLogicalVolume()->GetFieldManager())
|
|
{
|
|
|
|
G4cout << "Field Manager...OK: \n"
|
|
<< "particle name : " << p_name <<" ;\n "
|
|
<< "volume name : " << v_name <<" ;\n "
|
|
<< "position : " << aStep->GetPostStepPoint()->GetPosition() <<" ;\n "
|
|
<< "kenergy : " << aStep->GetTrack()->GetKineticEnergy()<<" ;\n "
|
|
<< "momentum direction: " << aStep->GetPostStepPoint()->GetMomentumDirection() <<" ;\n "
|
|
<< "polarization : " << aStep->GetPostStepPoint()->GetPolarization() <<" ;\n "
|
|
<< "time : " << aStep->GetPostStepPoint()->GetGlobalTime() <<" ;\n " ;
|
|
// field datas
|
|
// G4cout << "Field Step : " << aStep->GetTrack()->GetVolume()->GetLogicalVolume()->GetFieldManager()->GetDeltaOneStep() <<" ;\n " ;
|
|
double point[4];
|
|
point[0]= aStep->GetPostStepPoint()->GetPosition().x()/mm/100;
|
|
point[1]= aStep->GetPostStepPoint()->GetPosition().y()/mm/100;
|
|
point[2]= aStep->GetPostStepPoint()->GetPosition().z()/mm/100 + 5.67;
|
|
double field[6];
|
|
aStep->GetTrack()->GetVolume()->GetLogicalVolume()->GetFieldManager()->GetDetectorField()->GetFieldValue(point,field);
|
|
G4cout <<"pOSITION : "<< point[0] <<" "<<point[1]<<" "<<point[2] <<" ;\n "
|
|
<< "Field direction : " << field[0] <<" "<<field[1]<<" "<<field[2] <<" ;\n " ;
|
|
|
|
|
|
|
|
}
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
void LEMuSRSteppingAction:: ParticleInfo(const G4Step *aStep)
|
|
{
|
|
G4String p_name;
|
|
p_name = aStep->GetTrack()->GetDefinition()->GetParticleName(); // particle name
|
|
G4String v_name = aStep->GetTrack()->GetVolume()->GetLogicalVolume()->GetName();
|
|
|
|
// Get datas
|
|
//a
|
|
if(aStep->GetTrack()->GetDefinition())
|
|
{
|
|
p_name = aStep->GetTrack()->GetDynamicParticle()->GetDefinition()->GetParticleName(); // particle name
|
|
|
|
// if(p_name == "mu+")
|
|
// {
|
|
G4String nv_name;
|
|
G4String v_name = aStep->GetTrack()->GetVolume()->GetName(); // particle name
|
|
if(v_name!="pv_World")
|
|
{
|
|
nv_name = aStep->GetTrack()->GetNextVolume()->GetName(); // particle name
|
|
}
|
|
else if(v_name=="pv_World")
|
|
{
|
|
nv_name=v_name;
|
|
}
|
|
|
|
G4double spin= aStep->GetTrack()->GetDynamicParticle()->GetDefinition()->GetPDGSpin(); // spin in units of 1
|
|
//b
|
|
G4ThreeVector position = aStep->GetPostStepPoint()->GetPosition(); // position
|
|
G4ThreeVector momentum = aStep->GetPostStepPoint()->GetMomentum(); // momentum
|
|
//c
|
|
G4double tof = aStep->GetTrack()->GetLocalTime(); // time since track creation
|
|
G4double globaltime = aStep->GetTrack()->GetGlobalTime();// time since the event in which the track belongs is created.
|
|
G4double proptime = aStep->GetTrack()->GetProperTime(); // proper time of the particle
|
|
//d
|
|
G4double edep = aStep->GetTotalEnergyDeposit();
|
|
|
|
// Print Datas
|
|
|
|
if( aStep->GetTrack()->GetCreatorProcess())
|
|
{
|
|
G4cout << "NOT PRIMARY PARTICLE : created by : " << aStep->GetTrack()->GetCreatorProcess()->GetProcessName() <<" ;\n ";
|
|
// G4cout << "Parent particles are : " << aStep->GetTrack()->GetCreatorProcess()->GetProcessName() <<" ;\n ";
|
|
|
|
}
|
|
|
|
G4cout << "particle name : " << p_name <<" ;\n "
|
|
<< "volume name : " << v_name <<" ;\n "
|
|
<< "next volume name : " << nv_name <<" ;\n "
|
|
<< "spin : " << spin <<" ;\n "
|
|
<< "current energy : " << aStep->GetTrack()->GetDynamicParticle()->GetKineticEnergy()/keV <<" keV;\n "
|
|
<< "energy_deposition : " << G4BestUnit(edep,"Energy") <<" ;\n "
|
|
<< "time_of_flight : " << G4BestUnit(tof,"Time") <<" ; " << G4BestUnit(globaltime,"Time")<<" ; " << G4BestUnit(proptime,"Time") <<" ;\n "
|
|
<< "position : " << position <<" ;\n "
|
|
<< "momentum : " << momentum <<" ;\n "
|
|
<< "polarization : " << aStep->GetTrack()->GetDynamicParticle()->GetPolarization() <<" ;\n " ;
|
|
if(aStep->GetPostStepPoint()->GetProcessDefinedStep()){
|
|
G4cout << "process : " << aStep->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName() <<"\n"
|
|
<<G4endl;
|
|
}
|
|
|
|
|
|
}
|
|
// }
|
|
}
|