Getting the energy deposited by a specific gamma photon in a multiparticle source?

Hello all,

I have a multivolume(different physical copy numbers) detector.
The primary generator generates Multiple gammas simultaneously at the start of the event.

Is it possible to get the total energy deposited by each gamma individually on an event-by-event basis? Currently, I am getting the total energies deposited in the various volumes.
However, it ‘includes’ energies from two different gammas in the same volume and ‘excludes’ energy by the same gamma into multiple volumes (for example by Compton scattering).

I think if there is some attribute that remains the same for all the particles produced(e-, scattered gamma) from a given particle(gamma), I can access that and sum the energies accordingly.

I am a beginner in the concepts of Track, Step Stack, etc, if someone can please offer any assistance/hints it will be really helpful.

Please let me know if I am not able to clear the issue or more info.

Thanking You

If you do run/beamOn 10, then you are going to have 10 events in 1 run.

So you can implement some counting, recording, etc in your endofeventaction function. For instance, you would record the energy deposited from one event (one gamma), output that information to a file, and at the end of your run you can have recorded the energy depositions for each event (each gamma)

HI @loydms,

Thanks for your reply.

I am currently dealing with a specific case where I have more than one gamma per event, so doing

/run/beamOn 10 

Will have 10 events in one run, and each event can have multiple gammas.

I think I have the event-by-event information on the total energy deposited. However, this energy deposited can be due to ‘pileup’ of more than one gammas generated in each event. I am wondering how to scoop out the individual energies deposited by different gammas from the total energy deposited.

Please let me know fo further clarifications.

Thank You!

If you want to know the energy of individual gammas, why not just run 1 gamma per event?

You could probably use particle id’s to record the total energy from each gamma in stepping action, but I would guess that it would be easier to split by events, unless there was some other reason why you are running multiple gammas per event.

Once again thanks fo your time and efforts.

The physical reason for having multiple gamma s in each event is to understand pileup effects in detectors for high multiplicity sources. However, I am unsure if this would be the same as generating multiple events with single gammas. Please correct me if I am wrong here.

By particle ID’s you mean parentID >? I tried getting the parent ID however I think at present it’s complex for me to analyze the output as there are multiple parentID and trackIDs for each event.

Maybe you can distinguish each one of your generated from their TrackID?

Each primary gamma, I presume. You can do it, but it requires an understanding of how Geant4 tracking works, along with a fair bit of code (more below). I would strongly recommend that you save an N-tuple with the individual track/step energy deposits, and defer the extra work to your analysis code.

If you want to know the total energy delivered by one of your primaries, then you need to keep track of not only that primary’s own energy deposits, but also all of the secondaries that it produced, and the energy deposits of those secondaries (and the granddaughter secondaries, and so on).

For each track, Geant4 carries the track ID itself, along with the ID of the parent that created it (the parent track won’t be available once it finishes, but the ID is guaranteed to be unique within an event). Primary tracks always have parent ID=0.

If there were no “granddaughters” (secondaries of secondaries), then you could use the parent ID exactly as this attribute, and collect the energy deposit from each primary, along with the deposits from the daughters of each primary. But that’s not guaranteed (look up “electromagnetic shower” or “hadronic shower” for extreme examples).

In your output N-tuple, store all the following:

  • detector information (copy number, name, however you make your sensitive detectors identifiable);
  • track ID and parent ID;
  • particle type (either the PDG code, or the name string);
  • energy deposited in the step
  • the track’s kinetic energy, position, etc., if you think they might be useful (for fiducial cuts, energy-loss calculations, whatever)

Once you’ve got an N-tuple, you can use ROOT, or Pandas, or whatever analysis tools you like to follow the chain of trackID-parentID through the steps and sum up the total energy deposit from each primary. Plus, you can adjust your analysis approach without having to rerun the simulation over and over.

1 Like

Hello @AlexAerPhys
Thanks a lot for your suggestion.
Yes, I am also trying to follow the TrackIDs and sum energies considering their particleIDs.
But I am struggling to find a suitable data analysis procedure to old all the track energies into separate primary particle energies.

