Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Protected Attributes | Friends | List of all members
JDETECTOR::JPMTAnalogueSignalProcessor Struct Reference

PMT analogue signal processor. More...

#include <JPMTAnalogueSignalProcessor.hh>

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

Public Types

enum  JThresholdDomains { BELOW_THRESHOLDBAND = -1, IN_THRESHOLDBAND = 0, ABOVE_THRESHOLD = 2 }
 Threshold domain specifiers. More...
 

Public Member Functions

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

Static Public Member Functions

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

Public Attributes

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

Protected Attributes

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

Friends

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

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.

Definition at line 52 of file JPMTAnalogueSignalProcessor.hh.

Member Enumeration Documentation

Threshold domain specifiers.

Enumerator
BELOW_THRESHOLDBAND 

below threshold band

IN_THRESHOLDBAND 

inside threshold band

ABOVE_THRESHOLD 

above threshold band

Definition at line 89 of file JPMTSignalProcessorInterface.hh.

89  {
90  BELOW_THRESHOLDBAND = -1, //!< below threshold band
91  IN_THRESHOLDBAND = 0, //!< inside threshold band
92  ABOVE_THRESHOLD = 2 //!< above threshold band
93  };

Constructor & Destructor Documentation

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

Constructor.

Parameters
parametersPMT parameters

Definition at line 61 of file JPMTAnalogueSignalProcessor.hh.

61  :
64  decayTime_ns(0.0),
65  t1(0.0),
66  y1(0.0),
67  x1(std::numeric_limits<double>::max())
68  {
69  configure();
70  }
double x1
Transition point from a logarithmic to a linear relation between time-over-threshold and number of ph...
JPMTSignalProcessorInterface(const double max_charge=100.0)
Default constructor.
*fatal Wrong number of arguments esac JCookie sh typeset Z DETECTOR typeset Z SOURCE_RUN typeset Z TARGET_RUN set_variable PARAMETERS_FILE $WORKDIR parameters
Definition: diff-Tuna.sh:38
void configure()
Configure internal parameters.
JPMTParameters()
Default constructor.
double y1
amplitude at match point [npe]

Member Function Documentation

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 84 of file JPMTAnalogueSignalProcessor.hh.

85  {
86  static const int N = 100;
87  static const double precision = 1.0e-4;
88 
89  // check thresholdband
90 
91  if (threshold - thresholdBand < getTH0()) {
92  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::configure(): Invalid thresholdband [npe] " << thresholdBand);
93  }
94 
95  // check rise time
96 
98  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::configure(): Invalid rise time [ns] " << riseTime_ns);
99  }
100 
101  // decay time
102 
103  const double y = -log(threshold);
104 
105  const double a = y;
106  const double b = riseTime_ns * sqrt(2.0*y) - TIME_OVER_THRESHOLD_NS;
107  const double c = 0.5*riseTime_ns*riseTime_ns;
108 
109  const double Q = b*b - 4.0*a*c;
110 
111  if (Q > 0.0)
112  decayTime_ns = (-b + sqrt(Q)) / (2.0*a);
113  else
114  decayTime_ns = -b / (2.0*a);
115 
116  // fix matching of Gaussian and exponential
117 
118  const double x = riseTime_ns / decayTime_ns;
119 
120  t1 = riseTime_ns*x;
121  y1 = exp(-0.5*x*x);
122 
123  // determine transition point to linear dependence of time-over-threshold as a function of number of photo-electrons
124 
125  const double xs = saturation; // disable saturation
126 
127  saturation = 1.0e50;
128 
129  x1 = std::numeric_limits<double>::max(); // disable linearisation
130 
131  double xmin = 1.0;
132  double xmax = 1.0 / (getDerivative(1.0) * slope);
133 
134  for (int i = 0; i != N; ++i) {
135 
136  const double x = 0.5 * (xmin + xmax);
137  const double u = getDerivative(x) * slope;
138 
139  if (fabs(1.0 - u) < precision) {
140  break;
141  }
142 
143  if (u < 1.0)
144  xmin = x;
145  else
146  xmax = x;
147  }
148 
149  x1 = 0.5 * (xmin + xmax);
150 
151  saturation = xs;
152  }
double x1
Transition point from a logarithmic to a linear relation between time-over-threshold and number of ph...
const double TIME_OVER_THRESHOLD_NS
Specification for time-over-threshold corresponding to a one photo-electron pulse.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:670
double threshold
threshold [npe]
static double getTH0()
Get lower threshold for rise time evaluation.
then set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))*exp(-0.5 *(y-[1])*(y-[1])/([2]*[2]))" JF2 -o $WORKDIR/f2.root -F "$FORMULA" -@ "p0
double thresholdBand
threshold-band [npe]
double riseTime_ns
rise time of analogue pulse [ns]
virtual double getDerivative(const double npe) const override
Get derivative of number of photo-electrons to time-over-threshold.
static double getMaximalRiseTime(const double th)
Get maximal rise time for given threshold.
then $JPP_DIR software JCalibrate JCalibrateToT a
Definition: JTuneHV.sh:108
double saturation
saturation [ns]
double u[N+1]
Definition: JPolint.hh:739
then usage $script[input file[working directory[option]]] nWhere option can be N
Definition: JMuonPostfit.sh:37
double y1
amplitude at match point [npe]
double slope
slope [ns/npe]
double JDETECTOR::JPMTAnalogueSignalProcessor::getDecayTime ( ) const
inline

