Plot primary energy spectrum

Hi all,

Fairly new to Geant4.

I’m generating primaries using GeneralParticleSource from a macro as follows:

/gps/ene/type User
/gps/hist/type energy
/gps/hist/point 1.0 0.
/gps/hist/point 2.0 5.
/gps/hist/point 7.0 1.
/gps/hist/point 10. 1.

I’d like to be able to see a plot of the primary’ energy spectrum.

Is there a simple way of plotting the energy spectrum of the primaries?

Thanks!

In your UserEventAction, you can get the information about the generated primary in each event via G4Event::GetPrimaryVertex(G4int i) (use G4Event::GetNumberOfPrimaryVertex() to get the count and loop). From each vertex, you can get the individual primary particles from G4PrimaryVertex::GetPrmary(G4int i) and G4PrimaryVertex::GetNumberOfParticle().

If you used GPS to throw a single particle per event, you could skip the looping and just get index 0 in each case, but then you’d have to modify your code if you decided later to throw multiple primaries.

You could also make use of your SteppingAction or TrackingAction to look for tracks with ParentID==0, and get the initial energy that way. The method above has the advantage that it works any time during your event; tracks are deleted as soon as they finish tracking, so you can’t get back to them later.

1 Like

Hi all,
Thank you @DMH and @mkelsey for this topic.
I am a new user of Geant4 and I am working on the part of initializing the primary particles.
Your topic gives very good hints.
But I am having a hard time understanding how I can send just a single particle (single vertex) per event, as you (@mkelsey) mentioned.
I am following the instructions in https://indico.cern.ch/event/776050/contributions/3240644/attachments/1788889/2913530/PrimaryParticle.pdf slide nr 5.
And I am getting a different nr of primary particles per event.
Could you please point me to an example or documentation section that could help me with this?
Thank you very much,
Diana

Hi @LDAmorim

As you provide the slides, it looks you are using ParticleGun.
firstly, make sure if number of particle is one.

=>G4ParticleGun(n_particle);<=

where ‘n_particle’ should be “one”. It means “n_particle=1”.
secondly, I suggest you can use “/tracking/verbose 1” to run one event and debug your behavior of primary particle.

Hope these will helpful for you.

Ye

Hi,
Thank you for your reply and sorry for getting back to you so late. I am using n_particle=1.
I tested with verbose and I even tried to add fparticleGun-> SetNumberOfParticles(1);.
But I still get a nr of primary particles equal to 2.5 x the nr of events runManager->BeamOn(events).
Do you have any other idea on how to fix this?
Thank you for your help,
Diana

Hi again,
Thank you so much for all your suggestions. I discovered the error was coming from how I was storing the track information in each step in my code (without differentiating between event and step numbers). So thank you for helping me make it work.
Cheers,
Diana