Thanks a lot, @mkelsey, for your detailed explanation; it is very helpful for me.

Yes, correct. Sorry for my mistake.

Yes as you pointed out I think to sum the energies considering the direct and indirect parent-daughter associations is a bit tricky for me. I will try to follow some logic for that.

Here is a sample output for one event (number of primaries = 5)

ParentID: 1, TrackID: 7, t: 0.728001, Particle Name: e-, Energy: 0.0731583, Crystal: 4
ParentID: 1, TrackID: 7, t: 0.728148, Particle Name: e-, Energy: 0.156043, Crystal: 4
ParentID: 7, TrackID: 8, t: 0.728148, Particle Name: gamma, Energy: 9.298e-05, Crystal: 4
ParentID: 8, TrackID: 17, t: 0.72837, Particle Name: e-, Energy: 2.54e-05, Crystal: 4
ParentID: 8, TrackID: 16, t: 0.72837, Particle Name: e-, Energy: 9.26e-06, Crystal: 4
ParentID: 8, TrackID: 15, t: 0.72837, Particle Name: e-, Energy: 1.712e-05, Crystal: 4
ParentID: 8, TrackID: 14, t: 0.72837, Particle Name: e-, Energy: 8.63e-06, Crystal: 4
ParentID: 8, TrackID: 13, t: 0.72837, Particle Name: e-, Energy: 4.797e-05, Crystal: 4
ParentID: 8, TrackID: 12, t: 0.72837, Particle Name: e-, Energy: 4.86e-05, Crystal: 4
ParentID: 8, TrackID: 11, t: 0.72837, Particle Name: e-, Energy: 0.00097284, Crystal: 4
ParentID: 8, TrackID: 10, t: 0.72837, Particle Name: gamma, Energy: 4.864e-05, Crystal: 4
ParentID: 10, TrackID: 22, t: 0.728519, Particle Name: e-, Energy: 1.712e-05, Crystal: 4
ParentID: 10, TrackID: 21, t: 0.728519, Particle Name: e-, Energy: 1.712e-05, Crystal: 4
ParentID: 10, TrackID: 20, t: 0.728519, Particle Name: e-, Energy: 8.029e-05, Crystal: 4
ParentID: 10, TrackID: 19, t: 0.728519, Particle Name: e-, Energy: 1.608e-05, Crystal: 4
ParentID: 10, TrackID: 18, t: 0.728519, Particle Name: e-, Energy: 0.00966495, Crystal: 4
ParentID: 8, TrackID: 9, t: 0.72837, Particle Name: e-, Energy: 0.0172548, Crystal: 4
ParentID: 1, TrackID: 6, t: 0.64038, Particle Name: e-, Energy: 0.413284, Crystal: 4

I acknowledge that this question is not related to geant4 but mainly algorithmic. Do you have any suggestions/hints for doing this kind of summing? It will be very helpful.

At the moment, all I can think of is for each event

1.start with the smallest parentID, extract its event ID
2. Loop through the other parents ID and nested loop through the Track ID’s having same parent ID.

But I am clearly confused as how to control (exit) these statements for the loop.

Thank You

Also, I have a general query regarding the particle tracking, if you could offer some insight ?
The event output is not ordered in terms of parentID. Is it in some way that geant4 tracks? For example, the secondary data from a given parent (ParentID =1) are not tracked consecutively.

Naively, I expect it to track all the secondaries from a given parent first and then move to the next parent ID. However, looking at the output, it seems the tracking of ParentID 7 starts before the ending of ParentID 1.

Thank You

You should only need a single loop. Set up a couple of lookup tables:

  • std::map<G4int, G4double> (dictionary in Python), with the key being the trackID of your primary, and the value being the total energy deposited.
  • std::map<G4int, G4int> with the key being a trackID, and the value being the original primary parentID. (more below).