Get decay time.

Returns
decay time [ns]

Definition at line 160 of file JPMTAnalogueSignalProcessor.hh.

161  {
162  return decayTime_ns;
163  }
double JDETECTOR::JPMTAnalogueSignalProcessor::getT1 ( ) const
inline

Get time at transition point from Gaussian to exponential.

Returns
time [ns]

Definition at line 171 of file JPMTAnalogueSignalProcessor.hh.

172  {
173  return t1;
174  }
double JDETECTOR::JPMTAnalogueSignalProcessor::getY1 ( ) const
inline

Get amplitude at transition point from Gaussian to exponential.

Returns
amplitude [npe]

Definition at line 182 of file JPMTAnalogueSignalProcessor.hh.

183  {
184  return y1;
185  }
double y1
amplitude at match point [npe]
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 193 of file JPMTAnalogueSignalProcessor.hh.

194  {
195  return x1;
196  }
double x1
Transition point from a logarithmic to a linear relation between time-over-threshold and number of ph...
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 205 of file JPMTAnalogueSignalProcessor.hh.

206  {
207  if (t1_ns < t1) {
208 
209  const double x = t1_ns / riseTime_ns;
210 
211  return exp(-0.5*x*x); // Gaussian
212 
213  } else {
214 
215  const double x = t1_ns / decayTime_ns;
216 
217  return exp(-x) / y1; // exponential
218  }
219  }
then set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))*exp(-0.5 *(y-[1])*(y-[1])/([2]*[2]))" JF2 -o $WORKDIR/f2.root -F "$FORMULA" -@ "p0
double riseTime_ns
rise time of analogue pulse [ns]
double y1
amplitude at match point [npe]
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 230 of file JPMTAnalogueSignalProcessor.hh.

231  {
232  return riseTime_ns * sqrt(2.0*log(npe/th)); // Gaussian
233  }
double riseTime_ns
rise time of analogue pulse [ns]
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 244 of file JPMTAnalogueSignalProcessor.hh.

245  {
246  if (npe*y1 > th)
247  return decayTime_ns * (log(npe/th) - log(y1)); // exponential
248  else
249  return riseTime_ns * sqrt(2.0*log(npe/th)); // Gaussian
250  }
double riseTime_ns
rise time of analogue pulse [ns]
double y1
amplitude at match point [npe]
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 261 of file JPMTAnalogueSignalProcessor.hh.

262  {
263  if (th > 0.0 && th < 1.0)
264  return 0.5 * TIME_OVER_THRESHOLD_NS / sqrt(-2.0*log(th));
265  else
266  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getMaximalRiseTime(): Invalid threshold " << th);
267  }
const double TIME_OVER_THRESHOLD_NS
Specification for time-over-threshold corresponding to a one photo-electron pulse.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:670
double JDETECTOR::JPMTAnalogueSignalProcessor::applySaturation ( const double  T) const
inline

Get time-over-threshold with saturation.

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

Definition at line 276 of file JPMTAnalogueSignalProcessor.hh.

277  {
278  return saturation / sqrt(T*T + saturation*saturation) * T;
279  }
do set_variable OUTPUT_DIRECTORY $WORKDIR T
double saturation
saturation [ns]
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 288 of file JPMTAnalogueSignalProcessor.hh.

289  {
290  return saturation / sqrt(saturation*saturation - tot_ns*tot_ns) * tot_ns;
291  }
double saturation
saturation [ns]
double JDETECTOR::JPMTAnalogueSignalProcessor::getDerivativeOfSaturation ( const double  T) const
inline

Get derivative of saturation factor.

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

Definition at line 300 of file JPMTAnalogueSignalProcessor.hh.

