Jpp  18.1.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JPMTSignalProcessorInterface.hh
Go to the documentation of this file.
1 #ifndef __JDETECTOR__JPMTSIGNALPROCESSORINTERFACE__
2 #define __JDETECTOR__JPMTSIGNALPROCESSORINTERFACE__
3 
4 #include <algorithm>
5 #include <cmath>
6 #include <limits>
7 
10 
11 
12 /**
13  * \author mdejong
14  */
15 
16 namespace JDETECTOR {}
17 namespace JPP { using namespace JDETECTOR; }
18 
19 namespace JDETECTOR {
20 
21  /**
22  * PMT signal processor interface.
23  *
24  * This class supports the implementation of the PMT simulator interface using
25  * an alternative set of virtual methods.
26  * These methods constitute a user interface to the expected performance of a PMT.
27  *
28  * Each photon is converted to a photo-electron using member method getRandomTime().
29  * For this, the data structure JPhotoElectron is used.
30  * Note that the quantum efficiency (QE) of the PMT actually is included in the simulation of the detector response.
31  * A relative QE is applied via method applyQE().
32  * All photo-electrons are then time sorted.
33  * The photo-electrons which simultaneously arrive are merged.
34  * The corresponding condition is defined by member method compare().
35  * The time of the combined signal is determined by the time of the first photo-electron and
36  * the rise time of the analogue pulse to the threshold of the discriminator via method getRiseTime().
37  * In this, the actual amplitude of the combined analogue signal and the calibration of the PMT are taken into account.
38  * The amplitude of the combined analogue signal is simulated using member method getRandomCharge().
39  * For this, the data structure JPMTSignal is used.
40  *
41  * The analogue signals of the PMT are processed by a discriminator.
42  * This discriminator produces a time-over-threshold signal for each analogue signal that passes a preset threshold.
43  * The output signal is described by the time of the leading edge and the length of the time-over-threshold signal.
44  * For this, the data structure JPMTPulse is used.
45  * The determination of the time of the leading edge and the length of the time-over-threshold signal
46  * require a designated model.
47  * The class JPMTAnalogueSignalProcessor provides for an implementation of such a model.
48  *
49  * Overlapping time-over-threshold pulses are merged.
50  * The time of the leading edge is then set to the earliest leading edge and
51  * the time-over-threshold to the difference between
52  * the latest trailing edge and the earliest leading edge.
53  * The merging of hits is implemented in member method merge().
54  *
55  * The default implementation of these virtual methods corresponds to no smearing
56  * and a time-over-threshold value equal to a fixed two photo-electron resolution times the number of photo-electrons.
57  * The width of the charge distribution and the two photo-electron resolution are implemented in methods
58  * getQmin() and getTmin(), respectively.
59  *
60  * For a realistic PMT simulation, the derived class should provide for
61  * an implementation of each of these virtual methods.
62  */
64  public:
65  /**
66  * Default constructor.
67  */
69  {}
70 
71 
72  /**
73  * Virtual destructor.
74  */
76  {}
77 
78 
79  /**
80  * Process hits.
81  *
82  * Two (or more) photo-electrons are combined if they are comparable according method compare.\n
83  * Two (or more) consecutive hits hits maybe merged (according method merge).
84  *
85  * \param calibration PMT calibration
86  * \param input PMT signals
87  * \param output PMT hits
88  */
89  void operator()(const JCalibration& calibration,
90  const JPMTData<JPMTSignal>& input,
91  JPMTData<JPMTPulse>& output) const
92  {
93  // apply transition time distribution to each photo-electron.
94 
95  buffer.clear();
96 
97  for (JPMTData<JPMTSignal>::const_iterator hit = input.begin(); hit != input.end(); ++hit) {
98  for (int i = 0; i != hit->npe; ++i) {
99  if (applyQE()) {
100  buffer.push_back(JPhotoElectron(getRandomTime(hit->t_ns)));
101  }
102  }
103  }
104 
105  if (!buffer.empty()) {
106 
108 
109  buffer.sort();
110 
111 
112  // generate PMT hits from time sequence of photo-electrons.
113 
114  for (JPMTData<JPhotoElectron>::const_iterator q = buffer.begin(), p = q++; q != buffer.end(); ++q) {
115 
116  while (compare(*p,*q)) {
117  ++q;
118  }
119 
120  const double npe = getRandomCharge(distance(p,q));
121 
122  if (applyThreshold(npe)) {
123  output.push_back(JPMTPulse(putTime(p->t_ns + getRiseTime(npe), calibration), getTimeOverThreshold(npe)));
124  }
125 
126  p = q;
127  }
128 
129  // merge overlapping PMT hits.
130 
131  merge(output);
132  }
133  }
134 
135 
136  /**
137  * Apply relative QE.
138  * The default implementation returns <tt>true</tt>.
139  *
140  * \return true if accepted; false if rejected
141  */
142  virtual bool applyQE() const
143  {
144  return true;
145  }
146 
147 
148  /**
149  * Get randomised time according transition time distribution.
150  * The default implementation returns the given value.
151  *
152  * \param t_ns time [ns]
153  * \return time [ns]
154  */
155  virtual double getRandomTime(const double t_ns) const
156  {
157  return t_ns;
158  }
159 
160 
161  /**
162  * Compare arrival times of photo-electrons.
163  * The default implementation uses method getTmin as two photo-electron resolution.
164  *
165  * \param first first photo-electron
166  * \param second second photo-electron
167  * \return true if arrival times of photo-electrons are within two photo-electron resolution; else false
168  */
169  virtual bool compare(const JPhotoElectron& first, const JPhotoElectron& second) const
170  {
171  return second.t_ns - first.t_ns <= getTmin();
172  }
173 
174 
175  /**
176  * Get randomised charge according to gain and gain spread.
177  * The default implementation returns the given value.
178  *
179  * \param NPE number of photo-electrons
180  * \return number of photo-electrons
181  */
182  virtual double getRandomCharge(const int NPE) const
183  {
184  return (double) NPE;
185  }
186 
187 
188  /**
189  * Get probability density for given charge.
190  *
191  * \param npe observed number of photo-electrons
192  * \param NPE true number of photo-electrons
193  * \return probability [npe^-1]
194  */
195  virtual double getChargeProbability(const double npe, const int NPE) const
196  {
197  return (fabs(npe - NPE) <= 0.5 * getQmin() ? 1.0 / getQmin() : 0.0);
198  }
199 
200 
201  /**
202  * Apply threshold.
203  * The default implementation returns <tt>true</tt>.
204  *
205  * \param npe number of photo-electrons
206  * \return true if pass; else false
207  */
208  virtual bool applyThreshold(const double npe) const
209  {
210  return (npe > 0.0);
211  }
212 
213 
214  /**
215  * Get time to reach threshold.
216  * The default implementation returns <tt>0</tt>.
217  *
218  * \param npe number of photo-electrons
219  * \return time [ns]
220  */
221  virtual double getRiseTime(const double npe) const
222  {
223  return 0.0;
224  }
225 
226 
227  /**
228  * Get time-over-threshold (ToT).
229  * The default implementation corresponds to a linear relation between the number of photo-electrons and the time-over-threshold.
230  *
231  * \param npe number of photo-electrons
232  * \return ToT [ns]
233  */
234  virtual double getTimeOverThreshold(const double npe) const
235  {
236  return TIME_OVER_THRESHOLD_NS + getTmin() * (round(npe) - 1.0);
237  }
238 
239 
240  /**
241  * Get derivative of number of photo-electrons to time-over-threshold.
242  *
243  * \param npe number of photo-electrons
244  * \return dnpe/dToT [ns^-1]
245  */
246  virtual double getDerivative(const double npe) const
247  {
248  return 1.0 / getTmin();
249  }
250 
251 
252  /**
253  * Probability that a hit survives the simulation of the PMT.
254  * The default implementation returns <tt>1</tt> if given value larger than <tt>0</tt>.
255  *
256  * \param NPE number of photo-electrons
257  * \return probability
258  */
259  virtual double getSurvivalProbability(const int NPE) const
260  {
261  if (NPE > 0)
262  return 1.0;
263  else
264  return 0.0;
265  }
266 
267 
268  /**
269  * Get number of photo-electrons.
270  * The default implementation corresponds to a linear relation between the number of photo-electrons and the time-over-threshold.
271  *
272  * \param tot_ns time over threshold [ns]
273  * \return number of photo-electrons
274  */
275  virtual double getNPE(const double tot_ns) const
276  {
277  return 1.0 + (tot_ns - TIME_OVER_THRESHOLD_NS) / getTmin();
278  }
279 
280 
281  /**
282  * Merging of PMT hits.
283  *
284  * Hits with overlapping time-over-threshold signals should -de facto- be combined.
285  * In this, the leading edge is maintained and the time-over-threshold is
286  * set to the difference between the overall trailing and leading edges.
287  * As a result, the number of PMT hits may be reduced.
288  *
289  * \param data PMT hits (I/O)
290  */
291  virtual void merge(JPMTData<JPMTPulse>& data) const
292  {
293  using namespace std;
294 
295  JPMTData<JPMTPulse>::iterator out = data.begin();
296 
297  for (JPMTData<JPMTPulse>::iterator i = data.begin(); i != data.end(); ) {
298 
299  double t1 = i->t_ns;
300  double t2 = i->t_ns + i->tot_ns;
301 
302  while (++i != data.end() && i->t_ns < t2 + getTmin()) {
303  t2 = max(t2, i->t_ns + i->tot_ns);
304  }
305 
306  out->t_ns = t1;
307  out->tot_ns = t2 - t1;
308 
309  ++out;
310  }
311 
312  data.resize(distance(data.begin(), out));
313  }
314 
315 
316  /**
317  * Get two photo-electron resolution for time-over-threshold
318  *
319  * \return minimal time [ns]
320  */
321  static double getTmin()
322  {
323  return 1.0;
324  }
325 
326 
327  /**
328  * Get width of charge distribution.
329  *
330  * \return width charge distribution [npe]
331  */
332  static double getQmin()
333  {
334  return 1.0e-3;
335  }
336 
337 
338  private:
340  };
341 
342 
343  /**
344  * Get charge probability.
345  *
346  * \param pmt PMT signal processor
347  * \param npe measured number of photo-electrons
348  * \param NPE expected number of photo-electrons
349  * \param precision precision
350  * \return probability
351  */
353  const double npe,
354  const double NPE,
355  const double precision = 1.0e-4)
356  {
357  int i = (int) (NPE - 5.0 * sqrt(NPE));
358 
359  if (i < 1) {
360  i = 1;
361  }
362 
363  double p = NPE * exp(-NPE) / (double) 1;
364 
365  for (int __i = 1; __i != i; ++__i) {
366  p *= NPE / __i;
367  }
368 
369  double P = 0.0;
370 
371  for (double p0 = 0.0; (p >= p0 || p > precision); ++i, p0 = p, p *= NPE / (double) i) {
372  P += pmt.getChargeProbability(npe, i) * p;
373  }
374 
375  return P;
376  }
377 }
378 
379 #endif
Data structure for PMT digital pulse.
static double getQmin()
Get width of charge distribution.
Time calibration (including definition of sign of time offset).
virtual double getRandomCharge(const int NPE) const
Get randomised charge according to gain and gain spread.
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
Data structure for single photo-electron.
const double TIME_OVER_THRESHOLD_NS
Specification for time-over-threshold corresponding to a one photo-electron pulse.
virtual void merge(JPMTData< JPMTPulse > &data) const
Merging of PMT hits.
Data structure for time calibration.
virtual bool applyQE() const
Apply relative QE.
static double getTmin()
Get two photo-electron resolution for time-over-threshold.
void operator()(const JCalibration &calibration, const JPMTData< JPMTSignal > &input, JPMTData< JPMTPulse > &output) const
Process hits.
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
virtual double getRiseTime(const double npe) const
Get time to reach threshold.
virtual double getSurvivalProbability(const int NPE) const
Probability that a hit survives the simulation of the PMT.
double putTime(const T &t1, const JCalibration &cal)
Get de-calibrated time.
virtual bool applyThreshold(const double npe) const
Apply threshold.
std::vector< JElement_t >::iterator iterator
virtual ~JPMTSignalProcessorInterface()
Virtual destructor.
virtual bool compare(const JPhotoElectron &first, const JPhotoElectron &second) const
Compare arrival times of photo-electrons.
virtual double getNPE(const double tot_ns) const
Get number of photo-electrons.
virtual double getChargeProbability(const double npe, const int NPE) const
Get probability density for given charge.
static JPhotoElectron getEndMarker()
Get end marker.
virtual double getRandomTime(const double t_ns) const
Get randomised time according transition time distribution.
virtual double getTimeOverThreshold(const double npe) const
Get time-over-threshold (ToT).
Template data structure for PMT I/O.
virtual double getDerivative(const double npe) const
Get derivative of number of photo-electrons to time-over-threshold.
std::vector< JElement_t >::const_iterator const_iterator
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable STRING $argv[2] set_array QUANTILES set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable MODULE `getModule -a $DETECTOR -L "$STRING 0"` source JAcousticsToolkit.sh typeset -A TRIPODS get_tripods $WORKDIR/tripod.txt TRIPODS XMEAN
then $DIR JPlotNPE PDG P
Definition: JPlotNPE-PDG.sh:62