Loop over all the entries in your N-tuple for a single event (you’ll want to include EventID in your N-tuple; sorry I forgot to mention that!).

  • For each track, use the second map to see if you already know which primary it came from. If you do, great! If not, use that second map with this track’s parentID to see which primary the parent came from. Then add this track to that map.
  • Now for this same track, add its energy deposit to the total for the primary in the first map.
1 Like

That may depend on the process that created the secondary. One example I know of is that G4Optical processes (especially G4Cherenkov) pause the track, process the created G4OpticalPhoton secondaries first, then go back to the paused track to take another step. You may have other processes in your physics list that are configured that way.

Also, you can see that the track indexing is not sequential in the output. Geant4 puts each created track onto a LIFO stack, so the most recently created secondary is tracked first. If it creates more secondaries, then those will get tracked, before G4 eventually gets back to one of the other primaries.

1 Like

Thanks for that detailed explanation of analysis logic. I’ll try following it through.

Just another query is it safe to say that for every event, the parentID of the first track (produced anywhere, not limited to sensitive volume) corresponds to track id of the primary particle?

That should be true, but you shouldn’t need to rely on it for anything. The primary tracks will all have parentID == 0. That’s the canonical test.

Okay, thanks for your reply.

I think I am missing something here.
Hypothetically, suppose the data is ordered like this.

ParentID: 7, TrackID: 8, t: 0.728148, Particle Name: gamma, Energy: 9.298e-05, Crystal: 4
ParentID: 1, TrackID: 7, t: 0.728148, Particle Name: e-, Energy: 0.156043, Crystal: 4

Then (at least in my program), if it finds that both map 1 and map 2 are empty, it inserts the current trackID-parentID in map 2 and parentID-edep in map 1.

So, in this way, it will ‘mistake’ the current track’s parent as the original primary parent.
If it make things a bit clear, here the program

/****************** Looping across all track for each event ***************************************/
if (map1.count(tid) <= 0) {      //current track id is not in the track-parent map1
			if (map1.count(pid) <= 0){ //see if current track's parent is in map1
				map1[tid] = pid;        //if not add it as a key-value pair in map1
				map2[pid] = edep;       //enter(not sum) parentID-edep in map 2
			else {                      //its parent id is a track id in map1
				primary_pid = map1[pid];  //use map 1 to extrcat the primary parent id
				map1[tid] = primary_pid;    //add the track-primary parent ID in the map1
				if (map2.count(primary_pid) <= 0) //primary parent is in not in map2
						map2[primary_pid] = edep;       //enter it
						map2[primary_pid] += edep;     //primary parent id is in map2, add the track energy corresponding to that key
else {  //if the track id has primary parentID in the first map
			  primary_pid = map1[tid]; //extract primar parent id
			  if (map2.count(primary_pid) <= 0) //primary parent is in not in map2
				map2[primary_pid] = edep;       //enter it
	 			map2[primary_pid] += edep;     //primary parent id is in map2, add the track energy to that key

Can you please suggest if there is a more efficient way to deal with this as compared to the code above, basically for the very first track when both the emaps are empty.

Thanks a lot !

I’ll need to think about the problem and your code a bit. I see your point about the ordering; I think that can be solved in the parentage map, but I’ll have to think about how.

One important question – In your example output, none of the hits you list show, for example, a direct energy deposit from trackID 1, or from any other track with parentID=0. I assumed in my discussion that you would actually have hits in your SD from the primaries, so you’d be able to identify which ones they are.

If you are not getting hits from the primaries, that does make the problem much more complex (you have to have some way to know where the parentage chain terminates in a primary).

Thanks for your reply.

Sorry, I think I stated criteria that were too strict. However, it would be good if I could sum up the energy deposited in the SD by ‘unique’ parentIDs. I assumed all unique parentIDs correspond to a unique gamma primary parent, but now I can see this is false.

I will scan through the other output to see how often this complex case of not getting energy deposited by original primaries is happening.

Yes, I think that requires me to have a thorough understanding of the physics processes and probably filter out based on the first occurrence of a given process which leads to the termination of the parent chain into primary.

Thank You for your time and consideration

1 Like