//////////////////////////////////////////////////////////////////// /// \class RAT::PCAProc /// /// \brief Processor for generating constants from the PCA /// /// \author Freija Descamps (fbdescamps@lbl.gov) /// /// REVISION HISTORY: /// /// 01 May 2013 : F. Descamps - First version, laserball data only \n /// 01 Jun 2013 : F. Descamps - Added TELLIE first try \n /// 14 Oct 2014 : F. Descamps - Small cppcheck fix /// 11 Jan 2016 : F. Descamps - Add TELLIE DB access and Dijkstra algoritm /// - Major uncluttering using uncrustigy /// 20 Jul 2016 : M. Mottram - Update for ROOT6 compatibility /// 10 Oct 2016 : N. Barros - Completed the ROOT6 compatibility modifications /// 10 Nov 2016 : F. Descamps - 'Smarter' DB access /// /// \details This class generates PCA constants from the /// laserball or TELLIE calibration runs /// /// //////////////////////////////////////////////////////////////////// #ifndef __RAT_PCAProc__ #define __RAT_PCAProc__ // RAT stuff #include #include #include #include #include #include // C++ stuff #include #include class TFile; namespace RAT { class PCAProc : public Processor { public: PCAProc(); virtual ~PCAProc(); virtual void BeginOfRun( DS::Run& run ); virtual Processor::Result DSEvent( DS::Run& run, DS::Entry& ds ); protected: virtual Processor::Result Event( DS::Entry& ds, DS::EV& ev ); // PCA magic: // 1) /// Extract the interpolation points and fit to high charge points from t_ECA vs QHS data int WalkFitPoints(const int PMTNum, ///< Input PMT ID const int LEDNum, ///< Input LED ID std::vector &interX, ///< Return X-values of interpolation points std::vector &interY, ///< Return Y-values of interpolation points std::vector &interXErr, ///< Return Error on X-values of interpolation points std::vector &interYErr, ///< Return Error on Y-values of interpolation points std::vector &interFit ///< Return gradient, intercept and errors for high Q fit ); // 2) /// Extract charge spectrum QHS and QHL peak, threshold and hhp values. int GainFitPoints(const int pmt_num, ///< Input PMT ID const int chargeFlag, ///< Input charge flag 0=QHS, 1=QHL double &peak_bin_val, ///< Return the position of the peak double &thresh_val, ///< Return the position of the threshold double &hhp_val ///< Return the position of the high-half point ); // 3) /// Extract the time offsets between different LEDs in order to bring all LEDs to a common time-offset. int Bootstrap(); /// Tools /// Dijkstra implementation std::vector Dijkstra(std::vector > &LEDErrmatrix, ///< Input matric holding the error on the LED time-offsets std::vector > &LEDmatrix, ///< Input matric holding the LED time-offsets int startLED, ///< The start LED std::vector &inter, ///< The resulting offset array std::vector &interErr); ///< The resulting offset-error array /// Determine if the file exists static bool FileExists(const std::string &filename); /// Calculate median of input vectors, with errors. Also determine the interpolation point position and rms. void MultiCalc(const int type, ///< Tells the function which type of input to expect const std::vector &arg1, ///< Vector containing Y-values const std::vector &arg2 ///< Vector containing X-values ); // Input/Output: some are adapted from ECAProc.cc /// Read in the user input from PCA_GENERATION.ratdb void ReadInUserInput(); /// Initialise the output rootfile void InitPCARoot(const std::string &rootfilename); /// Initialise, write and close the output .ratdb files void PCATitles(); /// Close the output rootfile void ClosePCARoot(); /// Check data and PCA calibration quality void CheckPCARun(); /// Check that the occupancy is not too high void CheckOccupancy(); /// Handle non-fatal errors: warning but don't exit void DontDie(std::string message, ///< Error message int flag_err, ///< Error number float info1, ///< More info to input in error message float info2 ///< More info to input in error message ); /// Handle fatal errors: error and exit void Die(std::string message, ///< Error message int flag_err, ///< Error number float info1, ///< More info to input in error message float info2 ///< More info to input in error message ); // Flags to define the type of calibration and the verbosity int fPCAflag; int fPCAtype; bool fVerbosePCA; bool fDebugPCA; bool fDoTWPCA; bool fDoGFPCA; // Source related information int fPCASource; ///< Type of light source int fNSource; ///< Number of light-sources int fMaxNumFiber; ///< Maximum number of fibers double fMaxDistFromCenter; ///< For laserball, maximum allowed distance from center double fLaserballWL; ///< Laserball wavelength double fSelectAngle; ///< Maximum angle for fiber // Run related information int fTrigType; ///< Trigger type, should be EXTASYNC double fExtaSyncDelay; ///< Delay between trigger to laser/tellie and trigger to DAQ // Parameters that define the interpolation points int fTotIntPoints; ///< Total number of interpolation points for TW int fPeakIntPoints; ///< Number of interpolation points in peak region int fDropIntPoints; ///< Number of interpolation points in the peak drop-off region int fTailIntPoints; ///< Number of interpolation points in the tail region int fMinPointsPerBin; ///< Minimum required hits per interpolation bin int fMinHitsPerChan; ///< Minimum required hits per channel for TW to be performed std::vector fTailGap; // Minimum and maximum teca double fTimeMinVal; double fTimeMaxVal; // For the lightpath function double fDistTolerance; // For the TELLIE case bootstrap function double fMaxCorr; double fMaxErrLEDTime; /// For the MultiCalc function float fMedian; float fWidth; float fWidth2; float fSig; float fVar; float fGrad; // Run-specific checks, needs to be changed when DS changes bool fIsFirstRun; ///< Determine if we are at the first BeginOfRun bool fIsValid; ///< Determines if PCA calibration should be attempted bool fRunCheckOK; ///< Determines if the run pre-flight checks were done bool fRunCalcTransitTimes; ///< Determines if the transittimes were calculated for the current run bool fIsMC; int fDate; int fRunID; int fStartRunID; int fStopRunID; int fLEDnr; TVector3 fLEDPos; TVector3 fLEDDir; // Number of channels int fTotPMTs; /// Position and direction of the PMTs std::vector fTubePositionX; std::vector fTubePositionY; std::vector fTubePositionZ; std::vector fTubeDirX; std::vector fTubeDirY; std::vector fTubeDirZ; /// Status of the PMT std::vector fTubeStatus; /// Position and direction of the LEDs std::vector fLEDPositionX; std::vector fLEDPositionY; std::vector fLEDPositionZ; std::vector fLEDDirX; std::vector fLEDDirY; std::vector fLEDDirZ; std::vector fLEDWavelength; /// transittimes std::vector< std::vector > fTransitTimes; /// viewing angles std::vector< std::vector > fViewingAngles; /// Vectors to store the PMT accumulated hit information std::vector< std::vector< std::vector > > fQHSArray; std::vector< std::vector< std::vector > > fQHLArray; std::vector< std::vector< std::vector > > fTArray; // Store minimum and maximum number of hits per event int fNHitsMax; int fNHitsMin; // For Gain Fit (GF) results std::vector fPMTstatus; std::vector fGFpeaQHS; std::vector fGFthrQHS; std::vector fGFhhpQHS; std::vector fGFpeaQHL; std::vector fGFthrQHL; std::vector fGFhhpQHL; float fMeanhhpdb; float fMinhhpdb; float fMaxhhpdb; float fMeanhhp; float fMinhhp; float fMaxhhp; // For Gain Fit (GF) database values (later: check also if status of PMT has changed?) std::vector fPMTstatusdb; std::vector fGFpeaQHSdb; std::vector fGFthrQHSdb; std::vector fGFhhpQHSdb; std::vector fGFpeaQHLdb; std::vector fGFthrQHLdb; std::vector fGFhhpQHLdb; // For Gain Fit (GF) fitted values std::vector fGFhhpQHSfit; // For intermediate TimeWalk (TW) results, these are 3D arrays (yuk): Length: NPMTs, width=nSources, Depth= number of interpolation points std::vector< std::vector< std::vector > > fTWinterXQHS; std::vector< std::vector< std::vector > > fTWinterYQHS; // The error is the RMS for that interpolation bin std::vector< std::vector< std::vector > > fTWinterXErrQHS; std::vector< std::vector< std::vector > > fTWinterYErrQHS; // The fit to the high Q tail is a 2D array std::vector< std::vector > fTWinterI; std::vector< std::vector > fTWinterR; std::vector< std::vector > fTWinterIErr; std::vector< std::vector > fTWinterRErr; // A global RMS is returned for each PMT and used to judge the TW quality std::vector< std::vector > fTWRms; // Keep track of minimum and maximum for all charge vectors, this saves time afterwards!! std::vector< std::vector > fQHSMin; std::vector< std::vector > fQHSMax; std::vector< std::vector > fQHLMin; std::vector< std::vector > fQHLMax; // Keep track of how many hits are outside the TAC and QHS window std::vector< std::vector > fTQHSLess; std::vector< std::vector > fTMore; // Keep track of which LED is giving results for which PMT std::vector fSourceSelected; // Status flags BitManip *fBits; PCABits fPCABits; int fOverallMask; std::vector fPMTTWMask; // PMT TW status: 1 int per PMT std::vector fPMTGFMask; // PMT GF status: 1 int per PMT // Output histograms PCAHists fHists; // For the Bootstrap method (TELLIE) std::vector< std::vector > fLEDCheckRATid; // nPMTs x nLEDs std::vector fTOffsetArray; // nLEDs - final result std::vector fTOffsetErrArray; // nLEDs - errors final result std::vector fUsedLEDs; // the LEDs that were used in this round int fMinNPointsBoot; // Minimum required points for the bootstrap double fSelectTTransit; // Selection on transittime // For output std::string fPCALogName; std::string fPCATWDBName; std::string fPCARootName; std::string fPCAGFDBName; TFile *fPCARootFile; // The limits for validity of results int fPeakQHSDiffLim; int fPeakQHLDiffLim; int fThreshQHSDiffLim; int fThreshQHLDiffLim; int fHhpQHSDiffLim; int fHhpQHLDiffLim; int fThreshQHSLimhigh; int fThreshQHLLimhigh; int fThreshQHSLimlow; int fThreshQHLLimlow; int fMaxHhpQHSDiff; int fMaxPeakQHSDiff; int fMaxThreshQHSDiff; int fMaxHhpQHLDiff; int fMaxPeakQHLDiff; int fMaxThreshQHLDiff; int fMaxThreshQHSHigh; int fMaxThreshQHLHigh; int fMaxThreshQHSLow; int fMaxThreshQHLLow; int fMaxNPWLow; int fMaxFailed; int fMaxOffline; int fMaxZeroOcc; int fMaxLowOcc; int fMaxRmsHigh; int fMaxRmsVeryHigh; int fMaxTWnoHighQFit; int fNMaxTStepCut; int fNMaxFLateCut; int fNMaxFOutCut; int fNMaxIntPointMiss; int fNTWMaxGrad; int fNTWMinGrad; double fRmsHigh; double fRmsVeryHigh; double fMaxTStep; double fMaxFLate; double fMaxFOut; double fMaxGrad; double fMinGrad; }; } // namespace RAT #endif