Jpp test-rotations-new
the software that should make you happy
Loading...
Searching...
No Matches
JDETECTOR::JPMTAnalogueSignalProcessor Struct Reference

PMT analogue signal processor. More...

#include <JPMTAnalogueSignalProcessor.hh>

Inheritance diagram for JDETECTOR::JPMTAnalogueSignalProcessor:
JDETECTOR::JPMTSignalProcessorInterface JDETECTOR::JPMTParameters

Public Types

enum  JThresholdDomain { BELOW_THRESHOLD = -1 , THRESHOLDBAND = 0 , ABOVE_THRESHOLD = 2 }
 Threshold domain specifiers. More...
 

Public Member Functions

 JPMTAnalogueSignalProcessor (const JPMTParameters &parameters=JPMTParameters())
 Constructor.
 
void configure ()
 Configure internal parameters.
 
double getDecayTime () const
 Get decay time.
 
double getT1 () const
 Get time at transition point from Gaussian to exponential.
 
double getY1 () const
 Get amplitude at transition point from Gaussian to exponential.
 
double getStartOfLinearisation () const
 Get transition point from a model-dependent to linear relation between time-over-threshold and number of photo-electrons.
 
double getAmplitude (const double t1_ns) const
 Get amplitude at given time for a one photo-electron pulse.
 
double getRiseTime (const double npe, const double th) const
 Get time to pass from threshold to top of analogue pulse.
 
double getDecayTime (const double npe, const double th) const
 Get time to pass from top of analogue pulse to threshold.
 
double applySaturation (const double tot_ns) const
 Get time-over-threshold with saturation.
 
double removeSaturation (const double tot_ns) const
 Get time-over-threshold without saturation.
 
double getDerivativeOfSaturation (const double tot_ns) const
 Get derivative of saturation factor.
 
double getGainSpread (int NPE) const
 Get gain spread for given number of photo-electrons.
 
double getIntegralOfChargeProbability (const double xmin, const double xmax, const int NPE) const
 Get integral of probability.
 
double getIntegralOfChargeProbability (const JThresholdDomain domain, const int NPE) const
 Get integral of probability in specific threshold domain.
 
void setPMTParameters (const JPMTParameters &parameters)
 Set PMT parameters.
 
JThresholdDomain getThresholdDomain (const double npe) const
 Get threshold domain.
 
virtual bool applyQE () const override
 Apply relative QE.
 
virtual double getRandomTime (const double t_ns) const override
 Get randomised time according transition time distribution.
 
virtual bool compare (const JPhotoElectron &first, const JPhotoElectron &second) const override
 Compare arrival times of photo-electrons.
 
virtual double getRandomCharge (const int NPE) const override
 Get randomised charge according to gain and gain spread.
 
virtual double getChargeProbability (const double npe, const int NPE) const override
 Get probability density for given charge.
 
virtual bool applyThreshold (const double npe) const override
 Apply threshold.
 
virtual double getRiseTime (const double npe) const override
 Get time to reach threshold.
 
virtual double getTimeOverThreshold (const double npe) const override
 Get time-over-threshold (ToT).
 
virtual double getDerivative (const double npe) const override
 Get derivative of number of photo-electrons to time-over-threshold.
 
virtual double getSurvivalProbability (const int NPE) const override
 Probability that a hit survives the simulation of the PMT.
 
virtual double getNPE (const double tot_ns) const override
 Get number of photo-electrons.
 
double getTimeOverThresholdProbability (const double tot_ns, const int NPE) const
 Get probability of having a pulse with specific time-over-threshold.
 
double getIntegralOfTimeOverThresholdProbability (const double Tmin, const double Tmax, const int NPE) const
 Get cumulative probability of time-over-threshold distribution.
 
void operator() (const JCalibration &calibration, const JPMTData< JPMTSignal > &input, JPMTData< JPMTPulse > &output) const
 Process hits.
 
virtual void merge (JPMTData< JPMTPulse > &data) const
 Merging of PMT hits.
 
const JPMTParametersgetPMTParameters () const
 Get PMT parameters.
 
bool is_valid () const
 Check validity of PMT parameters.
 
JProperties getProperties (const JEquationParameters &equation=JPMTParameters::getEquationParameters())
 Get properties of this class.
 
JProperties getProperties (const JEquationParameters &equation=JPMTParameters::getEquationParameters()) const
 Get properties of this class.
 

Static Public Member Functions

static double getMaximalRiseTime (const double th)
 Get maximal rise time for given threshold.
 
static double getTH0 ()
 Get lower threshold for rise time evaluation.
 
static double getTH1 ()
 Get upper threshold for rise time evaluation.
 
static double getTmin ()
 Get two photo-electron resolution for time-over-threshold.
 
static double getQmin ()
 Get width of charge distribution.
 
static JEquationParametersgetEquationParameters ()
 Get equation parameters.
 
static void setEquationParameters (const JEquationParameters &equation)
 Set equation parameters.
 

Public Attributes

double QE
 relative quantum efficiency
 
double gain
 gain [unit]
 
double gainSpread
 gain spread [unit]
 
double riseTime_ns
 rise time of analogue pulse [ns]
 
double TTS_ns
 transition time spread [ns]
 
double threshold
 threshold [npe]
 
double PunderAmplified
 probability of underamplified hit
 
double thresholdBand
 threshold-band [npe]
 
double mean_ns
 mean time-over-threshold of threshold-band hits [ns]
 
double sigma_ns
 time-over-threshold standard deviation of threshold-band hits [ns]
 
double slope
 slope [ns/npe]
 
double saturation
 saturation [ns]
 
bool slewing
 time slewing of analogue signal
 

Protected Attributes

double decayTime_ns
 decay time [ns]
 
double t1
 time at match point [ns]
 
double y1
 amplitude at match point [npe]
 
double x1
 Transition point from a logarithmic to a linear relation between time-over-threshold and number of photo-electrons.
 

