Jpp
JPMTAnalogueSignalProcessor.hh
Go to the documentation of this file.
1 #ifndef __JDETECTOR__JPMTANALOGUESIGNALPROCESSOR__
2 #define __JDETECTOR__JPMTANALOGUESIGNALPROCESSOR__
3 
4 #include <istream>
5 #include <cmath>
6 #include <limits>
7 
8 #include "JLang/JException.hh"
12 
13 /**
14  * \file
15  *
16  * PMT analogue signal processor.
17  * \author mdejong
18  */
19 namespace JDETECTOR {
20 
22 
23 
24  /**
25  * PMT analogue signal processor.
26  *
27  * This class provides for an implementation of the JDETECTOR::JPMTSignalProcessorInterface
28  * using a specific model for the analogue pulse of the PMT.\n
29  * In this, the leading edge of the analogue pulse from the PMT is assumed to be a Gaussian and the tail an exponential.\n
30  * The width of the Gaussian is referred to as the rise time and
31  * the inverse slope of the exponential to the decay time.\n
32  * The two functions are matched at a point where the values and first derivatives are identical.\n
33  *
34  * Note that the decay time is related to the rise time via the specification JDETECTOR::TIME_OVER_THRESHOLD_NS.
35  *
36  * The charge distribution is assumed to be a Gaussian which is centered at the specified gain
37  * and truncated by the specified threshold.
38  *
39  * The transition times are generated according the specified spread as follows.\n
40  * - If the specified transition time spread (TTS) is negative,
41  * the transition times are generated according to the measurements (see method JDETECTOR::getTransitionTime).\n
42  * - If the specified TTS is positive,
43  * the transition times are generated according a Gaussian with a sigma equals to the TTS.\n
44  * - If the TTS is zero, the transition times are generated without any spread.
45  */
48  public JPMTParameters
49  {
50  /**
51  * Constructor.
52  *
53  * \param parameters PMT parameters
54  */
57  JPMTParameters(parameters),
58  decayTime_ns(0.0),
59  t1(0.0),
60  y1(0.0),
61  x1(std::numeric_limits<double>::max())
62  {
63  configure();
64  }
65 
66 
67  /**
68  * Configure internal parameters.
69  *
70  * This method provides the implementations for
71  * - matching of the leading edge of the analogue pulse (Gaussian) and the tail (exponential); and
72  * - determination of number of photo-electrons above which the time-over-threshold
73  * linearly depends on the number of photo-electrons (apart from saturation).
74  *
75  * Note that this method will throw an error if the value of the rise time (i.e. width of the Gaussian)
76  * is too large with respect to the specification JDETECTOR::TIME_OVER_THRESHOLD_NS.
77  */
78  void configure()
79  {
80  static const int N = 100;
81  static const double precision = 1.0e-3;
82 
83  // check rise time
84 
86  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::configure(): Invalid rise time [ns] " << riseTime_ns);
87  }
88 
89  // decay time
90 
91  const double y = -log(threshold);
92 
93  const double a = y;
94  const double b = riseTime_ns * sqrt(2.0*y) - TIME_OVER_THRESHOLD_NS;
95  const double c = 0.5*riseTime_ns*riseTime_ns;
96 
97  const double Q = b*b - 4.0*a*c;
98 
99  if (Q > 0.0)
100  decayTime_ns = (-b + sqrt(Q)) / (2.0*a);
101  else
102  decayTime_ns = -b / (2.0*a);
103 
104  // fix matching of Gaussian and exponential
105 
106  const double x = riseTime_ns / decayTime_ns;
107 
108  t1 = riseTime_ns*x;
109  y1 = exp(-0.5*x*x);
110 
111  // determine transition point to linear dependence of time-over-threshold as a function of number of photo-electrons
112 
113  const double v = slope / getDerivativeOfSaturation(getTimeOverThreshold(1.0));
114 
115  x1 = std::numeric_limits<double>::max(); // disable linearisation
116 
117  double xmin = 1.0;
118  double xmax = 1.0 / (getDerivative(1.0) * slope);
119 
120  for (int i = 0; i != N; ++i) {
121 
122  const double x = 0.5 * (xmin + xmax);
123  const double u = getDerivative(x) * v;
124 
125  if (fabs(1.0 - u) < precision) {
126  break;
127  }
128 
129  if (u < 1.0)
130  xmin = x;
131  else
132  xmax = x;
133  }
134 
135  x1 = 0.5 * (xmin + xmax);
136  }
137 
138 
139  /**
140  * Get decay time.
141  *
142  * \return decay time [ns]
143  */
144  double getDecayTime() const
145  {
146  return decayTime_ns;
147  }
148 
149 
150  /**
151  * Get time at transition point from Gaussian to exponential.
152  *
153  * \return time [ns]
154  */
155  double getT1() const
156  {
157  return t1;
158  }
159 
160 
161  /**
162  * Get amplitude at transition point from Gaussian to exponential.
163  *
164  * \return amplitude [npe]
165  */
166  double getY1() const
167  {
168  return y1;
169  }
170 
171 
172  /**
173  * Get transition point from a model dependent to linear relation between time-over-threshold and number of photo-electrons.
174  *
175  * \return number of photo-electrons [npe]
176  */
177  double getStartOfLinearisation() const
178  {
179  return x1;
180  }
181 
182 
183  /**
184  * Get amplitude at given time for a one photo-electron pulse.
185  *
186  * \param t1_ns time [ns]
187  * \return amplitude [npe]
188  */
189  double getAmplitude(const double t1_ns) const
190  {
191  if (t1_ns < t1) {
192 
193  const double x = t1_ns / riseTime_ns;
194 
195  return exp(-0.5*x*x); // Gaussian
196 
197  } else {
198 
199  const double x = t1_ns / decayTime_ns;
200 
201  return exp(-x) / y1; // exponential
202  }
203  }
204 
205 
206  /**
207  * Get time to pass from threshold to top of analogue pulse.\n
208  * In this, the leading edge of the analogue pulse is assumed to be Gaussian.
209  *
210  * \param npe number of photo-electrons
211  * \param th threshold [npe]
212  * \return time [ns]
213  */
214  double getRiseTime(const double npe, const double th) const
215  {
216  return riseTime_ns * sqrt(2.0*log(npe/th)); // Gaussian
217  }
218 
219 
220  /**
221  * Get time to pass from top of analogue pulse to threshold.\n
222  * In this, the trailing edge of the analogue pulse is assumed to be exponential.
223  *
224  * \param npe number of photo-electrons
225  * \param th threshold [npe]
226  * \return time [ns]
227  */
228  double getDecayTime(const double npe, const double th) const
229  {
230  if (npe*y1 > th)
231  return decayTime_ns * (log(npe/th) - log(y1)); // exponential
232  else
233  return riseTime_ns * sqrt(2.0*log(npe/th)); // Gaussian
234  }
235 
236 
237  /**
238  * Get maximal rise time for given threshold.
239  *
240  * Note that the rise time is entirely constraint by the specification JDETECTOR::TIME_OVER_THRESHOLD_NS.
241  *
242  * \param th threshold [npe]
243  * \return rise time [ns]
244  */
245  static double getMaximalRiseTime(const double th)
246  {
247  if (th > 0.0 && th < 1.0)
248  return 0.5 * TIME_OVER_THRESHOLD_NS / sqrt(-2.0*log(th));
249  else
250  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getMaximalRiseTime(): Invalid threshold " << th);
251  }
252 
253 
254  /**
255  * Get function value with saturation.
256  *
257  * \param y value
258  * \return value
259  */
260  double getSaturation(const double y) const
261  {
262  return saturation / sqrt(y*y + saturation*saturation);
263  }
264 
265 
266  /**
267  * Get derivative of function value with saturation.
268  *
269  * \param y value
270  * \return value
271  */
272  double getDerivativeOfSaturation(const double y) const
273  {
274  return (saturation*saturation + y*y) * sqrt(saturation*saturation + y*y) / (saturation * saturation * saturation);
275  }
276 
277 
278  /**
279  * Get gain spread for given number of photo-electrons.
280  *
281  * \param NPE number of photo-electrons
282  * \return gain spread
283  */
284  double getGainSpread(int NPE) const
285  {
286  return sqrt((double) NPE * gain) * gainSpread;
287  }
288 
289 
290  /**
291  * Get integral of probability.
292  *
293  * \param xmin minimum number of photo-electrons
294  * \param xmax maximum number of photo-electrons
295  * \param NPE true number of photo-electrons
296  * \return probability
297  */
298  double getIntegralOfProbability(const double xmin, const double xmax, const int NPE) const
299  {
300  double zmin = xmin;
301  double zmax = xmax;
302 
303  if (zmin < threshold) { zmin = threshold; }
304  if (zmax < threshold) { zmax = threshold; }
305 
306  const double mu = NPE * gain;
307  const double sigma = getGainSpread(NPE);
308 
309  return (0.5 * TMath::Erfc((zmin - mu) / sqrt(2.0) / sigma) -
310  0.5 * TMath::Erfc((zmax - mu) / sqrt(2.0) / sigma));
311  }
312 
313 
314  /**
315  * Set PMT parameters.
316  *
317  * \param parameters PMT parameters
318  */
319  void setPMTParameters(const JPMTParameters& parameters)
320  {
321  static_cast<JPMTParameters&>(*this).setPMTParameters(parameters);
322 
323  configure();
324  }
325 
326 
327  /**
328  * Read PMT signal from input.
329  *
330  * \param in input stream
331  * \param object PMT signal
332  * \return input stream
333  */
334  friend std::istream& operator>>(std::istream& in, JPMTAnalogueSignalProcessor& object)
335  {
336  in >> static_cast<JPMTParameters&>(object);
337 
338  object.configure();
339 
340  return in;
341  }
342 
343 
344  /**
345  * Apply relative QE.
346  *
347  * \return true if accepted; false if rejected
348  */
349  virtual bool applyQE() const
350  {
351  if (QE <= 0.0)
352  return false;
353  else if (QE < 1.0)
354  return gRandom->Rndm() < QE;
355  else
356  return true;
357  }
358 
359 
360  /**
361  * Get randomised time according transition time distribution.
362  *
363  * \param t_ns time [ns]
364  * \return time [ns]
365  */
366  virtual double getRandomTime(const double t_ns) const
367  {
368  if (TTS_ns < 0.0)
369  return t_ns + getTransitionTime(gRandom->Rndm());
370  else
371  return gRandom->Gaus(t_ns, TTS_ns);
372  }
373 
374 
375  /**
376  * Compare (arrival times of) photo-electrons.
377  *
378  * \param first first photo-electron
379  * \param second second photo-electron
380  * \return true if arrival times of photo-electrons are within rise time of analogue pulses; else false
381  */
382  virtual bool compare(const JPhotoElectron& first, const JPhotoElectron& second) const
383  {
384  return second.t_ns < first.t_ns + riseTime_ns;
385  }
386 
387 
388  /**
389  * Get randomised amplitude according gain and gain spread.
390  *
391  * \param NPE number of photo-electrons
392  * \return number of photo-electrons
393  */
394  virtual double getRandomAmplitude(const int NPE) const
395  {
396  double x;
397 
398  do {
399  x = gRandom->Gaus(NPE*gain, getGainSpread(NPE));
400  } while (x < 0.0);
401 
402  return x;
403  }
404 
405 
406  /**
407  * Get probability for given charge.
408  *
409  * \param npe observed number of photo-electrons
410  * \param NPE true number of photo-electrons
411  * \return probability [npe^-1]
412  */
413  virtual double getProbability(const double npe, const int NPE) const
414  {
415  if (npe >= threshold) {
416 
417  const double mu = NPE * gain;
418  const double sigma = getGainSpread(NPE);
419  const double V = 0.5 * TMath::Erfc((threshold - mu) / sqrt(2.0) / sigma);
420 
421  return TMath::Gaus(npe, mu, sigma, kTRUE) / V;
422  }
423 
424  return 0.0;
425  }
426 
427 
428  /**
429  * Apply threshold.
430  *
431  * \param npe number of photo-electrons
432  * \return true if pass; else false
433  */
434  virtual bool applyThreshold(const double npe) const
435  {
436  return npe >= threshold;
437  }
438 
439 
440  /**
441  * Get time to reach threshold.
442  *
443  * Note that the rise time is defined to be zero for a one photo-electron signal.
444  *
445  * \param npe number of photo-electrons
446  * \return time [ns]
447  */
448  virtual double getRiseTime(const double npe) const
449  {
450  if (slewing)
451  return ((getRiseTime(npe, getTH0()) - getRiseTime(npe, this->threshold)) -
452  (getRiseTime(1.0, getTH0()) - getRiseTime(1.0, this->threshold)));
453  else
454  return 0.0;
455  }
456 
457 
458  /**
459  * Get time-over-threshold (ToT).
460  *
461  * \param npe number of photo-electrons
462  * \return ToT [ns]
463  */
464  virtual double getTimeOverThreshold(const double npe) const
465  {
466  if (npe > threshold) {
467 
468  double tot = 0.0;
469 
470  if (npe*y1 > threshold) {
471 
472  if (npe <= getStartOfLinearisation()) {
473 
474  tot += getRiseTime (npe, threshold); // Gaussian
475  tot += getDecayTime(npe, threshold); // exponential
476 
477  } else {
478 
479  tot += getRiseTime (getStartOfLinearisation(), threshold); // Gaussian
480  tot += getDecayTime(getStartOfLinearisation(), threshold); // exponential
481 
482  tot += slope * (npe - getStartOfLinearisation()); // linear
483  }
484 
485  } else {
486 
487  tot += getRiseTime(npe, threshold); // Gaussian
488  tot += getRiseTime(npe, threshold); // Gaussian
489  }
490 
491  return tot * getSaturation(tot);
492 
493  } else {
494 
495  return 0.0;
496  }
497  }
498 
499 
500  /**
501  * Get derivative of number of photo-electrons to time-over-threshold.
502  *
503  * \param npe number of photo-electrons
504  * \return dnpe/dToT [ns^-1]
505  */
506  virtual double getDerivative(const double npe) const
507  {
508  if (npe >= threshold) {
509 
510  const double z = riseTime_ns / sqrt(2.0 * log(npe/threshold));
511 
512  double y = 0.0;
513 
514  if (npe*y1 > threshold) {
515 
516  if (npe <= getStartOfLinearisation())
517  y = npe / (z + decayTime_ns); // Gaussian + exponential
518  else
519  y = 1.0 / slope; // linear
520 
521  } else {
522 
523  y = npe / (2.0 * z); // Gaussian + Gaussian
524  }
525 
527 
528  } else {
529 
530  return 0.0;
531  }
532  }
533 
534 
535  /**
536  * Probability that a hit survives the simulation of the PMT.
537  *
538  * \param NPE number of photo-electrons
539  * \return probability
540  */
541  virtual double getSurvivalProbability(const int NPE) const
542  {
543  if (QE <= 0.0) {
544 
545  return 0.0;
546 
547  } else if (QE <= 1.0) {
548 
549  double P = 0.0;
550 
551  for (int i = 1; i <= NPE; ++i) {
552 
553  const double p = TMath::Binomial(NPE, i) * TMath::Power(QE, i) * TMath::Power(1.0 - QE, NPE - i);
554  const double mu = i * gain;
555  const double sigma = getGainSpread(i);
556  const double V = 0.5 * TMath::Erfc((threshold - mu) / sqrt(2.0) / sigma);
557  const double W = 0.5 * TMath::Erfc((0.0 - mu) / sqrt(2.0) / sigma);
558 
559  P += p*V/W;
560  }
561 
562  return P;
563 
564  } else {
565 
566  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getSurvivalProbability: Invalid QE " << QE);
567  }
568  }
569 
570 
571  /**
572  * Get lower threshold for rise time evaluation.
573  *
574  * \return threshold [npe]
575  */
576  static double getTH0()
577  {
578  return 0.1;
579  }
580 
581 
582  /**
583  * Get upper threshold for rise time evaluation.
584  *
585  * \return threshold [npe]
586  */
587  static double getTH1()
588  {
589  return 0.9;
590  }
591 
592  protected:
593  double decayTime_ns; //!< decay time [ns]
594  double t1; //!< time at match point [ns]
595  double y1; //!< amplitude at match point [npe]
596  /**
597  * Transition point from a logarithmic to a linear relation
598  * between time-over-threshold and number of photo-electrons.\n
599  * Measurements by B. Schermer and R. Bruijn at Nikhef.
600  */
601  double x1;
602  };
603 }
604 
605 #endif
JException.hh
JDETECTOR::JPMTAnalogueSignalProcessor::getIntegralOfProbability
double getIntegralOfProbability(const double xmin, const double xmax, const int NPE) const
Get integral of probability.
Definition: JPMTAnalogueSignalProcessor.hh:298
JDETECTOR::JPhotoElectron::t_ns
double t_ns
time [ns]
Definition: JPMTSimulator.hh:56
JDETECTOR::JPMTAnalogueSignalProcessor::t1
double t1
time at match point [ns]
Definition: JPMTAnalogueSignalProcessor.hh:594
JDETECTOR::JPMTParameters::saturation
double saturation
saturation [ns]
Definition: JPMTParameters.hh:223
JDETECTOR::JPMTAnalogueSignalProcessor::getSurvivalProbability
virtual double getSurvivalProbability(const int NPE) const
Probability that a hit survives the simulation of the PMT.
Definition: JPMTAnalogueSignalProcessor.hh:541
JDETECTOR::JPMTAnalogueSignalProcessor::x1
double x1
Transition point from a logarithmic to a linear relation between time-over-threshold and number of ph...
Definition: JPMTAnalogueSignalProcessor.hh:601
JDETECTOR::JPMTParameters::threshold
double threshold
threshold [npe]
Definition: JPMTParameters.hh:221
JDETECTOR::JPMTAnalogueSignalProcessor::getTimeOverThreshold
virtual double getTimeOverThreshold(const double npe) const
Get time-over-threshold (ToT).
Definition: JPMTAnalogueSignalProcessor.hh:464
JDETECTOR::JPMTAnalogueSignalProcessor::getAmplitude
double getAmplitude(const double t1_ns) const
Get amplitude at given time for a one photo-electron pulse.
Definition: JPMTAnalogueSignalProcessor.hh:189
JDETECTOR::TIME_OVER_THRESHOLD_NS
const double TIME_OVER_THRESHOLD_NS
Specification for time-over-threshold corresponding to a one photo-electron pulse.
Definition: JDetector/JCalibration.hh:29
JDETECTOR::JPMTAnalogueSignalProcessor::getY1
double getY1() const
Get amplitude at transition point from Gaussian to exponential.
Definition: JPMTAnalogueSignalProcessor.hh:166
JDETECTOR::JPMTParameters::TTS_ns
double TTS_ns
transition time spread [ns]
Definition: JPMTParameters.hh:220
JPMTSignalProcessorInterface.hh
JDETECTOR::getTransitionTime
static const JTransitionTimeGenerator getTransitionTime
Function object to generate transition time.
Definition: JPMTTransitTimeGenerator.hh:1071
JDETECTOR::JPMTAnalogueSignalProcessor::getDerivativeOfSaturation
double getDerivativeOfSaturation(const double y) const
Get derivative of function value with saturation.
Definition: JPMTAnalogueSignalProcessor.hh:272
JTOOLS::u
double u[N+1]
Definition: JPolint.hh:706
JDETECTOR::JPMTAnalogueSignalProcessor::decayTime_ns
double decayTime_ns
decay time [ns]
Definition: JPMTAnalogueSignalProcessor.hh:593
JDETECTOR::JPMTAnalogueSignalProcessor::getTH0
static double getTH0()
Get lower threshold for rise time evaluation.
Definition: JPMTAnalogueSignalProcessor.hh:576
JDETECTOR::JPMTAnalogueSignalProcessor::setPMTParameters
void setPMTParameters(const JPMTParameters &parameters)
Set PMT parameters.
Definition: JPMTAnalogueSignalProcessor.hh:319
JDETECTOR::JPMTAnalogueSignalProcessor::getRandomAmplitude
virtual double getRandomAmplitude(const int NPE) const
Get randomised amplitude according gain and gain spread.
Definition: JPMTAnalogueSignalProcessor.hh:394
JDETECTOR::JPMTAnalogueSignalProcessor::getGainSpread
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.
Definition: JPMTAnalogueSignalProcessor.hh:284
JDETECTOR::JPMTAnalogueSignalProcessor::getDecayTime
double getDecayTime(const double npe, const double th) const
Get time to pass from top of analogue pulse to threshold.
Definition: JPMTAnalogueSignalProcessor.hh:228
JDETECTOR::JPMTParameters::gainSpread
double gainSpread
gain spread [unit]
Definition: JPMTParameters.hh:218
JDETECTOR::JPMTAnalogueSignalProcessor::getStartOfLinearisation
double getStartOfLinearisation() const
Get transition point from a model dependent to linear relation between time-over-threshold and number...
Definition: JPMTAnalogueSignalProcessor.hh:177
JDETECTOR::JPhotoElectron
Data structure for single photo-electron.
Definition: JPMTSimulator.hh:26
JDETECTOR::JPMTParameters::slewing
bool slewing
time slewing of analogue signal
Definition: JPMTParameters.hh:224
JDETECTOR::JPMTAnalogueSignalProcessor::configure
void configure()
Configure internal parameters.
Definition: JPMTAnalogueSignalProcessor.hh:78
JDETECTOR::JPMTParameters::gain
double gain
gain [unit]
Definition: JPMTParameters.hh:217
JDETECTOR::JPMTAnalogueSignalProcessor::JPMTAnalogueSignalProcessor
JPMTAnalogueSignalProcessor(const JPMTParameters &parameters=JPMTParameters())
Constructor.
Definition: JPMTAnalogueSignalProcessor.hh:55
JDETECTOR::JPMTAnalogueSignalProcessor::getMaximalRiseTime
static double getMaximalRiseTime(const double th)
Get maximal rise time for given threshold.
Definition: JPMTAnalogueSignalProcessor.hh:245
JDETECTOR::JPMTAnalogueSignalProcessor::getT1
double getT1() const
Get time at transition point from Gaussian to exponential.
Definition: JPMTAnalogueSignalProcessor.hh:155
JDETECTOR::JPMTAnalogueSignalProcessor::operator>>
friend std::istream & operator>>(std::istream &in, JPMTAnalogueSignalProcessor &object)
Read PMT signal from input.
Definition: JPMTAnalogueSignalProcessor.hh:334
THROW
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:669
JDETECTOR::JPMTAnalogueSignalProcessor::applyThreshold
virtual bool applyThreshold(const double npe) const
Apply threshold.
Definition: JPMTAnalogueSignalProcessor.hh:434
JLANG::JValueOutOfRange
Exception for accessing a value in a collection that is outside of its range.
Definition: JException.hh:162
JDETECTOR::JPMTAnalogueSignalProcessor::getRiseTime
double getRiseTime(const double npe, const double th) const
Get time to pass from threshold to top of analogue pulse.
Definition: JPMTAnalogueSignalProcessor.hh:214
JDETECTOR::JPMTAnalogueSignalProcessor::getDecayTime
double getDecayTime() const
Get decay time.
Definition: JPMTAnalogueSignalProcessor.hh:144
JCalibration.hh
JDETECTOR::JPMTSignalProcessorInterface
PMT signal processor interface.
Definition: JPMTSignalProcessorInterface.hh:69
JDETECTOR::JPMTAnalogueSignalProcessor::getProbability
virtual double getProbability(const double npe, const int NPE) const
Get probability for given charge.
Definition: JPMTAnalogueSignalProcessor.hh:413
JTOOLS::v
data_type v[N+1][M+1]
Definition: JPolint.hh:707
JDETECTOR::JPMTAnalogueSignalProcessor::getSaturation
double getSaturation(const double y) const
Get function value with saturation.
Definition: JPMTAnalogueSignalProcessor.hh:260
JDETECTOR::JPMTAnalogueSignalProcessor::compare
virtual bool compare(const JPhotoElectron &first, const JPhotoElectron &second) const
Compare (arrival times of) photo-electrons.
Definition: JPMTAnalogueSignalProcessor.hh:382
JDETECTOR::JPMTAnalogueSignalProcessor::applyQE
virtual bool applyQE() const
Apply relative QE.
Definition: JPMTAnalogueSignalProcessor.hh:349
JDETECTOR::JPMTAnalogueSignalProcessor::getRandomTime
virtual double getRandomTime(const double t_ns) const
Get randomised time according transition time distribution.
Definition: JPMTAnalogueSignalProcessor.hh:366
std
Definition: jaanetDictionary.h:36
JDETECTOR::JPMTAnalogueSignalProcessor::getRiseTime
virtual double getRiseTime(const double npe) const
Get time to reach threshold.
Definition: JPMTAnalogueSignalProcessor.hh:448
JDETECTOR::JPMTParameters::slope
double slope
slope [ns/npe]
Definition: JPMTParameters.hh:222
JDETECTOR::JPMTAnalogueSignalProcessor::getDerivative
virtual double getDerivative(const double npe) const
Get derivative of number of photo-electrons to time-over-threshold.
Definition: JPMTAnalogueSignalProcessor.hh:506
JDETECTOR::JPMTAnalogueSignalProcessor::getTH1
static double getTH1()
Get upper threshold for rise time evaluation.
Definition: JPMTAnalogueSignalProcessor.hh:587
JDETECTOR::JPMTAnalogueSignalProcessor::y1
double y1
amplitude at match point [npe]
Definition: JPMTAnalogueSignalProcessor.hh:595
JDETECTOR::JPMTParameters::QE
double QE
relative quantum efficiency
Definition: JPMTParameters.hh:216
JDETECTOR::JPMTParameters
Data structure for PMT parameters.
Definition: JPMTParameters.hh:29
JPMTParameters.hh
JDETECTOR
Auxiliary classes and methods for detector calibration.
Definition: JAnchor.hh:12
JDETECTOR::JPMTParameters::riseTime_ns
double riseTime_ns
rise time of analogue pulse [ns]
Definition: JPMTParameters.hh:219
JDETECTOR::JPMTParameters::JPMTParameters
JPMTParameters()
Default constructor.
Definition: JPMTParameters.hh:41
JDETECTOR::JPMTAnalogueSignalProcessor
PMT analogue signal processor.
Definition: JPMTAnalogueSignalProcessor.hh:46