Multiple detections of the same particle in sensitive detector

Hello dear experts!
I’m trying to simulate multiple plastic scintillators and i wanna record the position of impact when a particle(not a photon) hits each scintillator. I have the following code:

G4bool MyBarDetector::ProcessHits(G4Step *aStep, G4TouchableHistory *ROhist){
    
    
    G4AnalysisManager *man = G4AnalysisManager::Instance();

    G4VPhysicalVolume* physicalVolume = aStep->GetPreStepPoint()->GetTouchableHandle()->GetVolume();
    G4int copyNumber = physicalVolume->GetCopyNo();
    G4Track *track=aStep->GetTrack();
    //track->SetTrackStatus(fStopAndKill);
    G4StepPoint *preStepPoint=aStep->GetPreStepPoint();
    G4StepPoint *postStepPoint=aStep->GetPostStepPoint();
    const G4ParticleDefinition* particleDefinition = track->GetDefinition();
    G4String particleName = particleDefinition->GetParticleName();
   // G4cout<<"particula detectada"<<G4endl;
    G4ThreeVector pos = preStepPoint->GetPosition();
    if(particleName != "opticalphoton"){
        G4cout<<"Particula: "<<particleName<<" coordenadas: "<<pos<<G4endl;
    }

    return true;

}

However, it seems that geant4 records the position of the same particle when it passes through my detector since i only shooted one muon and my terminal prints the following for one scintillator:
“G4WT0 > Particula: mu- coordenadas: (-0.949802,1.09331,800)
G4WT0 > Particula: mu- coordenadas: (-0.952262,1.09892,806.405)
G4WT0 > Particula: mu- coordenadas: (-0.950821,1.10216,812.81)”

Can anyone help me, please ? I need to record the position of impact, not several interactions

Hello
You have to kill the track before checking the position of impact.
Use track->SetTrackStatus(fStopAndKill);
before getting a pre-step point. This should work in my opinion.

@deepachoudhary That’s close, but not quite right. Your suggestion would end up with the OP stopping the muons before they completed their journey. What they’re asking is how to identify and record only the final position (if any).

Are you certain that the muons are actually going to stop in your scintillator assembly? Muons are generally minimum-ionizing, and will travel through a large thickness of light material without stopping.

If you’re sure you have stopping muons, then you can test (not set) the track status in your SD:

if (track->GetTrackStatus() == fStopAndKill) {
    G4cout << "Hooray, the track stopped here! " << track->GetPosition() << G4endl;
}
1 Like

Sorry, i did not explain myself well enough. Indeed, the muons will not stop in my scintillator assembly, they will pass through my detectors and continue until they exit my world. What i want to do is record the position where they get in contact the very first time with each scintillator without killing the track. I also tried your suggestion but i got the following output:

G4WT0 > Hooray, the track stopped here! (-0.330647,1.03567,923.064)
G4WT0 > Particula: e- coordenadas: (-0.330647,1.03567,923.064)

It did not record primary particles(muons).

I tried saving the trackID of each particle and the copyNumber of the detector that it hits my, in a map of the form

std::map<G4int , std::set<G4int>> registeredParticlesMap; //The first argument is the copyNumber

and then, in each hit, look for the trackID in my map and only record when the trackID is not found, which seemed to work just fine. However, right now i am testing my code and my main goal is to shoot a beam of particles and record secondary particles that hit my detectors, so this would be highly inefficient. Is there any other way to do it? I don’t know if Geant4 has predetermined way of doing this since i am relatively new at using the program

So you want to record just the first hit (not the total energy deposit) in each scintillator layer? Your use of a std::map is quite clever, but as you suspect, G4 does have a way of probably doing what you want. G4Step::IsFirstStepInVolume() returns true when a particle enters the volume. Assuming that every step deposits “some” energy (a valid assumption for muons, which are charged), then that will give you what you want.

Unfortunately, that didn’t work either. I got the following output:

G4WT0 > Particula: mu- coordenadas: (-0.949802,1.09331,800) CopyNo: 2
G4WT0 > Particula: mu- coordenadas: (-0.952262,1.09892,806.405) CopyNo: 2
G4WT0 > Particula: mu- coordenadas: (-0.950821,1.10216,812.81) CopyNo: 2

Let me show you the piece of code i modified:

G4bool MyBarDetector::ProcessHits(G4Step *aStep, G4TouchableHistory *ROhist){
    
    

    G4Track *track=aStep->GetTrack();
    const G4ParticleDefinition* particleDefinition = track->GetDefinition();
    G4String particleName = particleDefinition->GetParticleName();
    if(particleName == "opticalphoton"){
        return true;
    }
    
    
    if (aStep->IsFirstStepInVolume()) {
        G4AnalysisManager *man = G4AnalysisManager::Instance();

        G4VPhysicalVolume* physicalVolume = aStep->GetPreStepPoint()->GetTouchableHandle()->GetVolume();
        G4int copyNumber = physicalVolume->GetCopyNo();
        //track->SetTrackStatus(fStopAndKill);
        G4StepPoint *preStepPoint=aStep->GetPreStepPoint();
        G4StepPoint *postStepPoint=aStep->GetPostStepPoint();
        
    // G4cout<<"particula detectada"<<G4endl;
        G4ThreeVector pos = preStepPoint->GetPosition();
        
        G4cout<<"Particula: "<<particleName<<" coordenadas: "<<pos<<" CopyNo: "<<copyNumber<<G4endl;
    }
    
    

    return true;

}


