MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
SplineMonolith.h
Go to the documentation of this file.
1#pragma once
2
4
5//KS: Joy of forward declaration https://gieseanw.wordpress.com/2018/02/25/the-joys-of-forward-declarations-results-from-the-real-world/
6class SMonolithGPU;
7
12class SMonolith : public SplineBase {
13 public:
18 SMonolith(std::vector<std::vector<TResponseFunction_red*> > &MasterSpline,
19 const std::vector<RespFuncType> &SplineType,
20 const bool SaveFlatTree = false);
23 SMonolith(const std::string& FileName);
25 virtual ~SMonolith();
26
28 void Evaluate() override;
29
31 inline std::string GetName() const {return "SplineMonolith";};
32
35
39 inline const float* retPointer(const int event) {return &cpu_total_weights[event];}
40
43 inline void setSplinePointers(std::vector< const double* > spline_ParsPointers) {
44 splineParsPointer = spline_ParsPointers;
45 for (M3::int_t i = 0; i < nParams; ++i) SplineInfoArray[i].splineParsPointer = spline_ParsPointers[i];
46 };
47
52
53 private:
55 inline void Initialise();
66 inline void ScanMasterSpline(std::vector<std::vector<TResponseFunction_red*> > & MasterSpline,
67 unsigned int &nEvents,
68 short int &MaxPoints,
69 short int &numParams,
70 int &nSplines,
71 unsigned int &NSplinesValid,
72 unsigned int &numKnots,
73 unsigned int &nTF1Valid,
74 unsigned int &nTF1_coeff,
75 const std::vector<RespFuncType> &SplineType);
78 inline void PrepareForGPU(std::vector<std::vector<TResponseFunction_red*> > &MasterSpline, const std::vector<RespFuncType> &SplineType);
80 inline void MoveToGPU();
81
83 inline void PrintInitialsiation();
84
91 inline void getSplineCoeff_SepMany(TSpline3_red* &spl, int &nPoints, float *&xArray, float *&manyArray);
92
94 inline void CalcSplineWeights() override;
96 inline void ModifyWeights() override;
98 inline void ModifyWeights_GPU();
99
101 inline void PrepareSplineFile();
104 inline void LoadSplineFile(std::string FileName);
105
107 std::vector< const double* > splineParsPointer;
108
110 unsigned int NEvents;
112 short int _max_knots;
114 std::vector<int> index_spline_cpu;
116 std::vector<int> index_TF1_cpu;
117
119 unsigned int NSplines_valid;
121 unsigned int NTF1_valid;
122
125
127 unsigned int nKnots;
129 unsigned int nTF1coeff;
130
135
137 std::vector<unsigned int> cpu_nParamPerEvent;
138
140 std::vector<unsigned int> cpu_nParamPerEvent_tf1;
141
144
147
149 std::vector<float> cpu_coeff_TF1_many;
150
152 std::vector<short int> cpu_nPoints_arr;
153
155 std::vector<short int> cpu_paramNo_TF1_arr;
156
159};
Class responsible for calculating spline weight on GPU.
Even-by-event class calculating response for spline parameters. It is possible to use GPU acceleratio...
SMonolithGPU * gpu_spline_handler
KS: Store info about Spline monolith, this allow to obtain better step time. As all necessary informa...
std::vector< unsigned int > cpu_nParamPerEvent
KS: CPU map keeping track how many parameters applies to each event, we keep two numbers here {number...
unsigned int NTF1_valid
Number of valid TF1.
std::vector< short int > cpu_nPoints_arr
CPU arrays to hold number of points.
void Evaluate() override
CW: This Eval should be used when using two separate x,{y,a,b,c,d} arrays to store the weights; proba...
unsigned int nTF1coeff
Sum of all coefficients over all TF1.
void setSplinePointers(std::vector< const double * > spline_ParsPointers)
KS: Set pointers to spline params.
std::vector< float > cpu_coeff_TF1_many
CPU arrays to hold TF1 coefficients.
void Initialise()
KS: Set everything to null etc.
SplineMonoStruct * cpu_spline_handler
KS: Store info about Spline monolith, this allow to obtain better step time. As all necessary informa...
float * cpu_weights_tf1_var
CPU arrays to hold weight for each TF1.
const float * retPointer(const int event)
KS: Get pointer to total weight to make fit faster wrooom!
void PrintInitialsiation()
KS: Print info about how much knots etc has been initialised.
std::vector< unsigned int > cpu_nParamPerEvent_tf1
KS: CPU map keeping track how many parameters applies to each event, we keep two numbers here {number...
void LoadSplineFile(std::string FileName)
KS: Load preprocessed spline file.
void SynchroniseMemTransfer()
KS: After calculations are done on GPU we copy memory to CPU. This operation is asynchronous meaning ...
std::vector< int > index_TF1_cpu
holds the index for good TF1; don't do unsigned since starts with negative value!
void ModifyWeights_GPU()
Conversion from valid splines to all.
void PrepareSplineFile()
KS: Prepare spline file that can be used for fast loading.
float * cpu_weights
The returned gpu weights, read by the GPU.
float * cpu_weights_spline_var
CPU arrays to hold weight for each spline.
std::vector< int > index_spline_cpu
holds the index for good splines; don't do unsigned since starts with negative value!
void PrepareForGPU(std::vector< std::vector< TResponseFunction_red * > > &MasterSpline, const std::vector< RespFuncType > &SplineType)
CW: Prepare the TSpline3_red objects for the GPU.
void CalcSplineWeights() override
CPU based code which eval weight for each spline.
void MoveToGPU()
CW: The shared initialiser from constructors of TResponseFunction_red.
bool SaveSplineFile
Flag telling whether we are saving spline monolith into handy root file.
float * cpu_total_weights
KS: This holds the total CPU weights that gets read in samplePDFND.
unsigned int NSplines_valid
Number of valid splines.
virtual ~SMonolith()
Destructor for SMonolith class.
void ModifyWeights() override
Calc total event weight.
unsigned int NEvents
Number of events.
void getSplineCoeff_SepMany(TSpline3_red *&spl, int &nPoints, float *&xArray, float *&manyArray)
CW: This loads up coefficients into two arrays: one x array and one yabcd array.
std::string GetName() const
Get class name.
unsigned int NSplines_total_large
Number of total splines if each event had every parameter's spline.
void ScanMasterSpline(std::vector< std::vector< TResponseFunction_red * > > &MasterSpline, unsigned int &nEvents, short int &MaxPoints, short int &numParams, int &nSplines, unsigned int &NSplinesValid, unsigned int &numKnots, unsigned int &nTF1Valid, unsigned int &nTF1_coeff, const std::vector< RespFuncType > &SplineType)
CW: Function to scan through the MasterSpline of TSpline3.
short int _max_knots
Max knots for production.
std::vector< short int > cpu_paramNo_TF1_arr
CW: CPU array with the number of points per spline (not per spline point!)
std::vector< const double * > splineParsPointer
This holds pointer to parameter position which we later copy paste it to GPU.
unsigned int nKnots
Sum of all knots over all splines.
Base class for calculating weight from spline.
Definition: SplineBase.h:25
short int nParams
Number of parameters that have splines.
Definition: SplineBase.h:64
std::vector< FastSplineInfo > SplineInfoArray
Definition: SplineBase.h:57
CW: Reduced TSpline3 class.
int int_t
Definition: Core.h:29
KS: Struct storing information for spline monolith.
Definition: SplineCommon.h:30