u-eff-gee
commented
5 years ago Dear jvavrek,
Thank you for pointing out this issue and for the thorough documentation.
Indeed, the readout of the kinetic energy was in contradiction to what is written in the README. I applied your proposed fix to all detector classes. Now, the value of Ekin should make more sense.
So far, I did not notice this issue, because I usually do not read out Ekin when using an EnergyDepositionSD. This is because the EnergyDepositionSD stores information only at the level of G4Events (analysisManager->FillNtupleDColumn() is only called at the time of EnergyDepositionSD::EndOfEvent. Compare this to ParticleSD, where the ROOT file is filled at the level of G4Steps). Of course, several different particles may contribute to the total energy deposition, which is the quantity of interest for the EnergyDepositionSD. However, it only stores information about the first step (at least now it does) of the first particle in the event, since the actual filling of the ROOT file is done at the end of the event. I figured that the Ekin, position, momentum or particle information that an EnergyDepositionSD gives you would not be very useful anyway.
This comes from my background in Nuclear Resonance Fluorescence with single-crystal HPGes, where the only information accessible at the end is the total energy deposition by all particles in an event.
The pileup effect of the EnergyDepositionSD which you mention is caused by the same limitations of the sensitive detector. It can also be reproduced by shooting multiple particles in the same event. If the AngularCorrelationGenerator is used with the following macro commands
##################################
# Define the first step
##################################
# /angcorr/particle initiates the creation of a new cascade step.
# The macro will fail if this is not the first command in a step definition.
/angcorr/particle e+
/angcorr/energy 0.511 MeV
/angcorr/direction 1. 0. 0.
/angcorr/polarization 1. 0. 0.
##################################
# Define the second step
##################################
/angcorr/particle e-
/angcorr/energy 0.511 MeV
/angcorr/relativeangle 180. deg
/angcorr/polarization 1. 0. 0., an electron and a positron can be emitted in opposite directions inside a detector. They will be part of the same event. Similar to your pileup effect, the EnergyDepositionSD records
* Row * edep * ekin * particle * volume *
************************************************************
* 0 * 1.5329989 * 0.511 * 11 * 0 *, because the electron is emitted first. Note that it records an energy deposition which is even higher than the sum of the kinetical energies of the two emitted particles, because the positron was annihilated.
Actually, no sensitive detector exists at the moment that can store all the information that tracking/verbose gives you at the G4Step level. If you want to simulate gamma-ray tracking inside a detector in the GRETINA/AGATA style, then this would not be possible. In principle, ParticleSD can read out information about every hit point, but at the moment it does that only for the first step of a particle inside the sensitive detector. Particles are discriminated by their track ID and event ID. To have a 'true' tracking detector, one would need to remove the
if (trackID != getCurrentTrackID() || eventID != getCurrentEventID())condition in ParticleSD::ProcessHits.
jvavrek
Jack Silano (LLNL) and I have found that in the EnergyDepositionSD class, the kinetic energy appears to be tallied incorrectly. Based on the documentation, Ekin should be the kinetic energy of the incident particle, but based on experimentation, the kinetic energy is that after the primary's first interaction.
For example, firing 2500 keV photons at the ZeroDegree detector, a verbose output of EnergyDepositionSD::ProcessHits() gives the following:
Here we can see that the final value of Ekin = 1783.58 comes from the outgoing kinetic energy after the first Compton interaction at 2500 keV. Note also that this (quite often) causes Edep > Ekin.
A simple fix might be to change
hit->SetKineticEnergy(track->GetKineticEnergy());tohit->SetKineticEnergy(aStep->GetPreStepPoint()->GetKineticEnergy());This gives a better picture of what the true incident kinetic energy is, but doesn't perfectly solve the Edep > Ekin problem, as it doesn't account for pileup within the same event—multiple secondary photons depositing energy in the detector, e.g.:
Here two photons of energies 759.784 keV and 72.2999 keV both deposit their full energies in the detector, but only 759.784 keV is set as the value of Ekin.
Either way, I think the change recommended above is still a more intuitive measure of the 'Ekin' despite the pileup issue.