301  {
303  }
do set_variable OUTPUT_DIRECTORY $WORKDIR T
double saturation
saturation [ns]
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 312 of file JPMTAnalogueSignalProcessor.hh.

313  {
314  return sqrt((double) NPE * gain) * gainSpread;
315  }
double gain
gain [unit]
double gainSpread
gain spread [unit]
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 325 of file JPMTAnalogueSignalProcessor.hh.

326  {
327  if (applyThreshold(npe) > BELOW_THRESHOLDBAND) {
328 
329  const double f = gainSpread * gainSpread;
330 
331  double norm = 0.0;
332  double prob = 0.0;
333  double weight = (PunderAmplified > 0.0 ? pow(PunderAmplified, NPE) : 1.0);
334 
335  for (int k = (PunderAmplified > 0.0 ? 0 : NPE); k <= NPE; ++k) {
336 
337  const double mu = ((NPE-k) * f * gain) + (k * gain);
338  const double sigma = sqrt((NPE-k) * f + k) * getGainSpread(1);
339 
340  norm += weight * (0.5 * TMath::Erfc(((threshold-thresholdBand) - mu) / sqrt(2.0) / sigma));
341  prob += weight * TMath::Gaus(npe, mu, sigma, kTRUE);
342 
343  weight *= ( (NPE - k) / ((double) (k+1)) *
344  (1.0 - PunderAmplified) / PunderAmplified );
345  }
346 
347  return prob / norm;
348  }
349 
350  return 0.0;
351  }
double gain
gain [unit]
double gainSpread
gain spread [unit]
then fatal No sound hydrophone file $HYDROPHONE_TXT fi JGraph f $HYDROPHONE_TXT o $HYDROPHONE_ROOT sort gr k
do set_array DAQHEADER JPrintDAQHeader f
Definition: JTuneHV.sh:79
double threshold
threshold [npe]
double PunderAmplified
probability of underamplified hit
double thresholdBand
threshold-band [npe]
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:98
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.
virtual JThresholdDomains applyThreshold(const double npe) const override
Apply threshold.
std::vector< double > weight
Definition: JAlgorithm.hh:428
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 362 of file JPMTAnalogueSignalProcessor.hh.

