the methods i write may print out wrong results, but i can’t find out where the wrong happened.i have point the fdetConstruction to the shape2 physic volume which i want.After running one million electrons whose energys are both 100MeV hit target Wolfram,the Histogram reveals a strange image ,especially in the kinetic energy of gamma.almost all gamma i collect around 0.511MeV,then the energy spectrum is approximately exponential fall.if this code below doesn’t have problems, may some questions about physics program . i just import QGSP_BERT at main function.
#include “S2SteppingAction.hh”
#include “S2Analysis.hh”
#include “S2EventAction.hh”
#include “S2DetectorConstruction.hh”
#include “G4ParticleTypes.hh”
#include “G4Step.hh”
#include “G4RunManager.hh”
#include
#include “G4SystemOfUnits.hh”
//////////////
S2SteppingAction::S2SteppingAction(
const S2DetectorConstruction* detectorConstruction,
S2EventAction* )
: G4UserSteppingAction(),
fDetConstruction(detectorConstruction)
{}
//////////////
S2SteppingAction::~S2SteppingAction()
{}
//////////////////////////////////////
void S2SteppingAction::UserSteppingAction( const G4Step* step)
{
auto analysisManager = G4AnalysisManager::Instance();
const G4StepPoint* preStepPoint = step->GetPreStepPoint();
auto volume= preStepPoint->GetTouchableHandle()->GetVolume();
const G4Track* track=step->GetTrack();
if (volume == fDetConstruction->GetPV()&&(preStepPoint->GetStepStatus() == fGeomBoundary)){
const G4ParticleDefinition* particle = track->GetPaticleDefinition();
const G4double Kenergy = preStepPoint->GetKineticEnergy();
if (particle==G4Electron::Electron()) {
analysisManager->FillH1(1,Kenergy);}
else if (particle==G4Positron::Positron()) {
analysisManager->FillH1(2,Kenergy);}
else if (particle==G4Gamma::Gamma()) {
analysisManager->FillH1(0,Kenergy);}
}
}
Your check on the preStep status is correct, but for the volume check, you want to see the new postStep volume being entered. What you’re doing is looking for tracks which have just left (preStep) your volume of interest.
Just change the auto volume assignment to use the post-step point.
if i remove the target, the result are all right.However,if i use gamma and target W,the histogram is below
Good morning(maybe not morning to you)!i am not sure it is suitable to call you mike.if not ,i will apologize .but first, thank you very much!
i want to know whether the way i use 'if(particle ==)'to distinguish particles is feasible or not?
However,i have chosen vacuum to be the materials of shape2 and put the shape2 next the target without interval.
Ps:when i am writing repayment above, i have changed the auto volume by’const G4StepPoint* postStepPoint = step->GetPostStepPoint();and auto volume = postStepPoint->GetTouchableHandle()->GetVolume(); ',but the result is still have the problems i said before .
i want to know whether the way i use 'if(particle ==)'to distinguish particles is feasible or not?
Yes. In fact, what you do there is the preferred way to match particles: each one is a global singleton, and comparing name strings is more CPU intensive than a simple pointer comparison.
I notice that you’re collecting the particle energy at the preStep point, which means before interacting. Is that correct? Or do you want to see those particles after they’ve interacted?
Since you’re shooting primary electrons, are you wanting to find out the energy of generated gammas are? In that case, you may want a more complex stepping action, where you look for a primary electron which interacted, then look at the secondaries in the postStep of that electron, and where the secondary is a gamma, you collect its initial energy. You can look at the G4StepPoint.hh to understand what to call to get that information.
Good afternoon,Mike! maybe i have found one of the questions .In steppingAction,i export the particle’s position and the name of volume when it meets the condition ‘if (volume == fDetConstruction->GetPV()&&(preStepPoint->GetStepStatus() == fGeomBoundary))’. i set the target whose thickness is 7mm at (0,0,0) and the volume shape2 whose thickness is 16cm at (-17.5 cm,0,0). the function GetPV lists below
inline const G4VPhysicalVolume S2DetectorConstruction::GetPV() const
{
return fShape2;
}
however, the output reveal :'gammaof7.8785 -0.35atShape2. ’
Shape2 is the name of physic volume shape2.if i change the thickness of target(not the shape2),the ‘0.35’ will become the half thickness of target.
if you need other section of my code ,please tell me. i will paste immediately .
###PS:
const G4StepPoint* preStepPoint = step->GetPreStepPoint();
const G4StepPoint* postStepPoint = step->GetPostStepPoint();
auto volume = postStepPoint->GetTouchableHandle()->GetVolume();
auto volume1= preStepPoint->GetTouchableHandle()->GetVolume();
const G4Track* track=step->GetTrack();
if (volume == fDetConstruction->GetPV()&&(preStepPoint->GetStepStatus() == fGeomBoundary))
{
const G4String& name = volume->GetName();
const G4String& name1 = volume1->GetName();
const G4ThreeVector& preStepPosition = preStepPoint->GetPosition();
const G4ThreeVector& postStepPosition = postStepPoint->GetPosition();
G4double Px1= preStepPosition.x();
G4double Px= postStepPosition.x();
const G4ParticleDefinition* particle = track->GetParticleDefinition();
G4String pname = particle->GetParticleName();
const G4double Kenergy = postStepPoint->GetKineticEnergy();
const G4double Eenergy = preStepPoint->GetKineticEnergy();
std::cout<<"=========>>>>>"<<pname<<“of”<<Kenergy/MeV<<" “<<Px/cm<<“at”<<name<<”==>"<<Eenergy/MeV<<Px1/cm<<“at”<<name1<<std::endl;}
the output list below:
=========>>>>>gammaof0.510999 -9.5atShape2==>0.510999-0.35atWorld
=========>>>>>gammaof10.4968 -9.5atShape2==>10.4968-0.35atWorld
=========>>>>>gammaof0.386142 -9.5atShape2==>0.386142-0.35atWorld
=========>>>>>gammaof8.00037 -9.5atShape2==>8.00037-0.35atWorld
=========>>>>>gammaof0.554786 -9.5atShape2==>0.554786-0.35atWorld
=========>>>>>gammaof5.264 -9.5atShape2==>5.264-0.35atWorld
=========>>>>>gammaof1.89981 -9.5atShape2==>1.89981-0.35atWorld
=========>>>>>gammaof3.12891 -9.5atShape2==>3.12891-0.35atWorld
=========>>>>>gammaof0.252364 -9.5atShape2==>0.252364-0.35atWorld
=========>>>>>e+of15.9042 -9.5atShape2==>15.9042-0.35atWorld
=========>>>>>gammaof4.96709 -9.5atShape2==>4.96709-0.35atWorld
=========>>>>>gammaof0.463188 -9.5atShape2==>0.463188-0.35atWorld
=========>>>>>gammaof0.50737 -9.5atShape2==>0.50737-0.35atWorld
=========>>>>>e-of25.0006 -9.5atShape2==>25.0008-0.35atWorld
=========>>>>>gammaof2.47443 -9.5atShape2==>2.47443-0.35atWorld
=========>>>>>gammaof7.8785 -9.5atShape2==>7.8785-0.35atWorld
=========>>>>>gammaof1.75571 -9.5atShape2==>1.75571-0.35atWorld
=========>>>>>gammaof1.96374 -9.5atShape2==>1.96374-0.35atWorld
=========>>>>>gammaof0.965453 -9.5atShape2==>0.965453-0.35atWorld
Yes, you’re getting exactly what you have asked for. You are printing out the preStep position. You have asked for the preStep position to be the volume boundary. Therefore you are printing out the volume boundary position. You have placed your target at (0,0,0). Therefore the X coordinate of the volume boundary is at half the dimension along the X axis.
Morning Mike!,if i set maximum step length <1mm,and then the preStepPiont will stay at Shape2 other than World in base of former code?
how can i set the boundary of shape2 other than target as the Gemboundary？
I’m afraid I don’t really understand your questions. The boundary of any shape is set by that shape itself: The boundary of a cube is the surface of that cube. The boundary of a sphere is the surface of that sphere.
sorry，I didn’t express my idea appropriately. the key questions I want to ask the particle attributes saved in track is equal to the attributes saved in poststeppoint？
Almost. The track will get the values from the postStep point after the step is completed. In your stepping action, it’s safer to get the quantities directly from the step. If you decide to write a TrackingAction, then in your PostUserTrackingAction the postStep values will have already been copied.