MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
SplineBase Class Referenceabstract

Base class for calculating weight from spline. More...

#include <Splines/SplineBase.h>

Inheritance diagram for SplineBase:
[legend]
Collaboration diagram for SplineBase:
[legend]

Public Member Functions

 SplineBase ()
 Constructor.
 
virtual ~SplineBase ()
 Destructor.
 
virtual void Evaluate ()=0
 CW: This Eval should be used when using two separate x,{y,a,b,c,d} arrays to store the weights; probably the best one here! Same thing but pass parameter spline segments instead of variations.
 
virtual std::string GetName () const
 Get class name.
 
short int GetNParams () const
 Get number of spline parameters.
 

Protected Member Functions

void FindSplineSegment ()
 CW:Code used in step by step reweighting, Find Spline Segment for each param.
 
virtual void CalcSplineWeights ()=0
 CPU based code which eval weight for each spline.
 
virtual void ModifyWeights ()=0
 Calc total event weight.
 
void getTF1Coeff (TF1_red *&spl, int &nPoints, float *&coeffs)
 CW: Gets the polynomial coefficients for TF1.
 

Protected Attributes

std::vector< FastSplineInfoSplineInfoArray
 
short int * SplineSegments
 
float * ParamValues
 Store parameter values they are not in FastSplineInfo as in case of GPU we need to copy paste it to GPU.
 
short int nParams
 Number of parameters that have splines.
 

Detailed Description

Base class for calculating weight from spline.

Definition at line 25 of file SplineBase.h.

Constructor & Destructor Documentation

◆ SplineBase()

SplineBase::SplineBase ( )

Constructor.

Definition at line 7 of file SplineBase.cpp.

7 {
8// *****************************************
9 nParams = 0;
10 SplineSegments = nullptr;
11 ParamValues = nullptr;
12}
short int nParams
Number of parameters that have splines.
Definition: SplineBase.h:64
short int * SplineSegments
Definition: SplineBase.h:60
float * ParamValues
Store parameter values they are not in FastSplineInfo as in case of GPU we need to copy paste it to G...
Definition: SplineBase.h:62

◆ ~SplineBase()

SplineBase::~SplineBase ( )
virtual

Destructor.

Definition at line 16 of file SplineBase.cpp.

16 {
17// *****************************************
18}

Member Function Documentation

◆ CalcSplineWeights()

virtual void SplineBase::CalcSplineWeights ( )
protectedpure virtual

CPU based code which eval weight for each spline.

Implemented in PySplineBase, BinnedSplineHandler, and SMonolith.

◆ Evaluate()

virtual void SplineBase::Evaluate ( )
pure virtual

CW: This Eval should be used when using two separate x,{y,a,b,c,d} arrays to store the weights; probably the best one here! Same thing but pass parameter spline segments instead of variations.

Implemented in BinnedSplineHandler, PySplineBase, and SMonolith.

◆ FindSplineSegment()

void SplineBase::FindSplineSegment ( )
protected

CW:Code used in step by step reweighting, Find Spline Segment for each param.

Definition at line 24 of file SplineBase.cpp.

24 {
25// *************************
26 // Loop over the splines
27 //KS: Tried multithreading here with 48 splines and it is faster with one thread, maybe in future multithreading will be worth revisiting
28 for (M3::int_t i = 0; i < nParams; ++i)
29 {
30 const M3::int_t nPoints = SplineInfoArray[i].nPts;
31 const std::vector<M3::float_t>& xArray = SplineInfoArray[i].xPts;
32
33 // Get the variation for this reconfigure for the ith parameter
34 const float xvar = float(*SplineInfoArray[i].splineParsPointer);
35 ParamValues[i] = xvar;
36
37 // EM: if we have a parameter that has no response for any event (i.e. all splines have just one knot), then skip it and avoid a seg fault here
38 // In principle, such parameters shouldn't really be included in the first place, but with new det syst splines this
39 // could happen if say you were just running on one FHC run, then all RHC parameters would be flat and the code below would break.
40 if(xArray.size() == 0) continue;
41
42 // The segment we're interested in (klow in ROOT code)
43 M3::int_t segment = 0;
44 M3::int_t kHigh = nPoints-1;
45 //KS: We expect new segment is very close to previous
46 const M3::int_t PreviousSegment = SplineInfoArray[i].CurrSegment;
47
48 // If the variation is below the lowest saved spline point
49 if (xvar <= xArray[0]) {
50 segment = 0;
51 // If the variation is above the highest saved spline point
52 } else if (xvar >= xArray[nPoints-1]) {
53 //CW: Yes, the -2 is indeed correct, see TSpline.cxx:814 and //see: https://savannah.cern.ch/bugs/?71651
54 segment = kHigh;
55 //KS: It is quite probable the new segment is same as in previous step so try to avoid binary search, first we have to check if it is in bounds to avoid seg fault
56 } else if( xArray[PreviousSegment+1] > xvar && xvar >= xArray[PreviousSegment] ) {
57 segment = PreviousSegment;
58 // If the variation is between the maximum and minimum, perform a binary search
59 } else {
60 // The top point we've got
61 M3::int_t kHalf = 0;
62 // While there is still a difference in the points (we haven't yet found the segment)
63 // This is a binary search, incrementing segment and decrementing kHalf until we've found the segment
64 while (kHigh - segment > 1) {
65 // Increment the half-step
66 kHalf = M3::int_t((segment + kHigh)/2);
67 // If our variation is above the kHalf, set the segment to kHalf
68 if (xvar > xArray[kHalf]) {
69 segment = kHalf;
70 // Else move kHigh down
71 } else {
72 kHigh = kHalf;
73 }
74 } // End the while: we've now done our binary search
75 } // End the else: we've now found our point
76
77 if (segment >= nPoints-1 && nPoints > 1) segment = nPoints-2;
78
79 //CW: This way we avoid doing 1.2M+ binary searches on the GPU
80 // and literally just multiply lots of numbers together on the GPU without any algorithm
81 // Update the values and which segment it belongs to
82 SplineInfoArray[i].CurrSegment = segment;
83 SplineSegments[i] = short(SplineInfoArray[i].CurrSegment);
84
85 #ifdef DEBUG
86 if (SplineInfoArray[i].xPts[segment] > xvar && segment != 0) {
87 MACH3LOG_ERROR("Found a segment which is _ABOVE_ the variation!");
88 MACH3LOG_ERROR("IT SHOULD ALWAYS BE BELOW! (except when segment 0)");
89 MACH3LOG_ERROR("Spline: {}", i);
90 MACH3LOG_ERROR("Found segment = {}", segment);
91 MACH3LOG_ERROR("Doing variation = {}", xvar);
92 MACH3LOG_ERROR("x in spline = {}", SplineInfoArray[i].xPts[segment]);
93 for (M3::int_t j = 0; j < SplineInfoArray[j].nPts; ++j) {
94 MACH3LOG_ERROR(" {} = {}", j, SplineInfoArray[i].xPts[j]);
95 }
96 throw MaCh3Exception(__FILE__ , __LINE__ );
97 }
98 #endif
99 } //end loop over params
100}
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:25
Custom exception class for MaCh3 errors.
std::vector< FastSplineInfo > SplineInfoArray
Definition: SplineBase.h:57
int int_t
Definition: Core.h:29

