Scintillation By Particle Type bug
Kazu Terao found/diagnosed a bug in the optical scintillation light calculation in larsim. Some additional details can be found in MicroBooNE DocDB 8175, but I will summarize the main issues here:
---A simple bug is in ISCalculationSeparate.cxx (see https://cdcvs.fnal.gov/redmine/projects/larsim/repository/revisions/develop/entry/larsim/LArG4/ISCalculationSeparate.cxx#L142). On line 148, a pointer to a material property vector (Scint_Yield_Vector ) is set to null. An if/else collection follows, setting a scintillation yield based on particle type (where per-particle scintillation yields are defined in the material properties table, which itself is defined from LArG4 parameters [ or LArProperties ... I forget offhand ] ). On line 191, the Scint_Yield_Vector is checked to see if it is still null, and if so defaults to using the electron scintillation yield. Scint_Yield_Vector is not modified anywhere in the if statement, and so remains null, thus always defaulting to the electron scintillation yield.
---Beyond this, there is a deeper problem.
First, casting of float to integer in the calculation of the number of scintillation photons (see https://cdcvs.fnal.gov/redmine/projects/larsim/repository/revisions/develop/entry/larsim/LArG4/ISCalculationSeparate.cxx#L197 and following lines). If the number photons is not large, then this will lead to a systematic underestimation of the amount of scintillation light. Because MicroBooNE, and perhaps others, usually include a quantum efficiency and potentially other effects as part of their scintillation yield, this can lead to small numbers of photons per voxel. Note that this casting to int is required by the ISCalculation base classes.
Later usage does not support the choice of int in ISCalculation for the number of scintillation photons. See https://cdcvs.fnal.gov/redmine/projects/larsim/repository/revisions/develop/entry/larsim/LArG4/OpFastScintillation.cxx#L320, which immediately casts the output back to double, and is then used as the mean of a poisson random number generator for determing the number of photons to use at the optical detectors.
I believe then that two things need to be done:
(1) Remove the practically-but-not-literally unused pointer in the ISCalculationSeparate determination of yield by particle type.
(2) Rethink the ISCalculation interface. I think this thinking is not complicated: it seems clear that the expectation of later code and the expectation of a user generally would be to return the expected average number of scintillation photons / ionization electrons given an amount of energy deposited. There is no reason this could not be non-integer, and so the return type should be float or double (float is likely fine). The name of the variable/accessor could be changed to make this more explicit too, if desired (e.g. MeanNumberScintillationPhotons() or ExpectedNumberScintillationPhotons() ).
Obviously this will require coordination with the experiments to make sure the above would be ok, but as it exists it's a bug for everyone, and it would be important to see fixed ASAP. (I think the fixes are simple.)
I'm adding Kazu to the watch list so he can provide any additional insight or details I've forgotten.
#2 Updated by Gianluca Petrillo over 2 years ago
I think the changes are fine.
I have pushed branch
larsim updated to
Since Kazuhiro has also converted the number of ionisation electrons to floating point, I have updated the downstream code accordingly.
Please note that NEST still uses integers in its internal computation.