Friends

std::istream & operator>> (std::istream &in, JPMTAnalogueSignalProcessor &object)
 Read PMT signal from input.
 

Detailed Description

PMT analogue signal processor.

This class provides for an implementation of the JDETECTOR::JPMTSignalProcessorInterface using a specific model for the analogue pulse of the PMT.
In this, the leading edge of the analogue pulse from the PMT is assumed to be a Gaussian and the tail an exponential.
The width of the Gaussian is referred to as the rise time and the inverse slope of the exponential to the decay time.
The two functions are matched at a point where the values and first derivatives are identical.
Note that the decay time is related to the rise time via the specification JDETECTOR::TIME_OVER_THRESHOLD_NS.

The charge distribution is assumed to be a Gaussian which is centered at the specified gain and truncated by the specified threshold.

The transition times are generated according the specified spread as follows.

  • If the specified transition-time spread (TTS) is negative, the transition times are generated according to measurements (see method JDETECTOR::getTransitionTime).
    In this, the negated integral value of the TTS corresponds to the option which in turn corresponds to the detector identifier of the measurements.
  • If the specified TTS is positive, the transition times are generated according a Gaussian with a sigma equals to the given TTS.
  • If the TTS is zero, the transition times are generated without any spread.

Definition at line 52 of file JPMTAnalogueSignalProcessor.hh.

Member Enumeration Documentation

◆ JThresholdDomain

Threshold domain specifiers.

Enumerator
BELOW_THRESHOLD 

below threshold

THRESHOLDBAND 

inside threshold band

ABOVE_THRESHOLD 

above threshold

Definition at line 59 of file JPMTAnalogueSignalProcessor.hh.

59 {
60 BELOW_THRESHOLD = -1, //!< below threshold
61 THRESHOLDBAND = 0, //!< inside threshold band
62 ABOVE_THRESHOLD = 2 //!< above threshold
63 };

Constructor & Destructor Documentation

◆ JPMTAnalogueSignalProcessor()

JDETECTOR::JPMTAnalogueSignalProcessor::JPMTAnalogueSignalProcessor ( const JPMTParameters & parameters = JPMTParameters())
inline

Constructor.

Parameters
parametersPMT parameters

Definition at line 71 of file JPMTAnalogueSignalProcessor.hh.

71 :
73 JPMTParameters(parameters),
74 decayTime_ns(0.0),
75 t1(0.0),
76 y1(0.0),
77 x1(std::numeric_limits<double>::max())
78 {
79 configure();
80 }
JPMTParameters()
Default constructor.
double x1
Transition point from a logarithmic to a linear relation between time-over-threshold and number of ph...
void configure()
Configure internal parameters.

Member Function Documentation

◆ configure()

void JDETECTOR::JPMTAnalogueSignalProcessor::configure ( )
inline

Configure internal parameters.

This method provides the implementations for

  • matching of the leading edge of the analogue pulse (Gaussian) and the tail (exponential); and
  • determination of number of photo-electrons above which the time-over-threshold linearly depends on the number of photo-electrons (apart from saturation).

Note that this method will throw an error if the value of the rise time (i.e. width of the Gaussian) is too large with respect to the specification JDETECTOR::TIME_OVER_THRESHOLD_NS.

Definition at line 94 of file JPMTAnalogueSignalProcessor.hh.

95 {
96 static const int N = 100;
97 static const double precision = 1.0e-4;
98
99 // check thresholdband
100
101 if (threshold - thresholdBand < getTH0()) {
102 THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::configure(): Invalid thresholdband [npe] " << thresholdBand);
103 }
104
105 // check rise time
106
108 THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::configure(): Invalid rise time [ns] " << riseTime_ns);
109 }
110
111 // decay time
112
113 const double y = -log(threshold);
114
115 const double a = y;
116 const double b = riseTime_ns * sqrt(2.0*y) - TIME_OVER_THRESHOLD_NS;
117 const double c = 0.5*riseTime_ns*riseTime_ns;
118 const double Q = b*b - 4.0*a*c;
119
120 if (Q > 0.0)
121 decayTime_ns = (-b + sqrt(Q)) / (2.0*a);
122 else
123 decayTime_ns = -b / (2.0*a);
124
125 // fix matching of Gaussian and exponential
126
127 const double x = riseTime_ns / decayTime_ns;
128
129 t1 = riseTime_ns*x;
130 y1 = exp(-0.5*x*x);
131
132 // determine transition point to linear dependence of time-over-threshold as a function of number of photo-electrons
133
134 const double xs = saturation; // disable saturation
135
136 saturation = 1.0e50;
137
138 x1 = std::numeric_limits<double>::max(); // disable linearisation
139
140 double xmin = 1.0;
141 double xmax = 1.0 / (getDerivative(1.0) * slope);
142
143 for (int i = 0; i != N; ++i) {
144
145 const double x = 0.5 * (xmin + xmax);
146 const double u = getDerivative(x) * slope;
147
148 if (fabs(1.0 - u) < precision) {
149 break;
150 }
151
152 if (u < 1.0)
153 xmin = x;
154 else
155 xmax = x;
156 }
157
158 x1 = 0.5 * (xmin + xmax);
159
160 saturation = xs; // restore saturation
161 }
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
double thresholdBand
threshold-band [npe]
double riseTime_ns
rise time of analogue pulse [ns]
double threshold
threshold [npe]
double slope
slope [ns/npe]
double saturation
saturation [ns]
const double a
const double xmax
const double xmin
const double TIME_OVER_THRESHOLD_NS
Specification for time-over-threshold corresponding to a one photo-electron pulse.
static double getMaximalRiseTime(const double th)
Get maximal rise time for given threshold.
static double getTH0()
Get lower threshold for rise time evaluation.
virtual double getDerivative(const double npe) const override
Get derivative of number of photo-electrons to time-over-threshold.