◆ GetName()

virtual std::string SplineBase::GetName ( ) const
inlinevirtual

Get class name.

Reimplemented in SMonolith, and PySplineBase.

Definition at line 36 of file SplineBase.h.

36{return "SplineBase";};

◆ GetNParams()

short int SplineBase::GetNParams ( ) const
inline

Get number of spline parameters.

Definition at line 39 of file SplineBase.h.

39{return nParams;};

◆ getTF1Coeff()

void SplineBase::getTF1Coeff ( TF1_red *&  spl,
int &  nPoints,
float *&  coeffs 
)
protected

CW: Gets the polynomial coefficients for TF1.

Parameters
splpointer to TF1_red that will be checked
nPointsnumber of knots
coeffsArray holding coefficients for each knot

Definition at line 106 of file SplineBase.cpp.

106 {
107// *****************************************
108 // Initialise all arrays to 1.0
109 for (int i = 0; i < _nTF1Coeff_; ++i) {
110 coeffs[i] = 0.0;
111 }
112
113 // Get number of points in spline
114 nPoints = spl->GetSize();
115
116 if(nPoints > _nTF1Coeff_)
117 {
118 MACH3LOG_ERROR("Too big number of points for TF1");
119 throw MaCh3Exception(__FILE__ , __LINE__ );
120 }
121 // TSpline3 can only take doubles, not floats
122 // But our GPU is slow with doubles, so need to cast to float
123 for (int i = 0; i < nPoints; i++) {
124 coeffs[i] = float(spl->GetParameter(M3::int_t(i)));
125 }
126 // The structure is now coeffs = {a,b,c,d,e}
127}
#define _nTF1Coeff_
KS: For TF1 we store at most 5 coefficients, we could make it more flexible but for now define it her...
Definition: SplineCommon.h:15
double GetParameter(M3::int_t Parameter)
Get a parameter value.
int GetSize()
Get the size.

◆ ModifyWeights()

virtual void SplineBase::ModifyWeights ( )
protectedpure virtual

Calc total event weight.

Implemented in PySplineBase, BinnedSplineHandler, and SMonolith.

Member Data Documentation

◆ nParams

short int SplineBase::nParams
protected

Number of parameters that have splines.

Definition at line 64 of file SplineBase.h.

◆ ParamValues

float* SplineBase::ParamValues
protected

Store parameter values they are not in FastSplineInfo as in case of GPU we need to copy paste it to GPU.

Definition at line 62 of file SplineBase.h.

◆ SplineInfoArray

std::vector<FastSplineInfo> SplineBase::SplineInfoArray
protected

Array of FastSplineInfo structs: keeps information on each xsec spline for fast evaluation Method identical to TSpline3::Eval(double) but faster because less operations

Definition at line 57 of file SplineBase.h.

◆ SplineSegments

short int* SplineBase::SplineSegments
protected

Store currently found segment they are not in FastSplineInfo as in case of GPU we need to copy paste it to GPU

Warning
this is being used sometimes by GPU, therefore must be raw pointer!

Definition at line 60 of file SplineBase.h.


The documentation for this class was generated from the following files: