7#include <Math/Minimizer.h>
8#include <Math/Factory.h>
9#include <Math/Functor.h>
10#include <TClonesArray.h>
11#include <TObjString.h>
16#include <TMultiGraph.h>
27#include <TEfficiency.h>
30#include <THistPainter.h>
32#include <TFitResultPtr.h>
33#include <TFitResult.h>
35#include <TPaveStats.h>
82 unsigned int Millisecond;
83 unsigned long long int TDCsamIndex;
87 struct channel_data_without_measurement
91 int FirstCellToPlotsamIndex;
94 struct channel_data_with_measurement
98 int FirstCellToPlotsamIndex;
99 float MeasuredBaseline;
100 float AmplitudeValue;
101 float ComputedCharge;
102 float RiseTimeInstant;
103 float FallTimeInstant;
104 float RawTriggerRate;
123 "Error: No data has been loaded yet!\n \
124 Please call ReadFile() before calling any functions which manipulate data.\n \
141 void PlotChannelSums(
bool =
false,
bool =
false,
double = 0.,
double = 0.,
int = 2);
162 void PrintWFProjection(
float = 0,
float = 320,
float = -50,
float = 50,
int = 200);
168 void GetTimingCFD(
float = .3,
float = 100,
float = 140,
double = 0.,
bool =
true,
int = 2,
bool =
false,
bool =
false);
174 void FilterAll(
double = .3,
double = .9,
double = .2);
183 void PrintChargeSpectrumWF(
float,
float,
float = 0,
float = 300,
int = 1,
float = 0.,
float = 0.,
float = 0.,
float = 0.);
184 TH1F*
ChargeSpectrum(
int,
float,
float,
float = 0,
float = 300,
float = -50,
float = 600,
int = 750);
185 void PrintChargeSpectrum(
float,
float,
float = 0,
float = 300,
float = -50,
float = 600,
int = 750,
float = 0.,
float = 0.,
int = 99,
int = 0,
bool =
false);
190 float*
ChargeList(
int,
float = 20,
float = 80,
float = 0,
float = 300,
bool = 1);
191 void SaveChargeLists(
float = 20,
float = 80,
float = 0,
float = 300,
bool = 1);
192 void ChargeCorrelation(
float,
float,
float,
float,
float,
float,
int,
int,
int,
bool =
false);
195 void PrintDCR(
float = 15,
float = 85,
float = 0,
float = 300,
double = 3);
198 TH1F*
TimeDist(
int,
float = 0,
float = 300,
float = 0,
float = 300,
int = 100,
int = 0,
float = .3);
199 void PrintTimeDist(
float = 0,
float = 300,
float = 0,
float = 300,
int = 100,
int = 0,
float = .3);
204 void Print_GetTimingCFD(
float = 100,
float = 140,
int = 0,
int = -999,
string =
"S",
bool =
true);
215 double*
gety(
int,
int);
241 virtual void ReadFile(
string,
bool =
false,
int = 9,
string =
"out.root",
bool =
false,
long long = -1);
int omp_get_max_threads()
Main class containing the file reader and most analysis functions.
int FirstBinFileToRead
Number first of .bin file to be read in.
static pair< float, bool > LinearInterpolation(float, float, float, float, float, bool=false)
Simple linear interpolation for x.
TH1F * Getwf(int)
Helper that returns the waveform histogram for a certain waveform number number.
int maxNWF
Maximum possible number of waveforms in data: number of active channels x number of events.
void Print_GetTimingCFD(float=100, float=140, int=0, int=-999, string="S", bool=true)
Plot results of GetTimingCFD()
int PlotChannelAverages_cnt
Index for multiple executions of the same plotting function.
void ChargeCorrelation(float, float, float, float, float, float, int, int, int, bool=false)
Plot correlation of integrals/amplitudes between two channels.
float tCutEndg
End of time window for baseline correction when using ReadRun::Using_BaselineCorrection_in_file_loop.
TH1F * ChargeSpectrum(int, float, float, float=0, float=300, float=-50, float=600, int=750)
Histogram of the "charge" spectrum for one channel.
float DAQ_factor
DAQ conversion factor for wavecatcher output to mV.
void CorrectBaselineMin(vector< float >, double=0, int=2)
Baseline correction using minimum sum ( mean) in range for correction.
int nevents
Number of triggered events in data.
float IndexToTime(int)
Convert the bin number of the waveform to the time of the left bin edge
vector< bool > skip_event
Stores the event numbers which should be skipped in the analysis.
double * gety(int)
Get array of y values for a certain waveform.
vector< int > active_channels
Stores the numbers of the active channels.
vector< vector< float > > rundata
Stores waveforms.
void CorrectBaseline_function(vector< float > &, float, float, int)
Helper function called by CorrectBaseline()
vector< bool > PolarityMap(bool, int)
Channel map of polarity changes during reading. For parameters see ReadFile().
int TimeToIndex(float)
Convert time to the bin number of the waveform.
float GetPeakIntegral(TH1F *, float, float, float, float, int=0)
Calculate the integral around a peak with several options explained in GetIntWindow().
int GetEventIndex(unsigned int)
Returns index of a certain event number (if data files are read in parallel threads)
TH1F * His_GetTimingCFD(int, float, float, int=-999)
Plot results of GetTimingCFD()
void FilterAll(double=.3, double=.9, double=.2)
Filter all waveforms.
void PlotWFHeatmaps(float=-20, float=200, int=880, string="", float=0, EColorPalette=kRainBow)
Plot stacks of all non-skipped waveforms for all active channels.
int GetChannelIndex(int)
Match channel number (wavecatcher input channel) to channel index.
float tCutg
Start of time window for baseline correction when using ReadRun::Using_BaselineCorrection_in_file_loo...
float * ChargeList(int, float=20, float=80, float=0, float=300, bool=1)
Returns array with the individual "charge"/amplitude for all events of one channel.
TFile * root_out
Stores results of analysis.
int nChannelsWC
Wavecatcher hardware max. number of channels (reduce/increase if not using the 64 channel crate)
void PrintBaselineCorrectionResults(float=-5, float=5, int=200)
Print histogram of the baseline correction values for all channels.
virtual void ReadFile(string, bool=false, int=9, string="out.root", bool=false, long long=-1)
Routine to read files created by the wavecatcher.
virtual int GetWaveformIndex(int, int)
Returns index of a certain event number (if data files are read in parallel threads)
vector< float > PrintChargeSpectrum_pars
Starting values of the fit parameters for PrintChargeSpectrum()
bool Using_BaselineCorrection_in_file_loop
Set true for baseline correction during data reading. Needs to be called before ReadFile().
void PrintWFProjection(float=0, float=320, float=-50, float=50, int=200)
Plots waveform projection histograms of all channels.
void PlotChannelAverages(bool=false)
Plot averages only of the good, corrected waveforms for each channel.
bool PlotChannel(int)
Check if a channel index should be plotted according to ReadRun::plot_active_channels.
array< int, 3 > GetIntWindow(TH1F *, float, float, float, float, int=0)
Determine indices for integration window for peaks.
int PrintWFProjection_cnt
Index for multiple executions of the same plotting function.
int start_read_at_channel
Do analysis only for limited range of channels to reduce memory usage.
vector< int > maxSumBin
Stores bin numbers where the sum of waveforms have their maximum.
void PrintChargeSpectrumWF(float, float, float=0, float=300, int=1, float=0., float=0., float=0., float=0.)
Plot waveforms of all channels for a given event number and add the determined integration windows to...
vector< float > mean_integral
Stores the mean integral/lightyield from PrintChargeSpectrum() for all channels.
bool discard_original_eventnr
Can be used to discard the original event numbering of the data.
virtual int Nevents_good(int=0)
Number of good events that are not skipped.
vector< vector< float > > amplValuessum
Collects sums of all waveforms for each channel.
vector< int > switch_polarity_for_channels
Stores the channel number where the polarity should be inverted. Example use to switch polarity for c...
TH1F * WFProjectionChannel(int, int=0, int=1024, float=-50, float=50, int=200)
Waveform projections for one channel.
TGraph2D * MaxDist(int, float=0, float=300)
Finds maximum amplitude for a given channel in time window ["from", "to"] and creates 3d map of wavef...
int nwf
Total number of waveforms read from data: number of active channels x number of events.
TH1F * TimeDist(int, float=0, float=300, float=0, float=300, int=100, int=0, float=.3)
Time distribution of maximum, CFD, or 10% - 90% rise time in a certain time window.
void GetTimingCFD(float=.3, float=100, float=140, double=0., bool=true, int=2, bool=false, bool=false)
Determine the timing of the maximum peak with constant fraction discrimination.
void IntegralFilter(vector< float >, vector< bool >, float, float, float=50, float=250, bool=false, bool=false)
Skip events with threshold on integral.
TH1F * BaselineCorrectionResults(int, int, float=-5, float=5, int=200)
Histograms of the contents of baseline_correction_result.
void PrintDCR(float=15, float=85, float=0, float=300, double=3)
Calculate (SiPM) dark count rate.
int PrintChargeSpectrum_cnt
Index for multiple executions of the same plotting function.
void SmoothAll(double, int)
Smoothing all waveforms which are not skipped (for testing, careful when using for analysis!...
void CorrectBaseline(float, float=-999)
Baseline correction constant window.
int binNumber
Number of bins (usually 1024 samples per waveform).
vector< vector< float > > PrintChargeSpectrum_cal
Calibration values to normalize charge spectrum to number of photoelectrons Chennels must be ordered ...
vector< vector< float > > timing_results
Matrix to store timing of peaks from GetTimingCFD()
void SaveChargeLists(float=20, float=80, float=0, float=300, bool=1)
Saves TGraphs to root file with the individual "charge"/amplitude for all events and all channels.
vector< vector< float > > baseline_correction_result
Stores baseline correction results for CorrectBaseline() and related functions.
TH1F * His_GetTimingCFD_diff(vector< int >, vector< int >, float, float, int=-999)
Plot timing difference between the mean timings of two channel ranges.
void SkipEventsTimeDiffCut(int, int, double, double, bool=false)
Skip events where the time difference between two channels is outside of specified range.
int PlotWFHeatmaps_cnt
Index for multiple executions of the same plotting function.
void SkipEventsPerChannel(vector< float >, float=0, float=0, bool=false)
Skip events above/below individual thresholds per channel.
string data_path
Path to data.
void PrintMaxDist(float=0, float=300)
Finds maximum amplitude for a given channel in time window ["from", "to"] and creates 3d map of wavef...
void PrintChargeSpectrum(float, float, float=0, float=300, float=-50, float=600, int=750, float=0., float=0., int=99, int=0, bool=false)
Plots the "charge" spectrums of all channels.
virtual void checkData(bool isBaselineCorrection=false) const
Primitive check to see if data has been loaded.
int end_read_at_channel
See ReadRun::start_read_at_channel.
int nchannels
Number of active channels in data.
void UnskipAll()
Sets skip_event flag to false for all events, removing any previous cuts.
float SP
Sampling: ns per bin of data, sampling rate 3.2 GS/s -> 0.3125 ns.
vector< TFitResultPtr > timing_fit_results
Stores the fit results of Print_GetTimingCFD() for all channels.
T * getx(double shift=0.)
Get array of x axis (time of the bin centers) for standard wavecatcher settings.
virtual int GetCurrentEvent(int)
Get the current event index for a certain waveform index.
virtual ~ReadRun()
Destructor.
virtual bool SkipEvent(int, int=-1)
Check if event should be skipped.
void PlotChannelSums(bool=false, bool=false, double=0., double=0., int=2)
Plot sums of all raw waveforms for each channel.
int CheckBoundsX(int)
Check if index exists in time of waveforms.
void PrintTimeDist(float=0, float=300, float=0, float=300, int=100, int=0, float=.3)
Time distribution of maximum, CFD, or 10% - 90% rise time in a certain time window.
vector< int > plot_active_channels
Stores the numbers of the active channels which should be plotted.
void ShiftAllToAverageCF()
This function shifts all waveforms to the average signal starting times for each channel.
vector< TFitResultPtr > fit_results
Stores the fit results of PrintChargeSpectrum() for all channels and all function calls in ascending ...
virtual int GetCurrentChannel(int)
Get the current channel index for a certain waveform index.
void PrintSkippedEvents()
Prints a list of all skipped events into the terminal for diagnostics.
int LastBinFileToRead
Number of last .bin file to be read in.
void CorrectBaselineMinSlopeRMS(vector< float >, double=0, int=2)
Baseline correction method searching for non-monotonic, rather constant regions.
TH2F * WFHeatmapChannel(int, float=-20, float=200, int=880)
2D histogram of all non-skipped waveforms for one channel
void Print_GetTimingCFD_diff(vector< int >, vector< int >, float=100, float=140, int=0, int=-999, float=-999, float=-999, string="RS", bool=true)
Plot timing difference between the mean timings of two channel ranges.
vector< unsigned int > eventnr_storage
Events will be stored here in the order they have been read.