1 #ifndef __JDETECTOR__JPMTANALOGUESIGNALPROCESSOR__
2 #define __JDETECTOR__JPMTANALOGUESIGNALPROCESSOR__
77 x1(std::numeric_limits<double>::max())
96 static const int N = 100;
97 static const double precision = 1.0e-4;
118 const double Q = b*b - 4.0*a*
c;
138 x1 = std::numeric_limits<double>::max();
143 for (
int i = 0;
i !=
N; ++
i) {
145 const double x = 0.5 * (xmin +
xmax);
148 if (fabs(1.0 - u) < precision) {
220 return exp(-0.5*x*x);
272 if (th > 0.0 && th < 1.0)
275 THROW(
JValueOutOfRange,
"JPMTAnalogueSignalProcessor::getMaximalRiseTime(): Invalid threshold " << th);
302 return std::numeric_limits<double>::max();
347 if (zmin < th) { zmin = th; }
348 if (zmax < th) { zmax = th; }
356 const double mu = ((NPE-
k) * f *
gain) + (
k *
gain);
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));
363 weight *= ( (NPE -
k) / ((
double) (
k+1)) *
367 return cumulP / norm;
457 return gRandom->Rndm() <
QE;
474 return gRandom->Gaus(t_ns,
TTS_ns);
513 const double X = gRandom->Uniform();
516 for (
double sum_p = 0.0; k > 0; --
k) {
519 if (sum_p > X) {
break; }
521 weight *= ( k / ((double) (NPE - (k-1))) *
528 const double mu = ((NPE-
k) * f *
gain) + (k *
gain);
531 q = gRandom->Gaus(mu,sigma);
558 const double mu = ((NPE-
k) * f *
gain) + (
k *
gain);
562 prob += weight * TMath::Gaus(npe, mu, sigma, kTRUE);
564 weight *= ( (NPE -
k) / ((
double) (
k+1)) *
662 THROW(
JValueOutOfRange,
"JPMTAnalogueSignalProcessor::getTimeOverThreshold: Invalid charge " << npe);
718 }
else if (
QE <= 1.0) {
724 for (
int i = 1;
i <= NPE; ++
i) {
726 const double p = (TMath::Binomial(NPE,
i) *
727 TMath::Power(
QE,
i) * TMath::Power(1.0 -
QE, NPE -
i));
738 Ptotal += weight * (0.5 * TMath::Erfc((0.0 - mu) / sqrt(2.0) / sigma));
741 weight *= ( (
i -
k) / ((
double) (
k+1)) *
745 P += p*Pabove/Ptotal;
763 virtual double getNPE(
const double tot_ns)
const override
766 return std::numeric_limits<double>::max();
777 }
else if (tot <= TOT) {
782 const double z = (-b + sqrt(b*b - 4*a*c)) / (2*a);
804 const double npe =
getNPE(tot_ns);
808 const double RthBand = PthBand * TMath::Gaus(tot_ns,
mean_ns,
sigma_ns, kTRUE);
809 const double RaboveTh = y *
v;
811 return RthBand + RaboveTh;
828 const double IthBand = PthBand * (0.5 * TMath::Erfc((Tmin -
mean_ns) / sqrt(2.0) /
sigma_ns) -
832 return IthBand + IaboveTh;
883 const double precision = 1.0e-4)
885 int i = (int) (NPE - 5.0 * sqrt(NPE) - 0.5);
886 int M = (int) (NPE + 5.0 * sqrt(NPE) + 0.5);
888 if (i < 1) { i = 1; }
889 if (M < i) { M =
i; }
891 double p = NPE *
exp(-NPE) / (double) 1;
893 for (
int __i = 1; __i !=
i; ++__i) {
899 for (
double p0 = 0.0; (p >= p0 || p > precision) && i != M; ++
i, p0 = p, p *= NPE / (double) i) {
virtual double getRandomCharge(const int NPE) const override
Get randomised charge according to gain and gain spread.
double getDerivativeOfSaturation(const double tot_ns) const
Get derivative of saturation factor.
then fatal No hydrophone data file $HYDROPHONE_TXT fi sort gr k
Q(UTCMax_s-UTCMin_s)-livetime_s
virtual double getRiseTime(const double npe) const override
Get time to reach threshold.
double getIntegralOfChargeProbability(const JThresholdDomain domain, const int NPE) const
Get integral of probability in specific threshold domain.
double x1
Transition point from a logarithmic to a linear relation between time-over-threshold and number of ph...
void setPMTParameters(const JPMTParameters ¶meters)
Set PMT parameters.
double getTimeOverThresholdProbability(const JPMTAnalogueSignalProcessor &pmt, const double tot_ns, const double NPE, const double precision=1.0e-4)
Get time-over-threshold probability.
Time calibration (including definition of sign of time offset).
double getIntegralOfChargeProbability(const double xmin, const double xmax, const int NPE) const
Get integral of probability.
double decayTime_ns
decay time [ns]
double gainSpread
gain spread [unit]
Data structure for single photo-electron.
const double TIME_OVER_THRESHOLD_NS
Specification for time-over-threshold corresponding to a one photo-electron pulse.
void setPMTParameters(const JPMTParameters ¶meters)
Set PMT parameters.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
o $QUALITY_ROOT d $DEBUG!CHECK_EXIT_CODE JPlot1D f
*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
double threshold
threshold [npe]
static double getTH0()
Get lower threshold for rise time evaluation.
double getY1() const
Get amplitude at transition point from Gaussian to exponential.
double PunderAmplified
probability of underamplified hit
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
double thresholdBand
threshold-band [npe]
JDETECTOR::getTransitionTime getTransitionTime
Function object to generate transition time.
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.
bool slewing
time slewing of analogue signal
JThresholdDomain getThresholdDomain(const double npe) const
Get threshold domain.
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 applySaturation(const double tot_ns) const
Get time-over-threshold with saturation.
static double getMaximalRiseTime(const double th)
Get maximal rise time for given threshold.
void configure()
Configure internal parameters.
JPMTAnalogueSignalProcessor(const JPMTParameters ¶meters=JPMTParameters())
Constructor.
double sigma_ns
time-over-threshold standard deviation of threshold-band hits [ns]
friend std::istream & operator>>(std::istream &in, JPMTAnalogueSignalProcessor &object)
Read PMT signal from input.
JPMTParameters()
Default constructor.
T pow(const T &x, const double y)
Power .
virtual double getNPE(const double tot_ns) const override
Get number of photo-electrons.
PMT signal processor interface.
double getAmplitude(const double t1_ns) const
Get amplitude at given time for a one photo-electron pulse.
virtual bool applyThreshold(const double npe) const
Apply threshold.
double getGainSpread(int NPE) const
Get gain spread for given number of photo-electrons.
double getDecayTime() const
Get decay time.
PMT analogue signal processor.
$WORKDIR ev_configure_dqsimulator txt echo process $DQ_SIMULATOR $i $SOURCE_HOST[$index] csh c(setenv ROOTSYS $ROOTSYS &&source $JPP_DIR/setenv.csh $JPP_DIR &&($DQ_SIMULATOR\-u\$NAME\$\-H\$SERVER\$\-M\$LOGGER\$\-d $DEBUG</dev/null > &/dev/null &))'
then set_variable DIR else fatal Wrong number of arguments fi for INPUT_FILE in ls rt $DIR stage * log
virtual double getTimeOverThreshold(const double npe) const override
Get time-over-threshold (ToT).
then usage $script< input file >[option[primary[working directory]]] nWhere option can be N
virtual double getRandomTime(const double t_ns) const override
Get randomised time according transition time distribution.
double TTS_ns
transition time spread [ns]
double getIntegralOfTimeOverThresholdProbability(const double Tmin, const double Tmax, const int NPE) const
Get cumulative probability of time-over-threshold distribution.
double getTimeOverThresholdProbability(const double tot_ns, const int NPE) const
Get probability of having a pulse with specific time-over-threshold.
no fit printf nominal n $STRING awk v X
JThresholdDomain
Threshold domain specifiers.
Exception for accessing a value in a collection that is outside of its range.
double getT1() const
Get time at transition point from Gaussian to exponential.
Data structure for PMT parameters.
double saturation
saturation [ns]
double removeSaturation(const double tot_ns) const
Get time-over-threshold without saturation.
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 JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
static double getTH1()
Get upper threshold for rise time evaluation.
double QE
relative quantum efficiency
virtual bool compare(const JPhotoElectron &first, const JPhotoElectron &second) const override
Compare arrival times of photo-electrons.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable STRING $argv[2] set_array QUANTILES set_variable FORMULA *[0] exp(-0.5 *(x-[1])*(x-[1])/([2]*[2]))" set_variable MODULE `getModule -a $DETECTOR -L "$STRING 0"` source JAcousticsToolkit.sh typeset -A TRIPODS get_tripods $WORKDIR/tripod.txt TRIPODS XMEAN
double getDecayTime(const double npe, const double th) const
Get time to pass from top of analogue pulse to threshold.
virtual bool applyQE() const override
Apply relative QE.
double t1
time at match point [ns]
double getStartOfLinearisation() const
Get transition point from a model-dependent to linear relation between time-over-threshold and number...
virtual double getSurvivalProbability(const int NPE) const override
Probability that a hit survives the simulation of the PMT.
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.