363  {
364  double zmin = xmin;
365  double zmax = xmax;
366 
367  const double th = threshold - thresholdBand;
368  const double f = gainSpread * gainSpread;
369 
370  if (zmin < th) { zmin = th; }
371  if (zmax < th) { zmax = th; }
372 
373  double norm = 0.0;
374  double cumulP = 0.0;
375  double weight = (PunderAmplified > 0.0 ? pow(PunderAmplified, NPE) : 1.0);
376 
377  for (int k = (PunderAmplified > 0.0 ? 0 : NPE); k <= NPE; ++k) {
378 
379  const double mu = ((NPE-k) * f * gain) + (k * gain);
380  const double sigma = sqrt((NPE-k) * f + k) * getGainSpread(1);
381 
382  norm += weight * (0.5 * TMath::Erfc((th - mu) / sqrt(2.0) / sigma));
383  cumulP += weight * (0.5 * TMath::Erfc((zmin - mu) / sqrt(2.0) / sigma) -
384  0.5 * TMath::Erfc((zmax - mu) / sqrt(2.0) / sigma));
385 
386  weight *= ( (NPE - k) / ((double) (k+1)) *
387  (1.0 - PunderAmplified) / PunderAmplified);
388  }
389 
390  return cumulP / norm;
391  }
double gain
gain [unit]
double gainSpread
gain spread [unit]
then fatal No sound hydrophone file $HYDROPHONE_TXT fi JGraph f $HYDROPHONE_TXT o $HYDROPHONE_ROOT sort gr k
do set_array DAQHEADER JPrintDAQHeader f
Definition: JTuneHV.sh:79
double threshold
threshold [npe]
double PunderAmplified
probability of underamplified hit
double thresholdBand
threshold-band [npe]
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:98
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.
std::vector< double > weight
Definition: JAlgorithm.hh:428
double JDETECTOR::JPMTAnalogueSignalProcessor::getIntegralOfChargeProbability ( JThresholdDomains  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 401 of file JPMTAnalogueSignalProcessor.hh.

402  {
403  switch (domain) {
404  case ABOVE_THRESHOLD:
406 
407  case IN_THRESHOLDBAND:
409 
410  default:
411  return 0.0;
412  }
413  }
double getIntegralOfChargeProbability(const double xmin, const double xmax, const int NPE) const
Get integral of probability.
double threshold
threshold [npe]
double thresholdBand
threshold-band [npe]
void JDETECTOR::JPMTAnalogueSignalProcessor::setPMTParameters ( const JPMTParameters parameters)
inline

Set PMT parameters.

Parameters
parametersPMT parameters

Definition at line 421 of file JPMTAnalogueSignalProcessor.hh.

422  {
423  static_cast<JPMTParameters&>(*this).setPMTParameters(parameters);
424 
425  configure();
426  }
*fatal Wrong number of arguments esac JCookie sh typeset Z DETECTOR typeset Z SOURCE_RUN typeset Z TARGET_RUN set_variable PARAMETERS_FILE $WORKDIR parameters
Definition: diff-Tuna.sh:38
void configure()
Configure internal parameters.
JPMTParameters()
Default constructor.
virtual bool JDETECTOR::JPMTAnalogueSignalProcessor::applyQE ( ) const
inlineoverridevirtual

Apply relative QE.

Returns
true if accepted; false if rejected

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 451 of file JPMTAnalogueSignalProcessor.hh.

452  {
453  if (QE <= 0.0)
454  return false;
455  else if (QE < 1.0)
456  return gRandom->Rndm() < QE;
457  else
458  return true;
459  }
double QE
relative quantum efficiency
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 468 of file JPMTAnalogueSignalProcessor.hh.

469  {
470  if (TTS_ns < 0.0)
471  return t_ns + getTransitionTime(gRandom->Rndm(), -lrint(TTS_ns));
472  else if (TTS_ns > 0.0)
473  return gRandom->Gaus(t_ns, TTS_ns);
474  else
475  return t_ns;
476  }
getTransitionTime
Function object to generate transition time.
double TTS_ns
transition time spread [ns]
virtual bool JDETECTOR::JPMTAnalogueSignalProcessor::compare ( const JPhotoElectron first,
const JPhotoElectron second 
) const
inlineoverridevirtual

Compare (arrival times of) photo-electrons.

Parameters
firstfirst photo-electron
secondsecond photo-electron
Returns
true if arrival times of photo-electrons are within rise time of analogue pulses; else false

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 486 of file JPMTAnalogueSignalProcessor.hh.

487  {
488  return second.t_ns < first.t_ns + riseTime_ns;
489  }
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
double riseTime_ns
rise time of analogue pulse [ns]
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 498 of file JPMTAnalogueSignalProcessor.hh.

499  {
500  double q;
501  do {
502 
503  // Determine which contribution to sample from
504  int k = NPE;
505  if (PunderAmplified > 0.0) {
506 
507  const double X = gRandom->Uniform();
508  double weight = pow(1.0 - PunderAmplified, NPE);
509 
510  for (double sum_p = 0.0; k > 0; --k) {
511 
512  sum_p += weight;
513  if (sum_p > X) { break; }
514 
515  weight *= ( k / ((double) (NPE - (k-1))) *
516  PunderAmplified / (1.0 - PunderAmplified) );
517  }
518  }
519 
520  // Sample from chosen contribution
521  const double f = gainSpread * gainSpread;
522  const double mu = ((NPE-k) * f * gain) + (k * gain);
523  const double sigma = sqrt((NPE-k) * f + k) * getGainSpread(1);
524 
525  q = gRandom->Gaus(mu,sigma);
526 
527  } while (q < 0.0);
528 
529  return q;
530  }
double gain
gain [unit]
double gainSpread
gain spread [unit]
then fatal No sound hydrophone file $HYDROPHONE_TXT fi JGraph f $HYDROPHONE_TXT o $HYDROPHONE_ROOT sort gr k
do set_array DAQHEADER JPrintDAQHeader f
Definition: JTuneHV.sh:79
double PunderAmplified
probability of underamplified hit
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:98
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.
std::vector< double > weight
Definition: JAlgorithm.hh:428
virtual double JDETECTOR::JPMTAnalogueSignalProcessor::getNPE ( const double  tot_ns,
const double  eps = 1.0e-3 
) const
inlinevirtual

Get number of photo-electrons.

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

Reimplemented from JDETECTOR::JPMTSignalProcessorInterface.

Definition at line 540 of file JPMTAnalogueSignalProcessor.hh.

542  {
543  const double t = removeSaturation(tot_ns);
544 
545  double npe;
546  if (t <= 2*getRiseTime(threshold/y1,threshold)) { // Gaussian + Gaussian
547 
548  npe = threshold * exp(t*t/riseTime_ns/riseTime_ns/8.0);
549 
550  } else if (t <= (getRiseTime(getStartOfLinearisation(), threshold) +
552  // Gaussian + Exponential
553 
554  const double A = decayTime_ns;
555  const double B = sqrt(2.0) * riseTime_ns;
556  const double C = -(decayTime_ns*log(y1) + t);
557  const double z = (-B + sqrt(B*B - 4*A*C)) / (2*A);
558 
559  npe = threshold * exp(z*z);
560 
561  } else { // linear
562 
563  const double T = (getRiseTime(getStartOfLinearisation(), threshold) +
565 
566  npe = getStartOfLinearisation() + (t - T) / slope;
567  }
568 
569  return (npe <= MAX_CHARGE ? npe : MAX_CHARGE);
570  }
double threshold
threshold [npe]
then set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))*exp(-0.5 *(y-[1])*(y-[1])/([2]*[2]))" JF2 -o $WORKDIR/f2.root -F "$FORMULA" -@ "p0
double riseTime_ns
rise time of analogue pulse [ns]
static const double C
Physics constants.
do set_variable OUTPUT_DIRECTORY $WORKDIR T
double removeSaturation(const double tot_ns) const
Get time-over-threshold without saturation.
double getStartOfLinearisation() const
Get transition point from a model-dependent to linear relation between time-over-threshold and number...
source $JPP_DIR setenv csh $JPP_DIR eval JShellParser o a A
double y1
amplitude at match point [npe]
double slope
slope [ns/npe]
double getRiseTime(const double npe, const double th) const
Get time to pass from threshold to top of analogue pulse.
virtual JThresholdDomains 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 579 of file JPMTAnalogueSignalProcessor.hh.

580  {
581  if (npe > threshold) {
582 
583  return ABOVE_THRESHOLD;
584 
585  } else if (npe > threshold - thresholdBand) {
586 
587  return IN_THRESHOLDBAND;
588 
589  } else {
590 
591  return BELOW_THRESHOLDBAND;
592  }
593  }
double threshold
threshold [npe]
double thresholdBand
threshold-band [npe]
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 604 of file JPMTAnalogueSignalProcessor.hh.

605  {
606  if (slewing) {
607 
608  switch (applyThreshold(npe)) {
609  case IN_THRESHOLDBAND:
610  return ((getRiseTime(npe, getTH0()) - getRiseTime(npe, threshold-thresholdBand)) -
611  (getRiseTime(1.0, getTH0()) - getRiseTime(1.0, threshold-thresholdBand))) + this->mean_ns;
612 
613  case ABOVE_THRESHOLD:
614  return ((getRiseTime(npe, getTH0()) - getRiseTime(npe, threshold)) -
615  (getRiseTime(1.0, getTH0()) - getRiseTime(1.0, threshold)));
616 
617  default:
618  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getRiseTime: Invalid charge " << npe);
619  }
620 
621  } else {
622 
623  return 0.0;
624  }
625  }
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:670
double threshold
threshold [npe]
static double getTH0()
Get lower threshold for rise time evaluation.
double thresholdBand
threshold-band [npe]
bool slewing
time slewing of analogue signal
double mean_ns
mean time-over-threshold of threshold-band hits [ns]
virtual JThresholdDomains applyThreshold(const double npe) const override
Apply threshold.
double getRiseTime(const double npe, const double th) const
Get time to pass from threshold to top of analogue pulse.
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 634 of file JPMTAnalogueSignalProcessor.hh.

635  {
636 
637  switch (applyThreshold(npe)) {
638  case IN_THRESHOLDBAND: {
639 
640  return gRandom->Gaus(mean_ns, sigma_ns);
641  }
642 
643  case ABOVE_THRESHOLD: {
644 
645  double tot = 0.0;
646 
647  if (npe*y1 <= threshold) {
648 
649  tot += getRiseTime(npe, threshold); // Gaussian
650  tot += getRiseTime(npe, threshold); // Gaussian
651 
652  } else if (npe <= getStartOfLinearisation()) {
653 
654  tot += getRiseTime (npe, threshold); // Gaussian
655  tot += getDecayTime(npe, threshold); // exponential
656 
657  } else {
658 
659  tot += getRiseTime (getStartOfLinearisation(), threshold); // Gaussian
660  tot += getDecayTime(getStartOfLinearisation(), threshold); // exponential
661 
662  tot += slope * (npe - getStartOfLinearisation()); // linear
663  }
664 
665  return applySaturation(tot);
666  }
667 
668  default: {
669 
670  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getTimeOverThreshold: Invalid charge " << npe);
671  }}
672  }
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:670
double threshold
threshold [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 applySaturation(const double T) 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...
virtual JThresholdDomains applyThreshold(const double npe) const override
Apply threshold.
double y1
amplitude at match point [npe]
double slope
slope [ns/npe]
double getRiseTime(const double npe, const double th) const
Get time to pass from threshold to top of analogue pulse.
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 681 of file JPMTAnalogueSignalProcessor.hh.

682  {
683  switch(applyThreshold(npe)) {
684  case ABOVE_THRESHOLD: {
685 
686  const double z = riseTime_ns / sqrt(2.0 * log(npe/threshold));
687 
688  double y = 0.0;
689 
690  if (npe*y1 > threshold) {
691 
692  if (npe <= getStartOfLinearisation())
693  y = npe / (z + decayTime_ns); // Gaussian + exponential
694  else
695  y = 1.0 / slope; // linear
696 
697  } else {
698 
699  y = npe / (2.0 * z); // Gaussian + Gaussian
700  }
701 
702  const double tot_ns = getTimeOverThreshold(npe);
703 
704  return y / getDerivativeOfSaturation(removeSaturation(tot_ns));
705 
706  }
707 
708  default:
709  return 0.0;
710  }
711  }
double getDerivativeOfSaturation(const double T) const
Get derivative of saturation factor.
double threshold
threshold [npe]
double riseTime_ns
rise time of analogue pulse [ns]
virtual double getTimeOverThreshold(const double npe) const override
Get time-over-threshold (ToT).
double removeSaturation(const double tot_ns) const
Get time-over-threshold without saturation.
double getStartOfLinearisation() const
Get transition point from a model-dependent to linear relation between time-over-threshold and number...
virtual JThresholdDomains applyThreshold(const double npe) const override
Apply threshold.
double y1
amplitude at match point [npe]
double slope
slope [ns/npe]
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 721 of file JPMTAnalogueSignalProcessor.hh.

722  {
723 
724  const double PthBand = getIntegralOfChargeProbability(IN_THRESHOLDBAND, NPE);
725 
726  const double npe = getNPE(tot_ns);
727  const double y = getChargeProbability(npe, NPE);
728  const double v = getDerivative(npe);
729 
730  const double RthBand = PthBand * TMath::Gaus(tot_ns, mean_ns, sigma_ns, kTRUE);
731  const double RaboveTh = y * v;
732 
733  return RthBand + RaboveTh;
734  }
double getIntegralOfChargeProbability(const double xmin, const double xmax, const int NPE) const
Get integral of probability.
virtual double getDerivative(const double npe) const override
Get derivative of number of photo-electrons to time-over-threshold.
double mean_ns
mean time-over-threshold of threshold-band hits [ns]
virtual double getChargeProbability(const double npe, const int NPE) const override
Get probability density for given charge.
double sigma_ns
time-over-threshold standard deviation of threshold-band hits [ns]
virtual double getNPE(const double tot_ns, const double eps=1.0e-3) const
Get number of photo-electrons.
data_type v[N+1][M+1]
Definition: JPolint.hh:740
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 745 of file JPMTAnalogueSignalProcessor.hh.

746  {
747 
748  const double PthBand = getIntegralOfChargeProbability(IN_THRESHOLDBAND, NPE);
749 
750  const double IthBand = PthBand * (0.5 * TMath::Erfc((Tmin - mean_ns) / sqrt(2.0) / sigma_ns) -
751  0.5 * TMath::Erfc((Tmax - mean_ns) / sqrt(2.0) / sigma_ns));
752  const double IaboveTh = getIntegralOfChargeProbability(getNPE(Tmin), getNPE(Tmax), NPE);
753 
754  return IthBand + IaboveTh;
755  }
double getIntegralOfChargeProbability(const double xmin, const double xmax, const int NPE) const
Get integral of probability.
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]
virtual double getNPE(const double tot_ns, const double eps=1.0e-3) const
Get number of photo-electrons.
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 764 of file JPMTAnalogueSignalProcessor.hh.

765  {
766  if (QE <= 0.0) {
767 
768  return 0.0;
769 
770  } else if (QE <= 1.0) {
771 
772  double P = 0.0;
773 
774  const double f = gainSpread * gainSpread;
775 
776  for (int i = 1; i <= NPE; ++i) {
777 
778  const double p = (TMath::Binomial(NPE, i) *
779  TMath::Power(QE, i) * TMath::Power(1.0 - QE, NPE - i));
780 
781  double Ptotal = 0.0;
782  double Pabove = 0.0;
783  double weight = (PunderAmplified > 0.0 ? pow(PunderAmplified, i) : 1.0);
784 
785  for (int k = (PunderAmplified > 0.0 ? 0 : NPE); k <= NPE; ++k) {
786 
787  const double mu = ((i-k) * f * gain) + (k * gain);
788  const double sigma = sqrt((i-k) * f + k) * getGainSpread(1);
789 
790  Ptotal += weight * (0.5 * TMath::Erfc(( 0.0 - mu) / sqrt(2.0) / sigma));
791  Pabove += weight * (0.5 * TMath::Erfc((threshold - mu) / sqrt(2.0) / sigma));
792 
793  weight *= ( (i - k) / ((double) (k+1)) *
794  (1.0 - PunderAmplified) / PunderAmplified );
795  }
796 
797  P += p*Pabove/Ptotal;
798  }
799 
800  return P;
801 
802  } else {
803 
804  THROW(JValueOutOfRange, "JPMTAnalogueSignalProcessor::getSurvivalProbability: Invalid QE " << QE);
805  }
806  }
double gain
gain [unit]
double gainSpread
gain spread [unit]
then fatal No sound hydrophone file $HYDROPHONE_TXT fi JGraph f $HYDROPHONE_TXT o $HYDROPHONE_ROOT sort gr k
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:670
do set_array DAQHEADER JPrintDAQHeader f
Definition: JTuneHV.sh:79
double threshold
threshold [npe]
double PunderAmplified
probability of underamplified hit
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:98
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.
double QE
relative quantum efficiency
std::vector< double > weight
Definition: JAlgorithm.hh:428
then $DIR JPlotNPE PDG P
Definition: JPlotNPE-PDG.sh:60
static double JDETECTOR::JPMTAnalogueSignalProcessor::getTH0 ( )
inlinestatic

Get lower threshold for rise time evaluation.

Returns
threshold [npe]

Definition at line 814 of file JPMTAnalogueSignalProcessor.hh.

815  {
816  return 0.1;
817  }
static double JDETECTOR::JPMTAnalogueSignalProcessor::getTH1 ( )
inlinestatic

Get upper threshold for rise time evaluation.

Returns
threshold [npe]

Definition at line 825 of file JPMTAnalogueSignalProcessor.hh.

826  {
827  return 0.9;
828  }
void JDETECTOR::JPMTSignalProcessorInterface::operator() ( const JCalibration calibration,
const JPMTData< JPMTSignal > &  input,
JPMTData< JPMTPulse > &  output 
) const
inlineinherited

Process hits.

Parameters
calibrationPMT calibration
inputPMT signals
outputPMT hits

Definition at line 103 of file JPMTSignalProcessorInterface.hh.

106  {
107  // apply transition time distribution to each photo-electron.
108 
109  buffer.clear();
110 
111  for (JPMTData<JPMTSignal>::const_iterator hit = input.begin(); hit != input.end(); ++hit) {
112  for (int i = 0; i != hit->npe; ++i) {
113  if (applyQE()) {
114  buffer.push_back(JPhotoElectron(getRandomTime(hit->t_ns)));
115  }
116  }
117  }
118 
119  if (!buffer.empty()) {
120 
122 
123  buffer.sort();
124 
125 
126  // generate PMT hits from time sequence of photo-electrons.
127 
128  for (JPMTData<JPhotoElectron>::const_iterator q = buffer.begin(), p = q++; q != buffer.end(); ++q) {
129 
130  while (compare(*p,*q)) {
131  ++q;
132  }
133 
134  const double npe = getRandomCharge(distance(p,q));
135  if (applyThreshold(npe) > BELOW_THRESHOLDBAND) {
136  output.push_back(JPMTPulse(putTime(p->t_ns + getRiseTime(npe), calibration), getTimeOverThreshold(npe)));
137  }
138 
139  p = q;
140  }
141 
142  // merge overlapping PMT hits.
143 
144  merge(output);
145  }
146  }
Data structure for PMT digital pulse.
virtual double getRandomCharge(const int NPE) const
Get randomised charge in accordance with 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.
virtual void merge(JPMTData< JPMTPulse > &data) const
Merging of PMT hits.
virtual bool applyQE() const
Apply relative QE.
virtual double getRiseTime(const double npe) const
Get time to pass threshold.
double putTime(const T &t1, const JCalibration &cal)
Get de-calibrated time.
virtual JThresholdDomains applyThreshold(const double npe) const
Apply threshold.
virtual bool compare(const JPhotoElectron &first, const JPhotoElectron &second) const
Compare arrival times of photo-electrons.
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).
std::vector< JElement_t >::const_iterator const_iterator
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

