From caebbbc0cda17d4fcd6fe13969645b9f51a5edba Mon Sep 17 00:00:00 2001
From: Jack Christopher Hutchinson Rolph <jack.rolph@desy.de>
Date: Mon, 22 May 2023 14:23:56 +0200
Subject: [PATCH] Update 2 files

- /PeakOTron.py
- /Model.py
---
 Model.py     | 30 +++++++++---------------------
 PeakOTron.py | 29 +++++++++++++++++------------
 2 files changed, 26 insertions(+), 33 deletions(-)

diff --git a/Model.py b/Model.py
index f4522ca..15bf6c4 100644
--- a/Model.py
+++ b/Model.py
@@ -165,36 +165,24 @@ def Multichoose(n,k):
 def Beta(tau_R, tau_Ap, t_gate):
     return (tau_Ap - np.exp(-t_gate/tau_Ap)*(tau_Ap + tau_R*(1 - np.exp(-t_gate/tau_R))))/(tau_Ap + tau_R)
 
+
+
+
 ########################################################################################################
-## P_Ap(tau_R, tau_Ap, t_gate):
+## Beta_Inf(tau_R, tau_Ap, t_gate):
 ########################################################################################################
 #
-# Calculates afterpulse discharge probability from scaling factor p_Ap and maximum discharge scaling factor Beta.
+# Normalisation constant for the afterpulse p.d.f to infinity
+# 
+# Beta_Inf(tau_R, tau_Ap, t_gate)  = ∫(1 - exp(-t/tauR)) x exp(-t/tauAp)/tauAp from 0 to inf = tau_Ap/(tau_R + tau_Ap)
 #
-# Afterpulses can only have a maximum 
 ########################################################################################################
 @njit
-def P_Ap(tau_R, tau_Ap, t_gate, p_Ap):
-    return p_Ap*Beta(tau_R, tau_Ap, t_gate)
+def Beta_Inf(tau_R, tau_Ap):
+    return tau_Ap/(tau_Ap + tau_R)
 
 
-########################################################################################################
-## dPAp(tau_R, tau_Ap, t_gate, p_Ap)
-########################################################################################################
-##
-## Returns error on total afterpulse probability, dPAp, occurring between 0 and t_gate.
-##
-########################################################################################################
 
-@njit
-def dP_Ap(tau_R, tau_Ap, t_gate, p_Ap, dtau_Ap, dp_Ap):
-    dPApdtauAp = p_Ap*((t_gate*(-tau_Ap - tau_R*(1 - np.exp(-t_gate/tau_R)))*np.exp(-t_gate/tau_Ap))/tau_Ap**2 + 1 - np.exp(-t_gate/tau_Ap))/(tau_Ap+tau_R) - p_Ap*(tau_Ap - (tau_Ap + tau_R*(1 - np.exp(-t_gate/tau_R)))*np.exp(-t_gate/tau_Ap))/(tau_Ap + tau_R)**2
-    
-    dPApdpAp = (tau_Ap - (tau_Ap + tau_R*(1 - np.exp(-t_gate/tau_R)))*np.exp(-t_gate/tau_Ap))/(tau_Ap + tau_R)
-    
-    dPAp =  np.sqrt((dPApdtauAp**2)*(dtau_Ap**2) + (dPApdpAp**2)*(dp_Ap**2))
-    
-    return dPAp
 
 
 ########################################################################################################
diff --git a/PeakOTron.py b/PeakOTron.py
index 53a2f36..6bb5d03 100644
--- a/PeakOTron.py
+++ b/PeakOTron.py
@@ -2,7 +2,7 @@ from HelperFunctions import GP_gain, GP_lbda, GP_mu, HistMean, HistStd, HistSkew
 from HelperFunctions import LatexFormat, Linear
 from AdditionalPDFs import gpoisson
 from LossFunctions import *
-from Model import DRM, P_Ap, dP_Ap, sigmak, Beta
+from Model import DRM, sigmak, Beta, Beta_Inf
 from iminuit import Minuit
 from scipy.interpolate import interp1d, UnivariateSpline, InterpolatedUnivariateSpline
 from scipy.stats import poisson, norm
@@ -274,7 +274,7 @@ class PeakOTron:
     
     
     
-    def GetFitResults(self, debug=True, bin_units=True):
+    def GetFitResults(self, bin_units=True):
         
 
         
@@ -289,17 +289,22 @@ class PeakOTron:
             else:
                 temp_values = deepcopy(self._fit_values)
                 temp_errors = deepcopy(self._fit_errors)
-            
-            if(not debug):
 
-                temp_values["p_Ap"] = P_Ap(self._fit_values["tau_R"], self._fit_values["tau_Ap"], self._fit_values["t_gate"], self._fit_values["p_Ap"])
-                temp_errors["p_Ap"] = dP_Ap(self._fit_values["tau_R"], self._fit_values["tau_Ap"], self._fit_values["t_gate"], self._fit_values["p_Ap"], self._fit_errors["tau_Ap"], self._fit_errors["p_Ap"])
-                
-            else:
-                   
-                temp_values["P_Ap"] = P_Ap(self._fit_values["tau_R"], self._fit_values["tau_Ap"], self._fit_values["t_gate"], self._fit_values["p_Ap"])
-                temp_errors["P_Ap"] = dP_Ap(self._fit_values["tau_R"], self._fit_values["tau_Ap"], self._fit_values["t_gate"], self._fit_values["p_Ap"], self._fit_errors["tau_Ap"], self._fit_errors["p_Ap"])
-                
+
+            ## This part rescales the afterpulse to an infinite gate length. 
+            ## Please note that the effectiveness of the scaling factor depends on the ability to determine tau_Ap. 
+            ## This may be challenging for some spectra where there is insufficient information between the peaks to calculate tau_Ap. 
+            ## In that case, it is suggested that an appropriate tau_Ap is assumed (i.e. from overvoltages allowing an accurate assessment of tau_Ap)
+            ## so that this factor can be reconstructed. 
+
+            p_Ap_inf_sf = Beta_Inf(temp_values["tau_R"], temp_values["tau_Ap"])/Beta(temp_values["tau_R"], temp_values["tau_Ap"], temp_values["t_gate"])
+
+
+            temp_values["p_Ap_inf"] = temp_values["p_Ap"]*p_Ap_inf_sf
+            temp_errors["p_Ap_inf"] =  temp_errors["p_Ap"]*p_Ap_inf_sf ## TO DO: propogate errors correctly. 
+
+
+
             t_fit = self.GetFitTime(prefit=False)
             t_prefit = self.GetFitTime(prefit=True)
                 
-- 
GitLab