Photon statistics in GLG4sim

From: Stan Seibert (volsung@physics.utexas.edu)
Date: Thu Dec 01 2005 - 14:21:08 CST


Hi all,

Just to re-explain what I was saying during the meeting (since it was
pretty unclear) about how scintillation photons are handled in GLG4sim:

Current Method ("track weighting" method)
----------------------------------------------------------
There is a parameter called meanPhotonsPerSecondary. This number
defaults to 3.0, and (if I understand the code correctly, is used
this way:

The scintillation process calculates how much energy the charged
particle deposits in the scintillator during a step. It then
computes the mean number of photons expected for this energy
deposit. To decide how many optical photon tracks to give to GEANT4
to propagate, it picks a random number from a distribution whose mean
is meanNumPhotons / meanPhotonsPerSecondary. So by default, the mean
number of actual, propagating photons is one-third the true number.

This is compensated for by attaching to each photon track a weight
(3.0 in our case) which it carries along. I'm not entirely clear how
GEANT4 handles weight tracks in other processes, but at the other
end, when the photon hits the photocathode of a PMT, two checks are done

First, it calculates the mean number of expected photoelectrons,
which is the track weight (in a sense, the multiplicity of the photon
bundle) times the collection efficiency (absorption prob * quantum
efficiency). The actual number of produced photoelectrons is this
mean plus a random number X (range [0,1)), rounded to the nearest
integer. So, the higher the mean_pe, the more likely mean_pe+X is
greater than 1, and at least one photoelectron is detected. (And for
very high weights, multiple photoelectrons could be produced.)

Second, it decides if the track is absorbed. Using the same random
number as before, it is checked if (1-X) is less than the absorption
probability. If so, the track is absorbed and discarded. So, if the
track produced photoelectrons, it is probably going to be absorbed,
but not always, since there is some wiggle room, especially if
mean_pe is very large. Glen has copious comments in the code
claiming that he gets the statistics correct for weight = 1, and
mostly correct for weight > 1, but I still don't quite get it.

As a side note, I have confirmed this code produces slightly
different results (nhit shifts of several percent) when
meanPhotonsPerSecondary is 1 (no weighting) and when it is 3 (default).

Proposed Method ("efficiency scaling" method)
----------------------------------------------------------------
The method I have been trying to implement (and not having much
success yet) is a little easier to understand. The basic idea is to
try to throw photons that you "know" you are going to detect. So,
you select a "photon factor", much like the meanTracksPerSecondary
before.

When you generate photons (like in the scintillator), you scale the
number down by the photon_factor, which results in less tracks to
propagate, just as before. Each track however carries a weight of 1
still. (i.e. represents just a single photon)

At the photocathode, you scale up the detection probability by the
photon factor. Thus, the photon factor cannot be greater than 1/
peak_collection_efficiency, or you will find yourself needing to
detect photons with a probability greater than 1 to compensate for
the reduced production of photon tracks.

This method is not hard, but it requires the photon factor to be
applied everywhere optical photons are being produced. I think I am
currently missing some of these points. However, it should give you
identical results with photon_factor = 1 vs. photon_factor = 3, once
I can implement it right.

So that's the state of things. Clearly, we will leave Glenn's method
in place for the Version 0.1 release. The question I have is whether
we should set the default meanPhotonsPerSecondary=1, effectively
disabling the statistical fanciness until we can understand it or
replace it. The drawback is that our code will now run nearly 3
times slower.

---
Stan Seibert


This archive was generated by hypermail 2.1.6 : Fri Dec 02 2005 - 00:01:02 CST