Definition at line 294 of file JPMTSignalProcessorInterface.hh.

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

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

Returns
minimal time [ns]

Definition at line 380 of file JPMTSignalProcessorInterface.hh.

381  {
382  return 1.0;
383  }
static double JDETECTOR::JPMTSignalProcessorInterface::getQmin ( )
inlinestaticinherited

Get minimal width of charge distribution.

Returns
minimal charge [npe]

Definition at line 391 of file JPMTSignalProcessorInterface.hh.

392  {
393  return 1.0e-3;
394  }
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  }
Data structure for PMT parameters.
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  }
double gain
gain [unit]
double gainSpread
gain spread [unit]
double threshold
threshold [npe]
double thresholdBand
threshold-band [npe]
double QE
relative quantum efficiency
static JEquationParameters& JDETECTOR::JPMTParameters::getEquationParameters ( )
inlinestaticinherited

Get equation parameters.

Returns
equation parameters

Definition at line 192 of file JPMTParameters.hh.

193  {
194  static JEquationParameters equation("=", ",", "./", "#");
195 
196  return equation;
197  }
Simple data structure to support I/O of equations (see class JLANG::JEquation).
static void JDETECTOR::JPMTParameters::setEquationParameters ( const JEquationParameters equation)
inlinestaticinherited