It keeps showing me multiple interactions in the same volume and i don’t know why. Do you have any other suggestion, please?

The code looks right to me. Are you sure those print statements are all for the same event and the same track? Try adding a couple of additional outputs to your diagnostic.

  • Print the track number, aStep->GetTrack()->GetTrackID().
  • Print the event number, G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID().

The relevant intermediate pointers are guaranteed to exist within your SD’s ProcessHits() function (so you don’t need to do null testing). You’ll need to add header files to your code and recompile,

#include "G4EventManager.hh"
#include "G4Event.hh"
#include "G4Track.hh"

Here goes the output after i implemented what you suggested. This time i am also including the information of the second scintillator.

G4WT0 > Particula: mu- coordenadas: (-0.949802,1.09331,800) CopyNo: 2 TrackID: 1 EventID: 0
G4WT0 > Particula: mu- coordenadas: (-0.952262,1.09892,806.405) CopyNo: 2 TrackID: 1 EventID: 0
G4WT0 > Particula: mu- coordenadas: (-0.950821,1.10216,812.81) CopyNo: 2 TrackID: 1 EventID: 0
G4WT0 > Particula: mu- coordenadas: (-1.13761,1.00753,910) CopyNo: 6 TrackID: 1 EventID: 0
G4WT0 > Particula: mu- coordenadas: (-1.14058,0.997202,916.405) CopyNo: 6 TrackID: 1 EventID: 0
G4WT0 > Particula: e- coordenadas: (-1.13519,0.977785,922.442) CopyNo: 6 TrackID: 38500 EventID: 0
G4WT0 > Particula: e- coordenadas: (-0.623116,1.03637,922.78) CopyNo: 6 TrackID: 38500 EventID: 0
G4WT0 > Particula: e- coordenadas: (-0.489711,0.982879,923.036) CopyNo: 6 TrackID: 38500 EventID: 0
G4WT0 > Particula: e- coordenadas: (-0.480754,0.985098,923.038) CopyNo: 6 TrackID: 38500 EventID: 0
G4WT0 > Particula: e- coordenadas: (-0.474705,0.986368,923.039) CopyNo: 6 TrackID: 38500 EventID: 0
G4WT0 > Particula: e- coordenadas: (-0.330647,1.03567,923.064) CopyNo: 6 TrackID: 38500 EventID: 0
G4WT0 > Particula: mu- coordenadas: (-1.13519,0.977785,922.442) CopyNo: 6 TrackID: 1 EventID: 0
G4WT0 > Particula: mu- coordenadas: (-1.13235,0.953112,928.847) CopyNo: 6 TrackID: 1 EventID: 0

The first three lines corresponds to what i’ve been showing you and, indeed, those print statements do correspond to the same event and track.
Is there any chance this is a bug in Geant4?

That’s quite peculiar, and not what I was expecting. You can see all the optical photons clearly – the electron’s track ID was assigned after all of them were generated and tracked.

I wouldn’t expect so, but “anything is possible.” You can test for that relatively easily, by manually doing what that function should be doing already.

The first “full step” in a new volume will have the preStepPoint at the boundary of the volume, and the postStepPoint either within the volume, or at the outgoing boundary. If you modify your code to move the two G4StepPoint accesses outside the if block, then you can write:

  G4StepPoint *preStepPoint=aStep->GetPreStepPoint();
  G4StepPoint *postStepPoint=aStep->GetPostStepPoint();

  if (preStepPoint->GetStepStatus() == fGeomBoundary) {
    // Do your stuff here, including printing things out
  }

… I do notice one thing. Why are you using preStepPoint rather than postStepPoint for your position?

I implemented your suggestion, however, this time my console printed nothing. Then, I’m going to try and test if this is a bug unless you have another suggestion, i ran out of ideas.

Regarding your question, i really never thought about it, for my purposes i assumed it really did not matter because as I had understood, one gives you the position at the begging of the step and the other at the end and the length of the step is short. Should i use postStepPoint instead?

If the length of the the step is short enough (and that is usually the case with charged tracks), then you’re right that it doesn’t really matter in a practical sense. When you have dE/dx, or other “along step” energy, that energy loss is assigned to the interval between the pre-step and post-step points, and secondaries are emitted from the post-step point. So I tend to use the post-step position in my work, because “that’s where stuff comes from.”