XENONnT / fuse

Refactor XENONnT epix and WFsim code
BSD 3-Clause "New" or "Revised" License
6 stars 1 forks source link

Exciton Misplaced? #239

Closed FaroutYLq closed 1 month ago

FaroutYLq commented 3 months ago

What's the problem?

It seems to me that we might be doing something unphysical regarding treatment of excitons.

  1. Excitons are first computed here or there, together with photons and electrons before photon/electron collection efficiency applied. That means, the excitons we got here are also before photon/electron collection efficiency applied.
  2. The next and also last time excitons show up in fuse, is in S1 photon propagation. It was called in computing scintillation delay timing. Together with it in the same function nestpy_calc.GetPhotonTimes, the corresponding photons are already after collection efficiency.

Note that, we never implement a binomial process representing the loss of those excitons' photons in s1_photon_hits here. That being said, it is not surprising that you will have a lot of events with actual exciton fraction greater than 1, which is completely unphysical and already enters the wild extrapolation region in NEST. (Such breaking of physical sanity is even allowed in NEST!! See my PR for NEST here).

What it affects?

Of course, all the S1 pulse shape dependent studies are affected. However, as long as we have validated our simulation with data and got acceptable matches, we are not really wrong. There is no reason to interpret all the processes from NEST as physical, and it is OK we just use it as a grey box, with unphysical parameters fed into it.

How to do better?

A fix will be straightforward, just instead of only considering the total photons, apply two binomial processes and then find a way to return both photons and excitons:

n_recombined = n_photons - n_excitons
n_photon_hits_excitons = self.rng.binomial(n= n_excitons, p=ly)
n_photon_hits_recombined = self.rng.binomial(n= n_recombined, p=ly)
n_photon_hits = n_photon_hits_recombined +  n_photon_hits_excitons
dachengx commented 3 months ago

then find a way to return both photons and excitons

Do you mean propagate the n_photon_hits_excitons to S1PhotonPropagation plugin? In other words: change https://github.com/XENONnT/fuse/blob/75429f426376d93276de8c1039ab9b18acdc8c18/fuse/plugins/detector_physics/s1_photon_propagation.py#L381 to n_photon_hits_excitons?

FaroutYLq commented 3 months ago

then find a way to return both photons and excitons

Do you mean propagate the n_photon_hits_excitons to S1PhotonPropagation plugin? In other words: change

https://github.com/XENONnT/fuse/blob/75429f426376d93276de8c1039ab9b18acdc8c18/fuse/plugins/detector_physics/s1_photon_propagation.py#L381

to n_photon_hits_excitons?

@dachengx Yes exactly!

ramirezdiego commented 1 month ago

Resolved by #251.