Set equation parameters.

Parameters
equationequation parameters

Definition at line 205 of file JPMTParameters.hh.

206  {
207  getEquationParameters() = equation;
208  }
static JEquationParameters & getEquationParameters()
Get equation parameters.
JProperties JDETECTOR::JPMTParameters::getProperties ( const JEquationParameters equation = JPMTParameters::getEquationParameters())
inlineinherited

Get properties of this class.

Parameters
equationequation parameters

Definition at line 216 of file JPMTParameters.hh.

217  {
218  return JPMTParametersHelper(*this, equation);
219  }
JProperties JDETECTOR::JPMTParameters::getProperties ( const JEquationParameters equation = JPMTParameters::getEquationParameters()) const
inlineinherited

Get properties of this class.

Parameters
equationequation parameters

Definition at line 227 of file JPMTParameters.hh.

228  {
229  return JPMTParametersHelper(*this, equation);
230  }

Friends And Related Function Documentation

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 436 of file JPMTAnalogueSignalProcessor.hh.

437  {
438  in >> static_cast<JPMTParameters&>(object);
439 
440  object.configure();
441 
442  return in;
443  }
JPMTParameters()
Default constructor.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY source JAcoustics sh $DETECTOR_ID typeset A TRIPODS get_tripods $WORKDIR tripod txt TRIPODS for EMITTER in
Definition: JCanberra.sh:36