◆ getDecayTime() [1/2]

double JDETECTOR::JPMTAnalogueSignalProcessor::getDecayTime ( ) const
inline

Get decay time.

Returns
decay time [ns]

Definition at line 169 of file JPMTAnalogueSignalProcessor.hh.

170 {
171 return decayTime_ns;
172 }

◆ getT1()

double JDETECTOR::JPMTAnalogueSignalProcessor::getT1 ( ) const
inline

Get time at transition point from Gaussian to exponential.

Returns
time [ns]

Definition at line 180 of file JPMTAnalogueSignalProcessor.hh.

181 {
182 return t1;
183 }

◆ getY1()

double JDETECTOR::JPMTAnalogueSignalProcessor::getY1 ( ) const
inline

Get amplitude at transition point from Gaussian to exponential.

Returns
amplitude [npe]

Definition at line 191 of file JPMTAnalogueSignalProcessor.hh.

192 {
193 return y1;
194 }

◆ getStartOfLinearisation()

double JDETECTOR::JPMTAnalogueSignalProcessor::getStartOfLinearisation ( ) const
inline

Get transition point from a model-dependent to linear relation between time-over-threshold and number of photo-electrons.

Returns
number of photo-electrons [npe]

Definition at line 202 of file JPMTAnalogueSignalProcessor.hh.

203 {
204 return x1;
205 }

◆ getAmplitude()

double JDETECTOR::JPMTAnalogueSignalProcessor::getAmplitude ( const double t1_ns) const
inline

Get amplitude at given time for a one photo-electron pulse.

Parameters
t1_nstime [ns]
Returns
amplitude [npe]

Definition at line 214 of file JPMTAnalogueSignalProcessor.hh.

215 {
216 if (t1_ns < t1) {
217
218 const double x = t1_ns / riseTime_ns;
219
220 return exp(-0.5*x*x); // Gaussian
221
222 } else {
223
224 const double x = t1_ns / decayTime_ns;
225
226 return exp(-x) / y1; // exponential
227 }
228 }

◆ getRiseTime() [1/2]

double JDETECTOR::JPMTAnalogueSignalProcessor::getRiseTime ( const double npe,
const double th ) const
inline

Get time to pass from threshold to top of analogue pulse.


In this, the leading edge of the analogue pulse is assumed to be Gaussian.

Parameters
npenumber of photo-electrons
ththreshold [npe]
Returns
time [ns]

Definition at line 239 of file JPMTAnalogueSignalProcessor.hh.

240 {
241 return riseTime_ns * sqrt(2.0*log(npe/th)); // Gaussian
242 }

◆ getDecayTime() [2/2]

double JDETECTOR::JPMTAnalogueSignalProcessor::getDecayTime ( const double npe,
const double th ) const
inline

Get time to pass from top of analogue pulse to threshold.


In this, the trailing edge of the analogue pulse is assumed to be exponential.

Parameters
npenumber of photo-electrons
ththreshold [npe]
Returns
time [ns]

Definition at line 253 of file JPMTAnalogueSignalProcessor.hh.

254 {
255 if (npe*y1 > th)
256 return decayTime_ns * (log(npe/th) - log(y1)); // exponential
257 else
258 return riseTime_ns * sqrt(2.0*log(npe/th)); // Gaussian
259 }

◆ getMaximalRiseTime()

static double JDETECTOR::JPMTAnalogueSignalProcessor::getMaximalRiseTime ( const double th)
inlinestatic

Get maximal rise time for given threshold.

Note that the rise time is entirely constrained by the specification JDETECTOR::TIME_OVER_THRESHOLD_NS.

Parameters
ththreshold [npe]
Returns
rise time [ns]

Definition at line 270 of file JPMTAnalogueSignalProcessor.hh.

271 {
272 if (th > 0.0 && th < 1.0)
273 return 0.5 * TIME_OVER_THRESHOLD_NS / sqrt(-2.0*log(th));
274 else
275 THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getMaximalRiseTime(): Invalid threshold " << th);
276 }

◆ applySaturation()

double JDETECTOR::JPMTAnalogueSignalProcessor::applySaturation ( const double tot_ns) const
inline

Get time-over-threshold with saturation.

Parameters
tot_nstime-over-threshold without saturation
Returns
time-over-threshold with saturation

Definition at line 285 of file JPMTAnalogueSignalProcessor.hh.

286 {
287 return saturation / sqrt(tot_ns*tot_ns + saturation*saturation) * tot_ns;
288 }

◆ removeSaturation()

double JDETECTOR::JPMTAnalogueSignalProcessor::removeSaturation ( const double tot_ns) const
inline

Get time-over-threshold without saturation.

Parameters
tot_nstime-over-threshold with saturation
Returns
time-over-threshold without saturation

Definition at line 297 of file JPMTAnalogueSignalProcessor.hh.

298 {
299 if (tot_ns < saturation)
300 return saturation / sqrt(saturation*saturation - tot_ns*tot_ns) * tot_ns;
301 else
302 return std::numeric_limits<double>::max();
303 //THROW(JValueOutOfRange, "Time-over-threshold exceeds saturation " << tot_ns << " >= " << saturation);
304 }

◆ getDerivativeOfSaturation()

double JDETECTOR::JPMTAnalogueSignalProcessor::getDerivativeOfSaturation ( const double tot_ns) const
inline

Get derivative of saturation factor.

Parameters
tot_nstime-over-threshold without saturation
Returns
derivative of saturation factor

Definition at line 313 of file JPMTAnalogueSignalProcessor.hh.

314 {
315 return saturation * saturation * saturation / ((saturation*saturation + tot_ns*tot_ns) * sqrt(saturation*saturation + tot_ns*tot_ns));
316 }

◆ getGainSpread()

double JDETECTOR::JPMTAnalogueSignalProcessor::getGainSpread ( int NPE) const
inline

Get gain spread for given number of photo-electrons.

