MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
SampleHandlerFD.h
Go to the documentation of this file.
1#pragma once
2
3//MaCh3 includes
5
7
11
12#include "THStack.h"
13#include "TLegend.h"
14
19{
20public:
21 //######################################### Functions #########################################
24 SampleHandlerFD(std::string ConfigFileName, ParameterHandlerGeneric* xsec_cov,
25 const std::shared_ptr<OscillationHandler>& OscillatorObj_ = nullptr);
27 virtual ~SampleHandlerFD();
28
31 int GetNDim() const { return nDimensions; }
33 std::string GetSampleName(int iSample = 0) const override;
35 std::string GetTitle() const override {return SampleTitle;}
36
38 std::string GetXBinVarName() {return XVarStr;}
40 std::string GetYBinVarName() {return YVarStr;}
41
42 void PrintIntegral(TString OutputName="/dev/null", int WeightStyle=0, TString OutputCSVName="/dev/null");
43
44 //===============================================================================
45 // DB Reweighting and Likelihood functions
46
47 //ETA - abstract these to SampleHandlerFDBase
48 //DB Require these four functions to allow conversion from TH1(2)D to array for multi-threaded GetLikelihood
49 void AddData(TH1D* Data);
50 void AddData(TH2D* Data);
51 void AddData(std::vector<double> &data);
52 void AddData(std::vector< std::vector <double> > &data);
53
55 double GetLikelihood() override;
56 //===============================================================================
57
60 TH1* GetMCHist(const int Dimension);
61
64 TH1* GetW2Hist(const int Dimension);
65
68 TH1* GetDataHist(const int Dimension);
69
70 void Reweight() override;
71 M3::float_t GetEventWeight(const int iEntry) const;
72
76
77 void ReadSampleConfig();
78
80 int GetNOscChannels() override {return static_cast<int>(OscChannels.size());}
81
83 std::string GetFlavourName(const int iChannel) {
84 if (iChannel < 0 || iChannel > static_cast<int>(OscChannels.size())) {
85 MACH3LOG_ERROR("Invalid Channel Requested: {}", iChannel);
86 throw MaCh3Exception(__FILE__ , __LINE__);
87 }
88 return OscChannels[iChannel].flavourName;
89 }
90
92 TH1* Get1DVarHist(const std::string& ProjectionVar, const std::vector< KinematicCut >& EventSelectionVec = std::vector< KinematicCut >(),
93 int WeightStyle=0, TAxis* Axis=nullptr, const std::vector< KinematicCut >& SubEventSelectionVec = std::vector< KinematicCut >());
94 TH2* Get2DVarHist(const std::string& ProjectionVarX, const std::string& ProjectionVarY,
95 const std::vector< KinematicCut >& EventSelectionVec = std::vector< KinematicCut >(),
96 int WeightStyle=0, TAxis* AxisX=nullptr, TAxis* AxisY=nullptr,
97 const std::vector< KinematicCut >& SubEventSelectionVec = std::vector< KinematicCut >());
98
100 void Fill1DSubEventHist(TH1D* _h1DVar, const std::string& ProjectionVar, const std::vector< KinematicCut >& SubEventSelectionVec = std::vector< KinematicCut >(),
101 int WeightStyle=0);
102 void Fill2DSubEventHist(TH2D* _h2DVar, const std::string& ProjectionVarX, const std::string& ProjectionVarY,
103 const std::vector< KinematicCut >& SubEventSelectionVec = std::vector< KinematicCut >(), int WeightStyle=0);
104
106 TH1* Get1DVarHistByModeAndChannel(const std::string& ProjectionVar_Str, int kModeToFill=-1, int kChannelToFill=-1, int WeightStyle=0, TAxis* Axis=nullptr);
108 TH2* Get2DVarHistByModeAndChannel(const std::string& ProjectionVar_StrX, const std::string& ProjectionVar_StrY, int kModeToFill=-1, int kChannelToFill=-1, int WeightStyle=0, TAxis* AxisX=nullptr, TAxis* AxisY=nullptr);
109
111 TH1 *GetModeHist1D(int s, int m, int style = 0) {
112 return Get1DVarHistByModeAndChannel(XVarStr,m,s,style);
113 }
115 TH2 *GetModeHist2D(int s, int m, int style = 0) {
117 }
118
120 std::vector<TH1*> ReturnHistsBySelection1D(std::string KinematicProjection, int Selection1,int Selection2=-1, int WeightStyle=0, TAxis* Axis=0);
122 std::vector<TH2*> ReturnHistsBySelection2D(std::string KinematicProjectionX, std::string KinematicProjectionY, int Selection1, int Selection2=-1, int WeightStyle=0, TAxis* XAxis=0, TAxis* YAxis=0);
124 THStack* ReturnStackedHistBySelection1D(std::string KinematicProjection, int Selection1, int Selection2=-1, int WeightStyle=0, TAxis* Axis=0);
127
130 int ReturnKinematicParameterFromString(const std::string& KinematicStr) const;
133 std::string ReturnStringFromKinematicParameter(const int KinematicVariable) const;
134
136 void SaveAdditionalInfo(TDirectory* Dir) override;
137
138 // === JM declare the same functions for kinematic vectors ===
139 int ReturnKinematicVectorFromString(const std::string& KinematicStr) const;
140 std::string ReturnStringFromKinematicVector(const int KinematicVariable) const;
141 // ===========================================================
143 bool IsSubEventVarString(const std::string& VarStr);
144
145 protected:
147 virtual void SetupWeightPointers() = 0;
148
150 void SetupKinematicMap();
151
154 virtual void SetupSplines() = 0;
155
156 //DB Require all objects to have a function which reads in the MC
158 virtual void Init() = 0;
159
161 virtual int SetupExperimentMC() = 0;
162
164 virtual void SetupFDMC() = 0;
165
167 void Initialise();
168
170 std::unique_ptr<BinnedSplineHandler> SplineHandler;
171
173 std::shared_ptr<OscillationHandler> Oscillator;
174 //===============================================================================
177 void FillSplineBins();
178
179 //Functions which find the nominal bin and bin edges
182
183 //DB Overrided functions of base class which calculate erec bin and boundaries for reweighting speedup in beam samples
184 //ETA - this can be done using the core info stored in the new fdmc_struct
188 void Set1DBinning(size_t nbins, double* boundaries);
194 void Set2DBinning(size_t nbins1, double* boundaries1, size_t nbins2, double* boundaries2);
197 void Set1DBinning(std::vector<double> &XVec) {Set1DBinning(XVec.size()-1, XVec.data());};
201 void Set2DBinning(std::vector<double> &XVec, std::vector<double> &YVec) {Set2DBinning(XVec.size()-1, XVec.data(), YVec.size()-1, YVec.data());};
203 void SetupSampleBinning();
205 void SetupReweightArrays(const size_t numberXBins, const size_t numberYBins);
206
208 std::string XVarStr, YVarStr;
209 //===============================================================================
210
211 // ----- Functional Parameters -----
213 virtual void SetupFunctionalParameters();
215 void RegisterIndividualFunctionalParameter(const std::string& fpName, int fpEnum, FuncParFuncType fpFunc);
219 virtual void PrepFunctionalParameters(){};
221 virtual void ApplyShifts(int iEvent);
222
224 bool IsEventSelected(const int iEvent);
226 bool IsSubEventSelected(const std::vector<KinematicCut> &SubEventCuts, const int iEvent, unsigned const int iSubEvent, size_t nsubevents);
228 virtual void resetShifts(int iEvent) {(void)iEvent;};
230 std::vector<FunctionalParameter> funcParsVec;
233 std::unordered_map<std::string, int> funcParsNamesMap;
238 std::vector<FunctionalParameter *> funcParsMap;
241 std::unordered_map<int, FuncParFuncType> funcParsFuncMap;
243 std::vector<std::vector<int>> funcParsGrid;
245 std::vector<std::string> funcParsNamesVec = {};
246
248 void CalcNormsBins(std::vector<NormParameter>& norm_parameters, std::vector< std::vector< int > >& xsec_norms_bins);
252 M3::float_t CalcWeightNorm(const FarDetectorCoreInfo* MCEvent) const;
253
259 virtual void CalcWeightFunc(int iEvent){return; (void)iEvent;};
260
262 virtual double ReturnKinematicParameter(std::string KinematicParamter, int iEvent) = 0;
263 virtual double ReturnKinematicParameter(int KinematicVariable, int iEvent) = 0;
264
265 // === JM declare the same functions for kinematic vectors ===
266 virtual std::vector<double> ReturnKinematicVector(std::string KinematicParameter, int iEvent) {return {}; (void)KinematicParameter; (void)iEvent;};
267 virtual std::vector<double> ReturnKinematicVector(int KinematicVariable, int iEvent) {return {}; (void)KinematicVariable; (void)iEvent;};
268 // ===========================================================
269
271 std::vector<double> ReturnKinematicParameterBinning(const std::string& KinematicParameter);
272 virtual const double* GetPointerToKinematicParameter(std::string KinematicParamter, int iEvent) = 0;
273 virtual const double* GetPointerToKinematicParameter(double KinematicVariable, int iEvent) = 0;
274
276 const double* GetPointerToOscChannel(const int iEvent) const;
278 void SetupNormParameters();
279
280 //===============================================================================
281 //DB Functions required for reweighting functions
282 //DB Replace previous implementation with reading bin contents from SampleHandlerFD_array
284 void Fill1DHist();
286 void Fill2DHist();
287
289#ifdef MULTITHREAD
291 void FillArray_MP();
292#endif
294 void FillArray();
295
297 inline void ResetHistograms();
298
299 //===============================================================================
300 //DB Variables required for GetLikelihood
303
310 //===============================================================================
311
312 //===============================================================================
313 //MC variables
314 std::vector<FarDetectorCoreInfo> MCSamples;
315 std::vector<OscChannelInfo> OscChannels;
316 //===============================================================================
317
318 //===============================================================================
319 //DB Covariance Objects
320 //ETA - All experiments will need an xsec, det and osc cov
321 //these should be added to SampleHandlerBase to be honest
323
324 //===============================================================================
325
329 std::string SampleName;
330
332 std::string SampleTitle;
333
334 //===========================================================================
335 //DB Vectors to store which kinematic cuts we apply
336 //like in XsecNorms but for events in sample. Read in from sample yaml file
337 //What gets used in IsEventSelected, which gets set equal to user input plus
338 //all the vectors in StoreSelection
339
342 std::vector< KinematicCut > StoredSelection;
345 std::vector< KinematicCut > Selection;
346 //===========================================================================
347
349 const std::unordered_map<std::string, int>* KinematicParameters;
351 const std::unordered_map<int, std::string>* ReversedKinematicParameters;
352
353 // === JM mapping between string and kinematic vector enum ===
354 const std::unordered_map<std::string, int>* KinematicVectors;
355 const std::unordered_map<int, std::string>* ReversedKinematicVectors;
356 // ===========================================================
357
359 std::unique_ptr<manager> SampleManager;
364
366 std::vector<std::string> mc_files;
368 std::vector<std::string> spline_files;
369
370 std::unordered_map<std::string, double> _modeNomWeightMap;
371
372 //===============================================================================
374 TLegend* THStackLeg = nullptr;
375 //===============================================================================
376
381
382 TH1D *dathist;
384
385 // binned PDFs
386 TH1D* _hPDF1D;
387 TH2D* _hPDF2D;
388
390 NuPDG GetInitPDGFromFileName(const std::string& FileName) const {return FileToInitPDGMap.at(FileName);}
392 NuPDG GetFinalPDGFromFileName(const std::string& FileName) const {return FileToFinalPDGMap.at(FileName);}
393
394private:
395 std::unordered_map<std::string, NuPDG> FileToInitPDGMap;
396 std::unordered_map<std::string, NuPDG> FileToFinalPDGMap;
397
398};
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:25
std::function< void(const double *, std::size_t)> FuncParFuncType
NuPDG
Enum to track the incoming neutrino species.
Custom exception class for MaCh3 errors.
Class responsible for handling of systematic error parameters with different types defined in the con...
Class responsible for handling implementation of samples used in analysis, reweighting and returning ...
Class responsible for handling implementation of samples used in analysis, reweighting and returning ...
std::vector< std::vector< int > > funcParsGrid
HH - a grid of vectors of enums for each sample and event.
void CalcNormsBins(std::vector< NormParameter > &norm_parameters, std::vector< std::vector< int > > &xsec_norms_bins)
Check whether a normalisation systematic affects an event or not.
std::vector< FunctionalParameter > funcParsVec
HH - a vector that stores all the FuncPars struct.
virtual void resetShifts(int iEvent)
HH - reset the shifted values to the original values.
std::vector< FunctionalParameter * > funcParsMap
HH - a map that relates the funcpar enum to pointer of FuncPars struct HH - Changed to a vector of po...
virtual void Init()=0
Initialise any variables that your experiment specific SampleHandler needs.
std::string ReturnStringFromKinematicVector(const int KinematicVariable) const
std::vector< double > ReturnKinematicParameterBinning(const std::string &KinematicParameter)
Return the binning used to draw a kinematic parameter.
void AddData(TH1D *Data)
std::unique_ptr< BinnedSplineHandler > SplineHandler
Contains all your binned splines and handles the setup and the returning of weights from spline evalu...
virtual void CalcWeightFunc(int iEvent)
Calculate weights for function parameters.
SampleBinningInfo Binning
KS: This stores binning information, in future could be come vector to store binning for every used s...
virtual void SetupSplines()=0
initialise your splineXX object and then use InitialiseSplineObject to conviently setup everything up
bool IsEventSelected(const int iEvent)
DB Function which determines if an event is selected, where Selection double looks like {{ND280Kinema...
ParameterHandlerGeneric * ParHandler
std::string SampleTitle
the name of this sample e.g."muon-like"
bool FirstTimeW2
KS:Super hacky to update W2 or not.
std::unordered_map< std::string, NuPDG > FileToInitPDGMap
NuPDG GetInitPDGFromFileName(const std::string &FileName) const
Retrieve the initial neutrino PDG code associated with a given input file name.
std::vector< KinematicCut > Selection
a way to store selection cuts which you may push back in the get1DVar functions most of the time this...
virtual int SetupExperimentMC()=0
Experiment specific setup, returns the number of events which were loaded.
const std::unordered_map< int, std::string > * ReversedKinematicVectors
double ** SampleHandlerFD_array_w2
KS Array used for MC stat.
std::vector< KinematicCut > StoredSelection
What gets pulled from config options, these are constant after loading in this is of length 3: 0th in...
virtual std::vector< double > ReturnKinematicVector(std::string KinematicParameter, int iEvent)
std::vector< std::string > spline_files
names of spline files associated associated with this object
double GetLikelihood() override
DB Multi-threaded GetLikelihood.
std::unordered_map< std::string, NuPDG > FileToFinalPDGMap
std::unordered_map< int, FuncParFuncType > funcParsFuncMap
HH - a map that relates the funcpar enum to pointer of the actual function.
const double * GetPointerToOscChannel(const int iEvent) const
Get pointer to oscillation channel associated with given event. Osc channel is const.
std::vector< OscChannelInfo > OscChannels
std::string XVarStr
the strings associated with the variables used for the binning e.g. "RecoNeutrinoEnergy"
const std::unordered_map< std::string, int > * KinematicVectors
void SetupNormParameters()
Setup the norm parameters by assigning each event with bin.
void SetupKinematicMap()
Ensure Kinematic Map is setup and make sure it is initialised correctly.
void FillSplineBins()
Finds the binned spline that an event should apply to and stored them in a a vector for easy evaluati...
virtual void PrepFunctionalParameters()
Update the functional parameter values to the latest proposed values. Needs to be called before every...
void Fill2DSubEventHist(TH2D *_h2DVar, const std::string &ProjectionVarX, const std::string &ProjectionVarY, const std::vector< KinematicCut > &SubEventSelectionVec=std::vector< KinematicCut >(), int WeightStyle=0)
void Initialise()
Function which does a lot of the lifting regarding the workflow in creating different MC objects.
virtual ~SampleHandlerFD()
destructor
void FindNominalBinAndEdges1D()
std::vector< std::string > funcParsNamesVec
HH - a vector of string names for each functional parameter.
virtual void SetupFDMC()=0
Function which translates experiment struct into core struct.
std::string YVarStr
void SaveAdditionalInfo(TDirectory *Dir) override
Store additional info in a chan.
std::unique_ptr< manager > SampleManager
The manager object used to read the sample yaml file.
M3::float_t GetEventWeight(const int iEntry) const
void Fill1DHist()
Fill a 1D histogram with the event-level information used in the fit.
void InitialiseNuOscillatorObjects()
including Dan's magic NuOscillator
void Fill2DHist()
Fill a 2D histogram with the event-level information used in the fit.
void RegisterIndividualFunctionalParameter(const std::string &fpName, int fpEnum, FuncParFuncType fpFunc)
HH - a helper function for RegisterFunctionalParameter.
M3::float_t CalcWeightNorm(const FarDetectorCoreInfo *MCEvent) const
Calculate the norm weight for a given event.
bool UpdateW2
KS:Super hacky to update W2 or not.
std::shared_ptr< OscillationHandler > Oscillator
Contains all your binned splines and handles the setup and the returning of weights from spline evalu...
double ** SampleHandlerFD_array
DB Array to be filled after reweighting.
bool IsSubEventSelected(const std::vector< KinematicCut > &SubEventCuts, const int iEvent, unsigned const int iSubEvent, size_t nsubevents)
JM Function which determines if a subevent is selected.
virtual const double * GetPointerToKinematicParameter(double KinematicVariable, int iEvent)=0
const std::unordered_map< std::string, int > * KinematicParameters
Mapping between string and kinematic enum.
void Set2DBinning(size_t nbins1, double *boundaries1, size_t nbins2, double *boundaries2)
set the binning for 2D sample used for the likelihood calculation
void InitialiseSingleFDMCObject()
function to create the member of the FarDetectorInfo struct so they are the appropriate size.
virtual std::vector< double > ReturnKinematicVector(int KinematicVariable, int iEvent)
TH2 * Get2DVarHist(const std::string &ProjectionVarX, const std::string &ProjectionVarY, const std::vector< KinematicCut > &EventSelectionVec=std::vector< KinematicCut >(), int WeightStyle=0, TAxis *AxisX=nullptr, TAxis *AxisY=nullptr, const std::vector< KinematicCut > &SubEventSelectionVec=std::vector< KinematicCut >())
virtual void RegisterFunctionalParameters()=0
HH - a experiment-specific function where the maps to actual functions are set up.
M3::float_t CalcWeightSpline(const FarDetectorCoreInfo *MCEvent) const
Calculate the spline weight for a given event.
NuPDG GetFinalPDGFromFileName(const std::string &FileName) const
Retrieve the final neutrino PDG code associated with a given input file name.
std::vector< std::string > mc_files
names of mc files associated associated with this object
void ResetHistograms()
Helper function to reset histograms.
virtual double ReturnKinematicParameter(std::string KinematicParamter, int iEvent)=0
Return the value of an assocaited kinematic parameter for an event.
bool IsSubEventVarString(const std::string &VarStr)
JM Check if a kinematic parameter string corresponds to a subevent-level variable.
virtual void ApplyShifts(int iEvent)
ETA - generic function applying shifts.
std::unordered_map< std::string, double > _modeNomWeightMap
const std::unordered_map< int, std::string > * ReversedKinematicParameters
Mapping between kinematic enum and string.
int ReturnKinematicVectorFromString(const std::string &KinematicStr) const
double ** SampleHandlerFD_data
DB Array to be filled in AddData.
virtual void SetupWeightPointers()=0
DB Function to determine which weights apply to which types of samples pure virtual!...
void SetupSampleBinning()
wrapper to call set binning functions based on sample config info
virtual double ReturnKinematicParameter(int KinematicVariable, int iEvent)=0
void Reweight() override
void SetupNuOscillatorPointers()
void Set2DBinning(std::vector< double > &XVec, std::vector< double > &YVec)
set the binning for 2D sample used for the likelihood calculation
virtual void SetupFunctionalParameters()
ETA - a function to setup and pass values to functional parameters where you need to pass a value to ...
void Set1DBinning(std::vector< double > &XVec)
set the binning for 1D sample used for the likelihood calculation
void FillArray()
DB Nice new multi-threaded function which calculates the event weights and fills the relevant bins of...
int nDimensions
Keep track of the dimensions of the sample binning.
TLegend * THStackLeg
DB Miscellaneous Variables.
void PrintIntegral(TString OutputName="/dev/null", int WeightStyle=0, TString OutputCSVName="/dev/null")
std::vector< FarDetectorCoreInfo > MCSamples
void SetupReweightArrays(const size_t numberXBins, const size_t numberYBins)
Initialise data, MC and W2 histograms.
std::unordered_map< std::string, int > funcParsNamesMap
HH - a map that relates the name of the functional parameter to funcpar enum.
void Set1DBinning(size_t nbins, double *boundaries)
sets the binning used for the likelihood calculation, used for both data and MC
std::string SampleName
A unique ID for each sample based on powers of two for quick binary operator comparisons.
virtual const double * GetPointerToKinematicParameter(std::string KinematicParamter, int iEvent)=0
std::string GetYBinVarName()
TH1 * GetModeHist1D(int s, int m, int style=0)
TLegend * ReturnStackHistLegend()
TH1 * Get1DVarHist(const std::string &ProjectionVar, const std::vector< KinematicCut > &EventSelectionVec=std::vector< KinematicCut >(), int WeightStyle=0, TAxis *Axis=nullptr, const std::vector< KinematicCut > &SubEventSelectionVec=std::vector< KinematicCut >())
std::string ReturnStringFromKinematicParameter(const int KinematicVariable) const
ETA function to generically convert a kinematic type from xsec cov to a string.
std::string GetSampleName(int iSample=0) const override
void Fill1DSubEventHist(TH1D *_h1DVar, const std::string &ProjectionVar, const std::vector< KinematicCut > &SubEventSelectionVec=std::vector< KinematicCut >(), int WeightStyle=0)
int GetNDim() const
DB Function to differentiate 1D or 2D binning.
TH2 * GetModeHist2D(int s, int m, int style=0)
std::vector< TH2 * > ReturnHistsBySelection2D(std::string KinematicProjectionX, std::string KinematicProjectionY, int Selection1, int Selection2=-1, int WeightStyle=0, TAxis *XAxis=0, TAxis *YAxis=0)
std::string GetTitle() const override
std::string GetXBinVarName()
std::string GetFlavourName(const int iChannel)
int GetNOscChannels() override
TH2 * Get2DVarHistByModeAndChannel(const std::string &ProjectionVar_StrX, const std::string &ProjectionVar_StrY, int kModeToFill=-1, int kChannelToFill=-1, int WeightStyle=0, TAxis *AxisX=nullptr, TAxis *AxisY=nullptr)
TH1 * GetW2Hist(const int Dimension)
Get W2 histogram.
TH1 * Get1DVarHistByModeAndChannel(const std::string &ProjectionVar_Str, int kModeToFill=-1, int kChannelToFill=-1, int WeightStyle=0, TAxis *Axis=nullptr)
TH1 * GetMCHist(const int Dimension)
Get MC histogram.
std::vector< TH1 * > ReturnHistsBySelection1D(std::string KinematicProjection, int Selection1, int Selection2=-1, int WeightStyle=0, TAxis *Axis=0)
int ReturnKinematicParameterFromString(const std::string &KinematicStr) const
ETA function to generically convert a string from xsec cov to a kinematic type.
TH1 * GetDataHist(const int Dimension)
Get Data histogram.
THStack * ReturnStackedHistBySelection1D(std::string KinematicProjection, int Selection1, int Selection2=-1, int WeightStyle=0, TAxis *Axis=0)
double float_t
Definition: Core.h:28
static constexpr const int _BAD_INT_
Default value used for int initialisation.
Definition: Core.h:45
constructors are same for all three so put in here
KS: Small struct storying info about used binning.