Member Data Documentation

double JDETECTOR::JPMTAnalogueSignalProcessor::decayTime_ns
protected

decay time [ns]

Definition at line 832 of file JPMTAnalogueSignalProcessor.hh.

double JDETECTOR::JPMTAnalogueSignalProcessor::t1
protected

time at match point [ns]

Definition at line 833 of file JPMTAnalogueSignalProcessor.hh.

double JDETECTOR::JPMTAnalogueSignalProcessor::y1
protected

amplitude at match point [npe]

Definition at line 834 of file JPMTAnalogueSignalProcessor.hh.

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 840 of file JPMTAnalogueSignalProcessor.hh.

double JDETECTOR::JPMTSignalProcessorInterface::MAX_CHARGE
protectedinherited

Maximum charge [npe].

Definition at line 398 of file JPMTSignalProcessorInterface.hh.

double JDETECTOR::JPMTParameters::QE
inherited

relative quantum efficiency

Definition at line 233 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::gain
inherited

gain [unit]

Definition at line 234 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::gainSpread
inherited

gain spread [unit]

Definition at line 235 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::riseTime_ns
inherited

rise time of analogue pulse [ns]

Definition at line 236 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::TTS_ns
inherited

transition time spread [ns]

Definition at line 237 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::threshold
inherited

threshold [npe]

Definition at line 238 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::PunderAmplified
inherited

probability of underamplified hit

Definition at line 239 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::thresholdBand
inherited

threshold-band [npe]

Definition at line 240 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::mean_ns
inherited

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

Definition at line 241 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::sigma_ns
inherited

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

Definition at line 242 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::slope
inherited

slope [ns/npe]

Definition at line 243 of file JPMTParameters.hh.

double JDETECTOR::JPMTParameters::saturation
inherited

saturation [ns]

Definition at line 244 of file JPMTParameters.hh.

bool JDETECTOR::JPMTParameters::slewing
inherited

time slewing of analogue signal

Definition at line 245 of file JPMTParameters.hh.


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