Parameters
NPEnumber of photo-electrons
Returns
gain spread

Definition at line 325 of file JPMTAnalogueSignalProcessor.hh.

326 {
327 return sqrt((double) NPE * gain) * gainSpread;
328 }
double gainSpread
gain spread [unit]

◆ getIntegralOfChargeProbability() [1/2]

double JDETECTOR::JPMTAnalogueSignalProcessor::getIntegralOfChargeProbability ( const double xmin,
const double xmax,
const int NPE ) const
inline

Get integral of probability.

Parameters
xminminimum number of photo-electrons
xmaxmaximum number of photo-electrons
NPEtrue number of photo-electrons
Returns
probability

Definition at line 339 of file JPMTAnalogueSignalProcessor.hh.

340 {
341 double zmin = xmin;
342 double zmax = xmax;
343
344 const double th = threshold - thresholdBand;
345 const double f = gainSpread * gainSpread;
346
347 if (zmin < th) { zmin = th; }
348 if (zmax < th) { zmax = th; }
349
350 double norm = 0.0;
351 double cumulP = 0.0;
352 double weight = (PunderAmplified > 0.0 ? pow(PunderAmplified, NPE) : 1.0);
353
354 for (int k = (PunderAmplified > 0.0 ? 0 : NPE); k <= NPE; ++k) {
355
356 const double mu = ((NPE-k) * f * gain) + (k * gain);
357 const double sigma = sqrt((NPE-k) * f + k) * getGainSpread(1);
358
359 norm += weight * (0.5 * TMath::Erfc((th - mu) / sqrt(2.0) / sigma));
360 cumulP += weight * (0.5 * TMath::Erfc((zmin - mu) / sqrt(2.0) / sigma) -
361 0.5 * TMath::Erfc((zmax - mu) / sqrt(2.0) / sigma));
362
363 weight *= ( (NPE - k) / ((double) (k+1)) *
365 }
366
367 return cumulP / norm;
368 }
double PunderAmplified
probability of underamplified hit
const double sigma[]
T pow(const T &x, const double y)
Power .
Definition JMath.hh:97
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.

◆ getIntegralOfChargeProbability() [2/2]

double JDETECTOR::JPMTAnalogueSignalProcessor::getIntegralOfChargeProbability ( const JThresholdDomain domain,
const int NPE ) const
inline

Get integral of probability in specific threshold domain.

Parameters
domainthreshold domain
NPEtrue number of photo-electrons
Returns
probability

Definition at line 378 of file JPMTAnalogueSignalProcessor.hh.

379 {
380 switch (domain) {
381
382 case ABOVE_THRESHOLD:
384
385 case THRESHOLDBAND:
387
388 default:
389 return 0.0;
390 }
391 }
double getIntegralOfChargeProbability(const double xmin, const double xmax, const int NPE) const
Get integral of probability.

◆ setPMTParameters()

void JDETECTOR::JPMTAnalogueSignalProcessor::setPMTParameters ( const JPMTParameters & parameters)
inline

Set PMT parameters.

Parameters
parametersPMT parameters

Definition at line 399 of file JPMTAnalogueSignalProcessor.hh.

400 {
401 static_cast<JPMTParameters&>(*this).setPMTParameters(parameters);
402
403 configure();
404 }

◆ getThresholdDomain()

JThresholdDomain JDETECTOR::JPMTAnalogueSignalProcessor::getThresholdDomain ( const double npe) const
inline

Get threshold domain.

Parameters
npenumber of photo-electrons
Returns
threshold domain

Definition at line 430 of file JPMTAnalogueSignalProcessor.hh.

431 {
432 if (npe > threshold) {
433
434 return ABOVE_THRESHOLD;
435
436 } else if (npe > threshold - thresholdBand) {
437
438 return THRESHOLDBAND;
439
440 } else {
441
442 return BELOW_THRESHOLD;
443 }
444 }

◆ applyQE()

virtual bool JDETECTOR::JPMTAnalogueSignalProcessor::applyQE ( ) const
inlineoverridevirtual

Apply relative QE.

Returns
true if accepted; false if rejected

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 452 of file JPMTAnalogueSignalProcessor.hh.

453 {
454 if (QE <= 0.0)
455 return false;
456 else if (QE < 1.0)
457 return gRandom->Rndm() < QE;
458 else
459 return true;
460 }
double QE
relative quantum efficiency

◆ getRandomTime()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getRandomTime ( const double t_ns) const
inlineoverridevirtual

Get randomised time according transition time distribution.

Parameters
t_nstime [ns]
Returns
time [ns]

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 469 of file JPMTAnalogueSignalProcessor.hh.

470 {
471 if (TTS_ns < 0.0)
472 return t_ns + getTransitionTime(gRandom->Rndm(), -lrint(TTS_ns));
473 else if (TTS_ns > 0.0)
474 return gRandom->Gaus(t_ns, TTS_ns);
475 else
476 return t_ns;
477 }
double TTS_ns
transition time spread [ns]
JDETECTOR::getTransitionTime getTransitionTime
Function object to generate transition time.

◆ compare()

virtual bool JDETECTOR::JPMTAnalogueSignalProcessor::compare ( const JPhotoElectron & first,
const JPhotoElectron & second ) const
inlineoverridevirtual

Compare arrival times of photo-electrons.

This implementation uses the internal rise time as two photo-electron resolution.

Two (or more) photo-electrons are merged if they are comparable.

Parameters
firstfirst photo-electron
secondsecond photo-electron
Returns
true if arrival times of photo-electrons are within two photo-electron resolution; else false

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 490 of file JPMTAnalogueSignalProcessor.hh.

491 {
492 return second.t_ns < first.t_ns + riseTime_ns;
493 }

◆ getRandomCharge()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getRandomCharge ( const int NPE) const
inlineoverridevirtual

Get randomised charge according to gain and gain spread.

Parameters
NPEnumber of photo-electrons
Returns
number of photo-electrons

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 502 of file JPMTAnalogueSignalProcessor.hh.

503 {
504 double q;
505
506 do {
507
508 // Determine which contribution to sample from
509 int k = NPE;
510
511 if (PunderAmplified > 0.0) {
512
513 const double X = gRandom->Uniform();
514 double weight = pow(1.0 - PunderAmplified, NPE);
515
516 for (double sum_p = 0.0; k > 0; --k) {
517
518 sum_p += weight;
519 if (sum_p > X) { break; }
520
521 weight *= ( k / ((double) (NPE - (k-1))) *
523 }
524 }
525
526 // Sample from chosen contribution
527 const double f = gainSpread * gainSpread;
528 const double mu = ((NPE-k) * f * gain) + (k * gain);
529 const double sigma = sqrt((NPE-k) * f + k) * getGainSpread(1);
530
531 q = gRandom->Gaus(mu,sigma);
532
533 } while (q < 0.0);
534
535 return q;
536 }

◆ getChargeProbability()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getChargeProbability ( const double npe,
const int NPE ) const
inlineoverridevirtual

Get probability density for given charge.

Parameters
npeobserved number of photo-electrons
NPEtrue number of photo-electrons
Returns
probability [npe^-1]

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 546 of file JPMTAnalogueSignalProcessor.hh.

547 {
549
550 const double f = gainSpread * gainSpread;
551
552 double norm = 0.0;
553 double prob = 0.0;
554 double weight = (PunderAmplified > 0.0 ? pow(PunderAmplified, NPE) : 1.0);
555
556 for (int k = (PunderAmplified > 0.0 ? 0 : NPE); k <= NPE; ++k) {
557
558 const double mu = ((NPE-k) * f * gain) + (k * gain);
559 const double sigma = sqrt((NPE-k) * f + k) * getGainSpread(1);
560
561 norm += weight * (0.5 * TMath::Erfc(((threshold-thresholdBand) - mu) / sqrt(2.0) / sigma));
562 prob += weight * TMath::Gaus(npe, mu, sigma, kTRUE);
563
564 weight *= ( (NPE - k) / ((double) (k+1)) *
566 }
567
568 return prob / norm;
569 }
570
571 return 0.0;
572 }
JThresholdDomain getThresholdDomain(const double npe) const
Get threshold domain.

◆ applyThreshold()

virtual bool JDETECTOR::JPMTAnalogueSignalProcessor::applyThreshold ( const double npe) const
inlineoverridevirtual

Apply threshold.

Parameters
npenumber of photo-electrons
Returns
true if pass; else false

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 581 of file JPMTAnalogueSignalProcessor.hh.

582 {
584 }

◆ getRiseTime() [2/2]

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getRiseTime ( const double npe) const
inlineoverridevirtual

Get time to reach threshold.

Note that the rise time is defined to be zero for a one photo-electron signal.

Parameters
npenumber of photo-electrons
Returns
time [ns]

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 595 of file JPMTAnalogueSignalProcessor.hh.

596 {
597 if (slewing) {
598
599 switch (getThresholdDomain(npe)) {
600
601 case THRESHOLDBAND:
602 return ((getRiseTime(npe, getTH0()) - getRiseTime(npe, threshold-thresholdBand)) -
603 (getRiseTime(1.0, getTH0()) - getRiseTime(1.0, threshold-thresholdBand))) + this->mean_ns;
604
605 case ABOVE_THRESHOLD:
606 return ((getRiseTime(npe, getTH0()) - getRiseTime(npe, threshold)) -
607 (getRiseTime(1.0, getTH0()) - getRiseTime(1.0, threshold)));
608
609 default:
610 THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getRiseTime: Invalid charge " << npe);
611 }
612
613 } else {
614
615 return 0.0;
616 }
617 }
double mean_ns
mean time-over-threshold of threshold-band hits [ns]
bool slewing
time slewing of analogue signal
double getRiseTime(const double npe, const double th) const
Get time to pass from threshold to top of analogue pulse.

◆ getTimeOverThreshold()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getTimeOverThreshold ( const double npe) const
inlineoverridevirtual

Get time-over-threshold (ToT).

Parameters
npenumber of photo-electrons
Returns
ToT [ns]

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 626 of file JPMTAnalogueSignalProcessor.hh.

627 {
628 switch (getThresholdDomain(npe)) {
629
630 case THRESHOLDBAND: {
631
632 return gRandom->Gaus(mean_ns, sigma_ns);
633 }
634
635 case ABOVE_THRESHOLD: {
636
637 double tot = 0.0;
638
639 if (npe*y1 <= threshold) {
640
641 tot += getRiseTime(npe, threshold); // Gaussian
642 tot += getRiseTime(npe, threshold); // Gaussian
643
644 } else if (npe <= getStartOfLinearisation()) {
645
646 tot += getRiseTime (npe, threshold); // Gaussian
647 tot += getDecayTime(npe, threshold); // exponential
648
649 } else {
650
651 tot += getRiseTime (getStartOfLinearisation(), threshold); // Gaussian
652 tot += getDecayTime(getStartOfLinearisation(), threshold); // exponential
653
654 tot += slope * (npe - getStartOfLinearisation()); // linear
655 }
656
657 return applySaturation(tot);
658 }
659
660 default: {
661
662 THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getTimeOverThreshold: Invalid charge " << npe);
663 }}
664 }
double sigma_ns
time-over-threshold standard deviation of threshold-band hits [ns]
double applySaturation(const double tot_ns) const
Get time-over-threshold with saturation.
double getStartOfLinearisation() const
Get transition point from a model-dependent to linear relation between time-over-threshold and number...

◆ getDerivative()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getDerivative ( const double npe) const
inlineoverridevirtual

Get derivative of number of photo-electrons to time-over-threshold.

Parameters
npenumber of photo-electrons
Returns
dnpe/dToT [ns^-1]

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 673 of file JPMTAnalogueSignalProcessor.hh.

674 {
675 switch (getThresholdDomain(npe)) {
676
677 case ABOVE_THRESHOLD: {
678
679 const double z = riseTime_ns / sqrt(2.0 * log(npe/threshold));
680
681 double y = 0.0;
682
683 if (npe*y1 > threshold) {
684
685 if (npe <= getStartOfLinearisation())
686 y = npe / (z + decayTime_ns); // Gaussian + exponential
687 else
688 y = 1.0 / slope; // linear
689
690 } else {
691
692 y = npe / (2.0 * z); // Gaussian + Gaussian
693 }
694
695 const double tot_ns = getTimeOverThreshold(npe);
696
698 }
699
700 default:
701 return 0.0;
702 }
703 }
double getDerivativeOfSaturation(const double tot_ns) const
Get derivative of saturation factor.
double removeSaturation(const double tot_ns) const
Get time-over-threshold without saturation.
virtual double getTimeOverThreshold(const double npe) const override
Get time-over-threshold (ToT).

◆ getSurvivalProbability()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getSurvivalProbability ( const int NPE) const
inlineoverridevirtual

Probability that a hit survives the simulation of the PMT.

Parameters
NPEnumber of photo-electrons
Returns
probability

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 712 of file JPMTAnalogueSignalProcessor.hh.

713 {
714 if (QE <= 0.0) {
715
716 return 0.0;
717
718 } else if (QE <= 1.0) {
719
720 double P = 0.0;
721
722 const double f = gainSpread * gainSpread;
723
724 for (int i = 1; i <= NPE; ++i) {
725
726 const double p = (TMath::Binomial(NPE, i) *
727 TMath::Power(QE, i) * TMath::Power(1.0 - QE, NPE - i));
728
729 double Ptotal = 0.0;
730 double Pabove = 0.0;
731 double weight = (PunderAmplified > 0.0 ? pow(PunderAmplified, i) : 1.0);
732
733 for (int k = (PunderAmplified > 0.0 ? 0 : NPE); k <= NPE; ++k) {
734
735 const double mu = ((i-k) * f * gain) + (k * gain);
736 const double sigma = sqrt((i-k) * f + k) * getGainSpread(1);
737
738 Ptotal += weight * (0.5 * TMath::Erfc((0.0 - mu) / sqrt(2.0) / sigma));
739 Pabove += weight * (0.5 * TMath::Erfc((threshold - thresholdBand - mu) / sqrt(2.0) / sigma));
740
741 weight *= ( (i - k) / ((double) (k+1)) *
743 }
744
745 P += p*Pabove/Ptotal;
746 }
747
748 return P;
749
750 } else {
751
752 THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getSurvivalProbability: Invalid QE " << QE);
753 }
754 }

◆ getNPE()

virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getNPE ( const double tot_ns) const
inlineoverridevirtual

Get number of photo-electrons.

Parameters
tot_nstime over threshold [ns]
Returns
number of photo-electrons

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 763 of file JPMTAnalogueSignalProcessor.hh.

764 {
765 if (tot_ns >= saturation) {
766 return std::numeric_limits<double>::max();
767 }
768
769 const double tot = removeSaturation(tot_ns);
770 const double TOT = (getRiseTime (getStartOfLinearisation(), threshold) +
772
773 if (tot <= 2*getRiseTime(threshold/y1,threshold)) { // Gaussian + Gaussian
774
775 return threshold * exp(tot*tot/riseTime_ns/riseTime_ns/8.0);
776
777 } else if (tot <= TOT) { // Gaussian + Exponential
778
779 const double a = decayTime_ns;
780 const double b = sqrt(2.0) * riseTime_ns;
781 const double c = -(decayTime_ns*log(y1) + tot);
782 const double z = (-b + sqrt(b*b - 4*a*c)) / (2*a);
783
784 return threshold * exp(z*z);
785
786 } else { // linear
787
788 return getStartOfLinearisation() + (tot - TOT) / slope;
789 }
790 }

◆ getTimeOverThresholdProbability()

double JDETECTOR::JPMTAnalogueSignalProcessor::getTimeOverThresholdProbability ( const double tot_ns,
const int NPE ) const
inline

Get probability of having a pulse with specific time-over-threshold.

Parameters
tot_nstime-over-threshold with saturation [ns]
NPEtrue number of photo-electrons
Returns
probability [ns^-1]

Definition at line 800 of file JPMTAnalogueSignalProcessor.hh.

801 {
802 const double PthBand = getIntegralOfChargeProbability(THRESHOLDBAND, NPE);
803
804 const double npe = getNPE(tot_ns);
805 const double y = getChargeProbability(npe, NPE);
806 const double v = getDerivative(npe);
807
808 const double RthBand = PthBand * TMath::Gaus(tot_ns, mean_ns, sigma_ns, kTRUE);
809 const double RaboveTh = y * v;
810
811 return RthBand + RaboveTh;
812 }
virtual double getChargeProbability(const double npe, const int NPE) const override
Get probability density for given charge.
virtual double getNPE(const double tot_ns) const override
Get number of photo-electrons.

◆ getIntegralOfTimeOverThresholdProbability()

double JDETECTOR::JPMTAnalogueSignalProcessor::getIntegralOfTimeOverThresholdProbability ( const double Tmin,
const double Tmax,
const int NPE ) const
inline

Get cumulative probability of time-over-threshold distribution.

Parameters
Tminminimum time-over-threshold (with saturation) [ns]
Tmaxmaximum time-over-threshold (with saturation) [ns]
NPEtrue number of photo-electrons
Returns
probability [ns^-1]

Definition at line 823 of file JPMTAnalogueSignalProcessor.hh.

824 {
825
826 const double PthBand = getIntegralOfChargeProbability(THRESHOLDBAND, NPE);
827
828 const double IthBand = PthBand * (0.5 * TMath::Erfc((Tmin - mean_ns) / sqrt(2.0) / sigma_ns) -
829 0.5 * TMath::Erfc((Tmax - mean_ns) / sqrt(2.0) / sigma_ns));
830 const double IaboveTh = getIntegralOfChargeProbability(getNPE(Tmin), getNPE(Tmax), NPE);
831
832 return IthBand + IaboveTh;
833 }

◆ getTH0()

static double JDETECTOR::JPMTAnalogueSignalProcessor::getTH0 ( )
inlinestatic

Get lower threshold for rise time evaluation.

Returns
threshold [npe]

Definition at line 841 of file JPMTAnalogueSignalProcessor.hh.

842 {
843 return 0.1;
844 }

◆ getTH1()

static double JDETECTOR::JPMTAnalogueSignalProcessor::getTH1 ( )
inlinestatic

Get upper threshold for rise time evaluation.

Returns
threshold [npe]

Definition at line 852 of file JPMTAnalogueSignalProcessor.hh.

853 {
854 return 0.9;
855 }

◆ operator()()

void JDETECTOR::JPMTSignalProcessorInterface::operator() ( const JCalibration & calibration,
const JPMTData< JPMTSignal > & input,
JPMTData< JPMTPulse > & output ) const
inlineinherited

Process hits.

Two (or more) photo-electrons are combined if they are comparable according method compare.
Two (or more) consecutive hits hits maybe merged (according method merge).

Parameters
calibrationPMT calibration
inputPMT signals
outputPMT hits

Definition at line 89 of file JPMTSignalProcessorInterface.hh.

92 {
93 // apply transition time distribution to each photo-electron.
94
95 JPMTData<JPhotoElectron> buffer;
96
97 for (JPMTData<JPMTSignal>::const_iterator hit = input.begin(); hit != input.end(); ++hit) {
98
99 for (int i = 0; i != hit->npe; ++i) {
100
101 if (applyQE()) {
102 buffer.push_back(JPhotoElectron(getRandomTime(hit->t_ns)));
103 }
104 }
105 }
106
107 if (!buffer.empty()) {
108
109 buffer.push_back(JPhotoElectron::getEndMarker());
110
111 buffer.sort();
112
113
114 // generate PMT hits from time sequence of photo-electrons.
115
116 for (JPMTData<JPhotoElectron>::const_iterator q = buffer.begin(), p = q++; q != buffer.end(); ++q) {
117
118 while (compare(*p,*q)) {
119 ++q;
120 }
121
122 const double npe = getRandomCharge(distance(p,q));
123
124 if (applyThreshold(npe)) {
125 output.push_back(JPMTPulse(putTime(p->t_ns + getRiseTime(npe), calibration), getTimeOverThreshold(npe)));
126 }
127
128 p = q;
129 }
130
131 // merge overlapping PMT hits.
132
133 merge(output);
134 }
135 }
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
std::vector< JElement_t >::const_iterator const_iterator
virtual bool applyQE() const
Apply relative QE.
virtual double getTimeOverThreshold(const double npe) const
Get time-over-threshold (ToT).
virtual double getRandomCharge(const int NPE) const
Get randomised charge according to gain and gain spread.
virtual bool applyThreshold(const double npe) const
Apply threshold.
virtual bool compare(const JPhotoElectron &first, const JPhotoElectron &second) const
Compare arrival times of photo-electrons.
virtual double getRiseTime(const double npe) const
Get time to reach threshold.
virtual double getRandomTime(const double t_ns) const
Get randomised time according transition time distribution.
virtual void merge(JPMTData< JPMTPulse > &data) const
Merging of PMT hits.
double putTime(const T &t1, const JCalibration &cal)
Get de-calibrated time.
static JPhotoElectron getEndMarker()
Get end marker.

◆ merge()

virtual void JDETECTOR::JPMTSignalProcessorInterface::merge ( JPMTData< JPMTPulse > & data) const
inlinevirtualinherited

Merging of PMT hits.

Hits with overlapping time-over-threshold signals should -de facto- be combined. In this, the leading edge is maintained and the time-over-threshold is set to the difference between the overall trailing and leading edges. As a result, the number of PMT hits may be reduced.

Parameters
dataPMT hits (I/O)

Definition at line 293 of file JPMTSignalProcessorInterface.hh.

294 {
295 using namespace std;
296
298
299 for (JPMTData<JPMTPulse>::iterator i = data.begin(); i != data.end(); ) {
300
301 double t1 = i->t_ns;
302 double t2 = i->t_ns + i->tot_ns;
303
304 while (++i != data.end() && i->t_ns < t2 + getTmin()) {
305 t2 = max(t2, i->t_ns + i->tot_ns);
306 }
307
308 out->t_ns = t1;
309 out->tot_ns = t2 - t1;
310
311 ++out;
312 }
313
314 data.resize(distance(data.begin(), out));
315 }
std::vector< JElement_t >::iterator iterator
static double getTmin()
Get two photo-electron resolution for time-over-threshold.

◆ getTmin()

static double JDETECTOR::JPMTSignalProcessorInterface::getTmin ( )
inlinestaticinherited

Get two photo-electron resolution for time-over-threshold.

Returns
minimal time [ns]

Definition at line 323 of file JPMTSignalProcessorInterface.hh.

324 {
325 return 1.0;
326 }

◆ getQmin()

static double JDETECTOR::JPMTSignalProcessorInterface::getQmin ( )
inlinestaticinherited

Get width of charge distribution.

Returns
width charge distribution [npe]

Definition at line 334 of file JPMTSignalProcessorInterface.hh.

335 {
336 return 1.0e-3;
337 }

◆ getPMTParameters()

const JPMTParameters & JDETECTOR::JPMTParameters::getPMTParameters ( ) const
inlineinherited

Get PMT parameters.

Returns
PMT parameters

Definition at line 111 of file JPMTParameters.hh.

112 {
113 return static_cast<const JPMTParameters&>(*this);
114 }

◆ is_valid()

bool JDETECTOR::JPMTParameters::is_valid ( ) const
inlineinherited

Check validity of PMT parameters.

Returns
true if valid; else false

Definition at line 133 of file JPMTParameters.hh.

134 {
135 if (this->QE < 0.0 ||
136 this->gain < 0.0 ||
137 this->gainSpread < 0.0 ||
138 this->threshold < 0.0 ||
139 this->thresholdBand < 0.0) {
140 return false;
141 }
142
143 return true;
144 }

◆ getEquationParameters()

static JEquationParameters & JDETECTOR::JPMTParameters::getEquationParameters ( )
inlinestaticinherited

Get equation parameters.

Returns
equation parameters

Definition at line 194 of file JPMTParameters.hh.

195 {
196 static JEquationParameters equation("=", ",", "./", "#");
197
198 return equation;
199 }
Simple data structure to support I/O of equations (see class JLANG::JEquation).

◆ setEquationParameters()

static void JDETECTOR::JPMTParameters::setEquationParameters ( const JEquationParameters & equation)
inlinestaticinherited

Set equation parameters.

Parameters
equationequation parameters

Definition at line 207 of file JPMTParameters.hh.

208 {
209 getEquationParameters() = equation;
210 }
static JEquationParameters & getEquationParameters()
Get equation parameters.

◆ getProperties() [1/2]

JProperties JDETECTOR::JPMTParameters::getProperties ( const JEquationParameters & equation = JPMTParameters::getEquationParameters())
inlineinherited

Get properties of this class.

Parameters
equationequation parameters

Definition at line 218 of file JPMTParameters.hh.

219 {
220 return JPMTParametersHelper(*this, equation);
221 }

◆ getProperties() [2/2]

JProperties JDETECTOR::JPMTParameters::getProperties ( const JEquationParameters & equation = JPMTParameters::getEquationParameters()) const
inlineinherited

Get properties of this class.

Parameters
equationequation parameters

Definition at line 229 of file JPMTParameters.hh.

230 {
231 return JPMTParametersHelper(*this, equation);
232 }

Friends And Related Symbol Documentation

◆ operator>>

std::istream & operator>> ( std::istream & in,
JPMTAnalogueSignalProcessor & object )
friend

Read PMT signal from input.

Parameters
ininput stream
objectPMT signal
Returns
input stream

Definition at line 414 of file JPMTAnalogueSignalProcessor.hh.

415 {
416 in >> static_cast<JPMTParameters&>(object);
417
418 object.configure();
419
420 return in;
421 }

Member Data Documentation

◆ decayTime_ns

double JDETECTOR::JPMTAnalogueSignalProcessor::decayTime_ns
protected

decay time [ns]

Definition at line 859 of file JPMTAnalogueSignalProcessor.hh.

◆ t1

double JDETECTOR::JPMTAnalogueSignalProcessor::t1
protected

time at match point [ns]

Definition at line 860 of file JPMTAnalogueSignalProcessor.hh.

◆ y1

double JDETECTOR::JPMTAnalogueSignalProcessor::y1
protected

amplitude at match point [npe]

Definition at line 861 of file JPMTAnalogueSignalProcessor.hh.

◆ x1

double JDETECTOR::JPMTAnalogueSignalProcessor::x1
protected

Transition point from a logarithmic to a linear relation between time-over-threshold and number of photo-electrons.


Measurements by B. Schermer and R. Bruijn at Nikhef.

Definition at line 867 of file JPMTAnalogueSignalProcessor.hh.

◆ QE

double JDETECTOR::JPMTParameters::QE
inherited

relative quantum efficiency

Definition at line 235 of file JPMTParameters.hh.

◆ gain

double JDETECTOR::JPMTParameters::gain
inherited

gain [unit]

Definition at line 236 of file JPMTParameters.hh.

◆ gainSpread

double JDETECTOR::JPMTParameters::gainSpread
inherited

gain spread [unit]

Definition at line 237 of file JPMTParameters.hh.

◆ riseTime_ns

double JDETECTOR::JPMTParameters::riseTime_ns
inherited

rise time of analogue pulse [ns]

Definition at line 238 of file JPMTParameters.hh.

◆ TTS_ns

double JDETECTOR::JPMTParameters::TTS_ns
inherited

transition time spread [ns]

Definition at line 239 of file JPMTParameters.hh.

◆ threshold

double JDETECTOR::JPMTParameters::threshold
inherited

threshold [npe]

Definition at line 240 of file JPMTParameters.hh.

◆ PunderAmplified

double JDETECTOR::JPMTParameters::PunderAmplified
inherited

probability of underamplified hit

Definition at line 241 of file JPMTParameters.hh.

◆ thresholdBand

double JDETECTOR::JPMTParameters::thresholdBand
inherited

threshold-band [npe]

Definition at line 242 of file JPMTParameters.hh.

◆ mean_ns

double JDETECTOR::JPMTParameters::mean_ns
inherited

mean time-over-threshold of threshold-band hits [ns]

Definition at line 243 of file JPMTParameters.hh.

◆ sigma_ns

double JDETECTOR::JPMTParameters::sigma_ns
inherited

time-over-threshold standard deviation of threshold-band hits [ns]

Definition at line 244 of file JPMTParameters.hh.

◆ slope

double JDETECTOR::JPMTParameters::slope
inherited

slope [ns/npe]

Definition at line 245 of file JPMTParameters.hh.

◆ saturation

double JDETECTOR::JPMTParameters::saturation
inherited

saturation [ns]

Definition at line 246 of file JPMTParameters.hh.

◆ slewing

bool JDETECTOR::JPMTParameters::slewing
inherited

time slewing of analogue signal

Definition at line 247 of file JPMTParameters.hh.


The documentation for this struct was generated from the following file: