////////////////////////////////////////////////////////////////////////////// // Class: RAT::Methods::PMTTimeDistns // // Brief: Calculate various probabilities related to an observed PMT hit time, // for 0 to fMaxPE possible generated photoelectrons // // Author: W. Heintzelman // // REVISION HISTORY: // 25/05/2016 : W. Heintzelman - New file // 30/04/2018 : W. Heintzelman - extensive revisions, addition of // calculation of distributions of order statistics // // Detail: // Given an input table that defines the pdf of the time a photon spends // in the PMT bucket (i.e., after crossing the PMT face) until a PE is // produced, the distribution of the PMT transit time (i.e., the time // interval between when a PE is produced and the recorded hit time) from // class PMTTransitTime, and a list of tabulation-area photon-crossing // times from the auxiliary simulation, this class computes values of the // distribution of the total time between the occurence of an event and // the time of the electonic pulse corresponding to a single photon. For // a single-PE PMT hit, this would be the recorded hit time. It uses the // fact that the characteristic function of a sum of independent random // variables is the product of the characteristic functions of the // variables, and uses Fast Fourier Transform (FFT) methods to calculate // and invert characteristic functions. // // Notes: Input times are relative either to the time of the auxiliary // simulated event or to the seed fitted event time. Internally in // the pdfs, times are offset by fTimeOffset to allow for the possibility // of somewhat negative recorded hit times. // // This is a singleton class // ////////////////////////////////////////////////////////////////////////////// #ifndef __RAT_PMTTimeDistns_ #define __RAT_PMTTimeDistns_ #include #include #include const int fNPts=4096; // no. of input points to discrete Fourier transforms const int fNOv2P1 = fNPts/2 +1; const double fBinWidth = 0.2; // bin width in ns of input and output time pdfs // Storage for discrete points of a pdf struct DiscretePDFT{ double p[fNPts]; }; // Storage for Fourier coefficients of a discrete characteristic function struct CharFnT{ double real[fNOv2P1]; double imag[fNOv2P1]; }; namespace RAT { class PMTTimeDistns{ public: static PMTTimeDistns* Instance(); // Return pointer to singleton void BeginOfRun(); void SetEventParameters (double estTriggerTime, size_t maxPE); // Following functions have as inputs: // nExpNoisePEs - The auxiliary simulation is executed without // noise. This is the expected number of noise PEs in the // PMT in question during the trigger window if noise had // been modelled in the auxiliary simulation, multiplied by // the expected ratio of the number of photons incident on // PMT tabulation area to the number of photoelectrons // produced by photons. This ratio is calculated in PEnergy // as equal to the tabulation area factor times the number of // initiating tracks in the auxiliary simulation divided by // the average PMT effectiveness factor for the incident // photons. // relHitTime = the hit time of the PMT relative to fitted event time // tCross = list of times photons crossed the tabulation area // of the PMT in the auxiliary simulation, relative to // the start of the MC event. // tWeight = list of weights for the forgoing crossing times // // Return value of the probability density function at the time double GetRelHitTimeProb( double nExpNoisePEs, std::vector& tCross, std::vector& tWeight, double relHitTime, UInt_t pmtIdx); // Return value of the cumulative distn function at the time double GetRelHitTimeCDF( double nExpNoisePEs, std::vector& tCross, std::vector& tWeight, double relHitTime, UInt_t pmtIdx); void CalcOrderStats( double nExpNoisePEs, std::vector& tCross, std::vector& tWeight, double relHitTime, UInt_t pmtIdx, double pWidth, double chgIntT, double tLimit, std::vector< std::vector >& g, std::vector< std::vector >& Gwid, std::vector< std::vector >& Gprm, std::vector< std::vector >& Glim); void CalcIntervalProbs(double nExpNoisePEs, std::vector& tCross, std::vector& tWeight, double relHitTime, UInt_t pmtIdx, double delta, std::vector< std::vector >& pdelta); private: static PMTTimeDistns* fPtrSelf; // Ptr to singleton of this class PMTTimeDistns(); // private because this is a singleton class TFFTRealComplex* fftFwd; // FFT from pdf to characteristic fn. TFFTComplexReal* fftRev; // FFT from char. fn. to pdf CharFnT fDelayCharFn; // characteristic fn of total time between // when photon enters concentrator and recorded hit on FEC. void CalculateTimeDFs( double nExpNoisePEs, std::vector& tCross, std::vector& tWeight); // The value of fEstTriggerTime for the event must have been set by // a call to SetEventParameters before this function is called. void CalcBinomCoeff(int maxn, std::vector< std::vector >& c); void MultiplyCharFn( CharFnT& a, CharFnT& b, CharFnT& result); double fTimeOffset; // See notes in header. double fTriggerWindow; // length of the trigger window (gate) in ns double fTriggerToFECDelay; // delay of trigger signal back to FEC in ns double fEstTriggerTime; // estimated trigger time of auxiliary event // relative to the event start time int fMaxPE; // maximum number of PEs that need be considered double fHighLimit; // Upper limit on the time pdf range double fNoiseLimit; // Upper time limit to which noise hits are added // to the hit time distribution int fLastBin; UInt_t fPrevPMTIdx; // Index of PMT for which class was previously called std::vector< std::vector > fBinomCoeff; // storage for // binomial coeff. DiscretePDFT fHitTimePDF, fHitTimeCDF; // fHitTimePDF.p[iBin] is the amount of probability in bin iBin // fHitTimeCDF.p[iBin] is the sum of the probability from bin 0 // up to and including bin iBin. }; } // namespace RAT #endif