Merge pull request #36 from DUNE/35-modifyfix-scintillation-time-func…

…tion

35 modifyfix scintillation time function

include/utils.h

@@ -165,6 +165,8 @@ namespace scintillation
 const double tscin = 3.08;
 const double tscex = 0.588;
 const double vlfb = 5.85;
+const double rise_time = 0.7; // ns
+const double decay_time = 3.0; // ns
 }  // namespace scintillation
 
 namespace photo_sensor

src/SANDDigitization.cpp

@@ -13,9 +13,26 @@ TRandom3 rand(0);
 
 namespace ecal
 {
-// simulate pe arrival time to pmt
-double photo_electron_time_to_pmt_arrival_time(double t0, double d)
+
+// simulate time of emission of the scintillation photons
+double photons_emission_time()
+{
+  // https://github.com/Geant4/geant4/blob/e58e650b32b961c8093f3dd6a2c3bc917b2552be/source/processes/electromagnetic/xrays/src/G4Scintillation.cc#L638
+  double t;
+
+  do
+  {
+    // The exponential distribution as an envelope function: very efficient
+    t  = -1.0 * sand_reco::ecal::scintillation::decay_time * TMath::Log(1.0 - rand.Uniform() );
+  }
+  while (rand.Uniform()  > (1.0 - TMath::Exp(-t/sand_reco::ecal::scintillation::rise_time)));
+
+  return t;
+}
+
+double photons_emission_time_as_kloe()
 {
+
   /*
      - For each photoelectron: Time for TDC simulation obtained from
 
@@ -34,9 +51,34 @@ C                      + 1ns  uncertainty
       TSCEX  0.588
   */
 
-  double tdec = sand_reco::ecal::scintillation::tscin *
+  return sand_reco::ecal::scintillation::tscin *
                 TMath::Power(1. / rand.Uniform() - 1.,
                              sand_reco::ecal::scintillation::tscex);
+}
+
+
+// simulate pe arrival time to pmt
+double photo_electron_time_to_pmt_arrival_time(double t0, double d)
+{
+  /*
+     - For each photoelectron: Time for TDC simulation obtained from
+
+C  PHOTOELECTRON TIME :  Particle TIME in the cell
+C                      + SCINTILLATION DECAY TIME +
+C                      + signal propagation to the cell
+C                      + 1ns  uncertainty
+
+               TPHE = Part_time+TSDEC+DPM1*VLFB+Gauss(1ns)
+
+      VLFB = 5.85 ns/m
+!!!! Input-TDC Scintillation time -
+               TSDEC = TSCIN*(1./RNDMPH(1)-1)**TSCEX  (ns)
+
+      TSCIN  3.08  ns
+      TSCEX  0.588
+  */
+
+  double tdec = photons_emission_time();
 
   double time = t0 + tdec +
                 sand_reco::ecal::scintillation::vlfb * d * conversion::mm_to_m +
@@ -149,4 +191,4 @@ void eval_adc_and_tdc_from_photo_electrons(
 }
 }  // namespace ecal
 
-}  // namespace digitization
+}  // namespace digitization