Bug #24433

Implicit casting double to integer on larsim/LegacyLArG4/OpFastScintillation

Added by Iker de Icaza Astiz about 2 months ago. Updated 1 day ago.

Under Discussion
Target version:
Start date:
Due date:
% Done:


Estimated time:
Occurs In:


Another issue I came across whilst refactoring OpFastScintillation and that is present on larsim/PhotonPropagation/

The functions VUVHits() and VISHits() return the number of photons a given Optical Detector is expected to detect from a Scintillation Point inside the detector. They compute the apparent solid angle of the optical detectors from the perspective of the Scintillation Point and do some other corrections. Their input, amongst others, is the number of photons created in the scintillation point. This variable is treated as a double everywhere until is passed to the aforementioned functions. The relevant lines are:

499: double MeanNumberOfPhotons = larg4::IonizationAndScintillation::Instance()->NumberScintillationPhotons();
582: double Num = 0;

      if (scnt ==1){
         if ( ExcitationRatio == 1.0 ) {
            Num = std::min(YieldRatio,1.0)*MeanNumberOfPhotons;
          else {
            Num = std::min(ExcitationRatio,1.0)*MeanNumberOfPhotons;
        Num = MeanNumberOfPhotons - Num;  
1530:      DetThis = VUVHits(Num, ScintPoint,
1531:              ,;

1592:      ReflDetThis = VISHits(Num, ScintPoint,
1593:                   ,,
1594:                             cathode_hits_rec, hotspot);

The current declarations of VUVHits() and VISHits() truncate Num to an integer when they get called. They also round their calculation inside to return an integer.

I understand that physically only an integer number of photons are produced and detected. It is not clear to me wether is meaningful in these simulations to truncate them to integers. I suppose there's arbitrary values being used by g4, that call for MeanNumberOfPhotons and Num to be doubles.

My expectation is that if the desired behaviour of Num was to be a double everywhere until it was used for VUVHits and VISHits, then it would have been explicitly truncated using floor(). Hence, I believe using them as integers when they are originally double was an omission.

cosmics_g4_icarus_volCryostat-20200504-00-G4profiling.png (62.1 KB) cosmics_g4_icarus_volCryostat-20200504-00-G4profiling.png Screenshot of ICARUS cosmic ray LegacyLArG4 simulation job (5 events) Gianluca Petrillo, 05/26/2020 06:44 PM


#1 Updated by Iker de Icaza Astiz about 2 months ago

On my fork I have tested the outcome of changing the declarations of VUVHits() and VISHits() to use Nphotons_created (ie Num) as doubles.

The commit:

In the context of SBND, this yields about ~5% more VUV and VIS photons for light generated across the whole TPC. Clearly this affects all experiments using this code, I suspect that in similar fashion.

Intuitively I would assume a 5% discrepancy would have been observed when this method was validated. If that's the case, I might be opening a completely bogus concern or there might be some other issue in action that's balancing this discrepancy.

Again, I felt this called for larger input.

#2 Updated by Kyle Knoepfel about 2 months ago

  • Status changed from New to Under Discussion

This will require a discussion among LArSoft stakeholders. Thank you for bringing it up. Wei Mu, do you have any comments to make regarding Iker's findings?

#3 Updated by Lynn Garren about 1 month ago

From Wei Mu:

I agree with Iker that it will be better to truncate Num explicitly. In another photon simulation module (, I use round() to truncate the variables (num of photons and num of ions). Since the number of photons is eventually integer in OpDetBacktrackerRecord and SimPhotonsLite, truncating it to integer makes more sense.

The reason why MeanNumberOfPhotons and Num are doubles, I think, is they are calculated from variables that are doubles. It might be an omission using num of photon as integers in current code since there are also other functions which arbitrarily use double to read the num of photons (integer), such as AddScintillationPhotons().

#4 Updated by Alexander Himmel about 1 month ago

Even if what's in the module now faithfully reproduces what was in the legacy implementation (that's what would have been seen in the validation), we still probably want to take a step back and think intentionally about where the move from double -> integer happens. I actually think all of the simulation at this stage should be in floating point numbers, and we should only be truncating to integers once we reach the digitizer, apply quantum efficiency, and need to start making a waveform out of discrete photons.

#5 Updated by Gianluca Petrillo about 1 month ago

In LegacyLArG4 quantum efficiency can be applied there already (there are OnePhoton enough without saving the 80-90% of them that won't make it).

Back to Iker's report...
The rounding close to Num in the "old" code is after an extraction from Poisson statistics, which returns an integral value but represented in floating point because of the random library output; although redundant, it (expensively) clarifies an intention.
The new code in VUVHits() applies the Poisson statistics later, so I think that Num number should not be converted into integral.
In VISHits it's hard to guess the intention of the author, since the parameter of Num is unused.

But either roads, Poisson statistics is used and an integral number of photons is obtained. If this is not desired, a continuous distribution must be used instead. Ah, did I mention that that Poisson extraction is already the single most CPU-demanding operation of ICARUS LArG4 job? (see attachment)

#6 Updated by Iker de Icaza Astiz about 1 month ago

Hi Gianluca,

It's actually my fault that VISHits no longer uses Num, you see I realised that I could take it out of the function because the portion that used it was repeatedly running the same operation over and over on a loop. Missed to clean that up. As a comment: it was because of this move that I realised the impact of the integer truncation.

However if you see the commit I linked above there I'm removing the cast to int:

      (solid_angle_cathode / (4.*CLHEP::pi)) * int(Num);
      (solid_angle_cathode / (4.*CLHEP::pi)) * Num;

I put this cast on the refactored code to be able to recreate the behaviour.

I'm rather surprised about how slow the G4Poisson() is, there are a few gRandom->Poisson(), which are getting called as often and don't show up on your profiler.

#7 Updated by Iker de Icaza Astiz 1 day ago

I have created PR to address this issue.

I left the return type of VUVHits() and VISHits() as integers, this is because the implementations that handle the timing rely on integers. Hence it would be a change with many more implications and much harder to do.

Also available in: Atom PDF