BmnRoot
Loading...
Searching...
No Matches
BmnHgndCalibration.h
Go to the documentation of this file.
1#ifndef BmnHgndCalibration_H
2#define BmnHgndCalibration_H
3
4#include <boost/math/special_functions/lambert_w.hpp>
5#include <cmath>
6
7namespace BmnHgndCalibration
8{
9
10struct TdcFromTot
11{
12 double offset, p, tau, RC;
13 double scale;
14
15 TdcFromTot(double offset_, double p_, double tau_, double RC_)
16 : offset(offset_)
17 , p(p_)
18 , tau(tau_)
19 , RC(RC_)
20 {
21 const double e = std::exp(1.0);
22 const double Tp = (2 * e - 1) * p * tau * RC / (e - 1) / (RC * (1 - p) - tau);
23 scale = Tp / 2. / RC * std::sqrt((e - 1) / p) * std::exp(Tp / 2. / RC);
24 }
25
26 double operator()(double tot) const
27 {
28 if (!std::isfinite(scale))
29 return offset;
30 double arg = scale * std::exp(-tot / 2. / RC);
31 double tdc = offset + 2 * RC * boost::math::lambert_w0(arg);
32 return std::isfinite(tdc) ? tdc : offset;
33 }
34};
35
36struct TdcFromTotApprox
37{
38 double offset, p, tau, RC;
39 double scale;
40
41 TdcFromTotApprox(double offset_, double p_, double tau_, double RC_)
42 : offset(offset_)
43 , p(p_)
44 , tau(tau_)
45 , RC(RC_)
46 {
47 const double e = std::exp(1.0);
48 const double Tp = (2 * e - 1) * p * tau * RC / (e - 1) / (RC * (1 - p) - tau);
49 scale = Tp * std::sqrt((e - 1) / p) * std::exp(Tp / 2. / RC);
50 }
51
52 double operator()(double tot) const
53 {
54 if (!std::isfinite(scale))
55 return offset;
56 double tdc = offset + scale * std::exp(-tot / 2. / RC);
57 return std::isfinite(tdc) ? tdc : offset;
58 }
59};
60
61struct QdcFromTot
62{
63 double p, tau, RC;
64 double theta; // threshold in mV
65 double Rload; // load in Ohm
66 double Tp, scale;
67 TdcFromTot tdcRaw; // TDC with zero offset
68
69 QdcFromTot(double p_, double tau_, double RC_, double theta_ = 20.0, double Rload_ = 33.0)
70 : p(p_)
71 , tau(tau_)
72 , RC(RC_)
73 , theta(theta_)
74 , Rload(Rload_)
75 , tdcRaw(0.0, p_, tau_, RC_) // zero-offset TDC for internal computation
76 {
77 const double e = std::exp(1.0);
78 Tp = (2 * e - 1) * p * tau * RC / (e - 1) / (RC * (1 - p) - tau);
79 double como = p * Tp / 2. / (e - 1);
80 scale = theta / Rload * (RC - (1 - p) * tau + como);
81 }
82
83 double operator()(double tot) const
84 {
85 if (!std::isfinite(scale))
86 return 0.0;
87
88 // Reuse TDC logic here
89 const double tdc = tdcRaw(tot);
90 const double qdc = scale * std::exp((tot + tdc - Tp) / RC);
91 return std::isfinite(qdc) ? qdc : 0.0;
92 }
93};
94
95struct QdcFromTotApprox
96{
97 double p, tau, RC;
98 double theta; // threshold in mV
99 double Rload; // load in Ohm
100 double Tp, scale;
101
102 QdcFromTotApprox(double p_, double tau_, double RC_, double theta_ = 20.0, double Rload_ = 33.0)
103 : p(p_)
104 , tau(tau_)
105 , RC(RC_)
106 , theta(theta_)
107 , Rload(Rload_)
108 {
109 const double e = std::exp(1.0);
110 Tp = (2 * e - 1) * p * tau * RC / (e - 1) / (RC * (1 - p) - tau);
111 double lamda = RC + Tp;
112 scale = theta * lamda / Rload; // pC
113 }
114
115 double operator()(double tot) const
116 {
117 if (!std::isfinite(scale))
118 return 0.0;
119
120 const double qdc = scale * (1 + exp((tot) / (RC + tau)));
121 return std::isfinite(qdc) ? qdc : 0.0;
122 }
123};
124
125} // namespace BmnHgndCalibration
126
127#endif // BmnHgndCalibration_H
friend F32vec4 exp(const F32vec4 &a)
Definition P4_F32vec4.h:122
QdcFromTotApprox(double p_, double tau_, double RC_, double theta_=20.0, double Rload_=33.0)
double operator()(double tot) const
QdcFromTot(double p_, double tau_, double RC_, double theta_=20.0, double Rload_=33.0)
TdcFromTotApprox(double offset_, double p_, double tau_, double RC_)
double operator()(double tot) const
TdcFromTot(double offset_, double p_, double tau_, double RC_)