MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
SplineMonolith.cpp
Go to the documentation of this file.
1#include "SplineMonolith.h"
2
3#ifdef MaCh3_CUDA
5#endif
6
7#pragma GCC diagnostic ignored "-Wuseless-cast"
8#pragma GCC diagnostic ignored "-Wfloat-conversion"
9
10// *****************************************
11//Set everything to NULL or 0
13// *****************************************
14#ifdef MaCh3_CUDA
15 MACH3LOG_INFO("Using GPU version event by event monolith");
16 gpu_spline_handler = nullptr;
17#endif
18
20
21 nKnots = 0;
22 nTF1coeff = 0;
23 NEvents = 0;
24 _max_knots = 0;
25
27 NTF1_valid = 0;
29
30 cpu_weights_spline_var = nullptr;
31 cpu_weights = nullptr;
32 cpu_weights_tf1_var = nullptr;
33
34 cpu_total_weights = nullptr;
35}
36
37// *****************************************
38SMonolith::SMonolith(std::vector<std::vector<TResponseFunction_red*> > &MasterSpline, const std::vector<RespFuncType> &SplineType, const bool SaveFlatTree)
39: SplineBase() {
40// *****************************************
41
42 //KS: If true it will save spline monolith into huge ROOT file
43 SaveSplineFile = SaveFlatTree;
44 Initialise();
45 MACH3LOG_INFO("-- GPUING WITH arrays and master spline containing TResponseFunction_red");
46
47 // Convert the TSpline3 pointers to the reduced form and call the reduced constructor
48 PrepareForGPU(MasterSpline, SplineType);
49}
50
51// *****************************************
52// The shared initialiser from constructors of TSpline3 and TSpline3_red
53void SMonolith::PrepareForGPU(std::vector<std::vector<TResponseFunction_red*> > &MasterSpline, const std::vector<RespFuncType> &SplineType) {
54// *****************************************
55
56 // Scan for the max number of knots, the number of events (number of splines), and number of parameters
57 int maxnSplines = 0;
58 ScanMasterSpline(MasterSpline,
59 NEvents,
61 nParams,
62 maxnSplines,
64 nKnots,
67 SplineType);
68
69 MACH3LOG_INFO("Found {} events", NEvents);
70 MACH3LOG_INFO("Found {} knots at max", _max_knots);
71 MACH3LOG_INFO("Found {} parameters", nParams);
72 MACH3LOG_INFO("Found {} maximum number of splines in an event", maxnSplines);
73 MACH3LOG_INFO("Found total {} knots in all splines", nKnots);
74 MACH3LOG_INFO("Number of splines = {}", NSplines_valid);
75 MACH3LOG_INFO("Found total {} coeffs in all TF1", nTF1coeff);
76 MACH3LOG_INFO("Number of TF1 = {}", NTF1_valid);
77
78 // Can pass the spline segments to the GPU instead of the values
79 // Make these here and only refill them for each loop, avoiding unnecessary new/delete on each reconfigure
80 //KS: Since we are going to copy it each step use fancy CUDA memory allocation
81 #ifdef MaCh3_CUDA
84 #else
85 SplineSegments = new short int[nParams];
86 ParamValues = new float[nParams];
87 #endif
88
89 for (M3::int_t j = 0; j < nParams; j++)
90 {
91 SplineSegments[j] = 0;
92 ParamValues[j] = -999;
93 }
94
95 // Number of objects we have in total if each event has *EVERY* spline. Needed for some arrays
97
98 unsigned int event_size_max = _max_knots * nParams;
99 // Declare the {x}, {y,b,c,d} arrays for all possible splines which the event has
100 // We'll filter off the flat and "disabled" (e.g. CCQE event should not have MARES spline) ones in the next for loop, but need to declare these beasts here
101
102 // Declare the {y,b,c,d} for each knot
103 // float because GPU precision (could change to double, but will incur significant speed reduction on GPU unless you're very rich!)
104 cpu_spline_handler->coeff_many.resize(nKnots*_nCoeff_); // *4 because we store y,b,c,d parameters in this array
105 //KS: For x coeff we assume that for given dial (MAQE) spacing is identical, here we are sloppy and assume each dial has the same number of knots, not a big problem
106 cpu_spline_handler->coeff_x.resize(event_size_max);
107
108 // Set all the big arrays to -999 to keep us safe...
109 for (unsigned int j = 0; j < event_size_max; j++) {
110 cpu_spline_handler->coeff_x[j] = -999;
111 }
112
113 //CW: With TF1 we only save the coefficients and the order of the polynomial
114 // Makes most sense to have one large monolithic array, but then it becomes impossible to tell apart a coefficient from a "number of points". So have two arrays: one of coefficients and one of number of points
115 // Let's first assume all are of _max_knots size
116 // Now declare the arrays for each point in the valid splines which the event actually has (i.e. include the splines that the event undergoes)
117 // Also make array with the number of points per spline (not per spline point!)
118 // float because GPU precision (could change to double, but will incur significant speed reduction on GPU unless you're very rich!)
120 cpu_coeff_TF1_many.resize(nTF1coeff); // *5 because this array holds a,b,c,d,e parameters
121
122 #ifdef Weight_On_SplineBySpline_Basis
123 // This holds the index of each spline
126
127 #ifdef MULTITHREAD
128 #pragma omp parallel for
129 #endif
130 for (unsigned int j = 0; j < NSplines_total_large; j++) {
131 index_spline_cpu[j] = -1;
132 index_TF1_cpu[j] = -1;
133 }
134 // This holds the total CPU weights that gets read in samplePDFND
136 #else
137 //KS: Map keeping track how many parameters applies to each event, we keep two numbers here {number of splines per event, index where splines start for a given event}
138 cpu_nParamPerEvent.resize(2*NEvents);
140 #ifdef MULTITHREAD
141 #pragma omp parallel for
142 #endif
143 for (unsigned int j = 0; j < 2*NEvents; j++) {
144 cpu_nParamPerEvent[j] = -1;
146 }
147 #endif
148
149 // Make array with the number of points per spline (not per spline point!)
151 //KS: And array which tells where each spline stars in a big monolith array, sort of knot map
154
155 // Temporary arrays to hold the coefficients for each spline
156 // We get one x, one y, one b,... for each point, so only need to be _max_knots big
157 //KS: Some params has less splines but this is all right main array will get proper number while this temp will be deleted
158 float *x_tmp = new float[_max_knots]();
159 float *many_tmp = new float[_max_knots*_nCoeff_]();
160 float *temp_coeffs = new float[_nTF1Coeff_]();
161
162 // Count the number of events
163 unsigned int KnotCounter = 0;
164 unsigned int TF1PointsCounter = 0;
165 unsigned int NSplinesCounter = 0;
166 unsigned int TF1sCounter = 0;
167 int ParamCounter = 0;
168 int ParamCounterGlobal = 0;
169 int ParamCounter_TF1 = 0;
170 int ParamCounterGlobalTF1 = 0;
171 // Loop over events and extract the spline coefficients
172 for(unsigned int EventCounter = 0; EventCounter < MasterSpline.size(); ++EventCounter) {
173 // Structure of MasterSpline is std::vector<std::vector<TSpline3*>>
174 // A conventional iterator to count which parameter a given spline should be applied to
175 for(unsigned int ParamNumber = 0; ParamNumber < MasterSpline[EventCounter].size(); ++ParamNumber) {
176
177 // If NULL we don't have this spline for the event, so move to next spline
178 if (MasterSpline[EventCounter][ParamNumber] == NULL) continue;
179
180 if(SplineType[ParamNumber] == kTSpline3_red)
181 {
182 //KS: how much knots each spline has
183 int nPoints_tmp = 0;
184 // Get a pointer to the current spline for this event
185 TResponseFunction_red* TespFunc = MasterSpline[EventCounter][ParamNumber];
186 TSpline3_red* CurrSpline = static_cast<TSpline3_red*>(TespFunc);
187
188 // If the number of knots are greater than 2 the spline is not a dummy and we should extract coefficients to load onto the GPU
189 getSplineCoeff_SepMany(CurrSpline, nPoints_tmp, x_tmp, many_tmp);
190
191 //KS: One knot means flat spline so ignore
192 if (nPoints_tmp == 1) continue;
193 for (int j = 0; j < _max_knots; ++j) {
194 cpu_spline_handler->coeff_x[ParamNumber*_max_knots + j] = x_tmp[j];
195 }
196 //KS: Contrary to X coeff we keep for other coeff only filled knots, there is no much gain for doing so for x coeff
197 for (int j = 0; j < nPoints_tmp; ++j) {
198 for (int k = 0; k < _nCoeff_; k++) {
199 cpu_spline_handler->coeff_many[KnotCounter*_nCoeff_ + j*_nCoeff_ + k] = many_tmp[j*_nCoeff_+k];
200 }
201 }
202 // Set the parameter number for this spline
203 cpu_spline_handler->paramNo_arr[NSplinesCounter] = short(ParamNumber);
204 //KS: Fill map when each spline starts
205 cpu_spline_handler->nKnots_arr[NSplinesCounter] = KnotCounter;
206 KnotCounter += nPoints_tmp;
207
208 #ifdef Weight_On_SplineBySpline_Basis
209 // Set the index of the spline so we can tell apart from flat splines
210 index_spline_cpu[EventCounter*nParams + ParamNumber] = NSplinesCounter;
211 #else
212 ++ParamCounter;
213 #endif
214 // Increment the counter for the number of good splines we have
215 ++NSplinesCounter;
216 }
217 else if (SplineType[ParamNumber] == kTF1_red)
218 {
219 // Don't actually use this ever -- we give each spline the maximum number of points found in all splines
220 int nPoints_tmp = 0;
221 // Get a pointer to the current spline for this event
222 TF1_red* CurrSpline = dynamic_cast<TF1_red*>(MasterSpline[EventCounter][ParamNumber]);
223
224 // If the number of knots are greater than 2 the spline is not a dummy and we should extract coefficients to load onto the GPU
225 getTF1Coeff(CurrSpline, nPoints_tmp, temp_coeffs);
226 for (int j = 0; j < _nTF1Coeff_; ++j) {
227 cpu_coeff_TF1_many[TF1PointsCounter+j] = temp_coeffs[j];
228 }
229 // Save the number of points for this spline
230 cpu_nPoints_arr[TF1sCounter] = short(nPoints_tmp);
231
232 TF1PointsCounter += nPoints_tmp;
233 // Set the parameter number for this spline
234 cpu_paramNo_TF1_arr[TF1sCounter] = short(ParamNumber);
235 #ifdef Weight_On_SplineBySpline_Basis
236 // Set the index of the spline so we can tell apart from flat splines
237 index_TF1_cpu[EventCounter*nParams + ParamNumber] = TF1sCounter;
238 #else
239 ++ParamCounter_TF1;
240 #endif
241 // Increment the counter for the number of good splines we have
242 ++TF1sCounter;
243 }
244 //KS: Don't delete in debug
245 #ifndef DEBUG
246 delete MasterSpline[EventCounter][ParamNumber];
247 MasterSpline[EventCounter][ParamNumber] = NULL;
248 #endif
249 } // End the loop over the parameters in the MasterSpline
250 #ifndef Weight_On_SplineBySpline_Basis
251 cpu_nParamPerEvent[2*EventCounter] = ParamCounter;
252 cpu_nParamPerEvent[2*EventCounter+1] = ParamCounterGlobal;
253 ParamCounterGlobal += ParamCounter;
254
255 cpu_nParamPerEvent_tf1[2*EventCounter] = ParamCounter_TF1;
256 cpu_nParamPerEvent_tf1[2*EventCounter+1] = ParamCounterGlobalTF1;
257 ParamCounterGlobalTF1 += ParamCounter_TF1;
258
259 ParamCounter = 0;
260 ParamCounter_TF1 = 0;
261 #endif
262 } // End the loop over the number of events
263 delete[] many_tmp;
264 delete[] x_tmp;
265 delete[] temp_coeffs;
266
267 int BadXCounter = 0;
268 for (unsigned int j = 0; j < event_size_max; j++) {
269 if (cpu_spline_handler->coeff_x[j] == -999) BadXCounter++;
270 // Perform checks that all entries have been modified from initial values
271 if (cpu_spline_handler->coeff_x[j] == -999 && BadXCounter < 5) {
272 MACH3LOG_WARN("***** BAD X !! *****");
273 MACH3LOG_WARN("Indicates some parameter doesn't have a single spline");
274 MACH3LOG_WARN("j = {}", j);
275 //throw MaCh3Exception(__FILE__ , __LINE__ );
276 }
277 if(BadXCounter == 5) MACH3LOG_WARN("There is more unutilised knots although I will stop spamming");
278 }
279
280 MACH3LOG_WARN("Found in total {} BAD X", BadXCounter);
281 #ifdef Weight_On_SplineBySpline_Basis
282 // Make the array that holds all the returned weights from the GPU to pass to the CPU
284 cpu_weights_tf1_var = new float[NTF1_valid]();
285 #else
286 //KS: This is tricky as this variable use both by CPU and GPU, however if use CUDA we use cudaMallocHost
287 #ifndef MaCh3_CUDA
288 cpu_total_weights = new float[NEvents]();
290 cpu_weights_tf1_var = new float[NTF1_valid]();
291 #endif
292 #endif
293
294 // Print some info; could probably make this to a separate function
296
298
299 MoveToGPU();
300}
301
302// *****************************************
303// The shared initialiser from constructors of TSpline3 and TSpline3_red
305// *****************************************
306 #ifdef MaCh3_CUDA
307 unsigned int event_size_max = _max_knots * nParams;
308 MACH3LOG_INFO("Total size = {:.2f} MB memory on CPU to move to GPU",
309 (double(sizeof(float) * nKnots * _nCoeff_) + double(sizeof(float) * event_size_max) / 1.E6 +
310 double(sizeof(short int) * NSplines_valid)) / 1.E6);
311 MACH3LOG_INFO("Total TF1 size = {:.2f} MB memory on CPU to move to GPU",
312 double(sizeof(float) * NTF1_valid * _nTF1Coeff_) / 1.E6);
313 MACH3LOG_INFO("GPU weight array (GPU->CPU every step) = {:.2f} MB", static_cast<double>(sizeof(float)) * (NSplines_valid + NTF1_valid) / 1.0e6);
314 #ifndef Weight_On_SplineBySpline_Basis
315 MACH3LOG_INFO("Since you are running Total event weight mode then GPU weight array (GPU->CPU every step) = {:.2f} MB",
316 double(sizeof(float) * NEvents) / 1.E6);
317 #endif
318 MACH3LOG_INFO("Parameter value array (CPU->GPU every step) = {:.4f} MB", double(sizeof(float) * nParams) / 1.E6);
319 //CW: With the new set-up we have: 1 coefficient array of size coeff_array_size, all same size
320 // 1 coefficient array of size coeff_array_size*4, holding y,b,c,d in order (y11,b11,c11,d11; y12,b12,c12,d12;...) where ynm is n = spline number, m = spline point. Should really make array so that order is (y11,b11,c11,d11; y21,b21,c21,d21;...) because it will optimise cache hits I think; try this if you have time
321 // return gpu_weights
322
324
325 // The gpu_XY arrays don't actually need initialising, since they are only placeholders for what we'll move onto the GPU. As long as we cudaMalloc the size of the arrays correctly there shouldn't be any problems
326 // Can probably make this a bit prettier but will do for now
327 // Could be a lot smaller of a function...
329 #ifndef Weight_On_SplineBySpline_Basis
331 NEvents,
332 #endif
333 nKnots, // How many entries in coefficient array (*4 for the "many" array)
334 NSplines_valid, // What's the number of splines we have (also number of entries in gpu_nPoints_arr)
336 event_size_max //Knots times event number of unique splines
337 );
338
339 // Move number of splines and spline size to constant GPU memory; every thread does not need a copy...
340 // The implementation lives in splines/gpuSplineUtils.cu
341 // The GPU splines don't actually need declaring but is good for demonstration, kind of
342 // fixed by passing const reference
345
346 // TFI related now
349 #ifndef Weight_On_SplineBySpline_Basis
350 NEvents,
353 #endif
354 nParams,
357 nKnots,
358 NTF1_valid);
359
360 // Delete all the coefficient arrays from the CPU once they are on the GPU
363 #ifndef Weight_On_SplineBySpline_Basis
366 #endif
367 delete cpu_spline_handler;
368 cpu_spline_handler = nullptr;
369 MACH3LOG_INFO("Good GPU loading");
370 #endif
371}
372
373// Need to specify template functions in header
374// *****************************************
375// Scan the master spline to get the maximum number of knots in any of the TSpline3*
376void SMonolith::ScanMasterSpline(std::vector<std::vector<TResponseFunction_red*> > & MasterSpline,
377 unsigned int &nEvents,
378 short int &MaxPoints,
379 short int &numParams,
380 int &nSplines,
381 unsigned int &NSplinesValid,
382 unsigned int &numKnots,
383 unsigned int &nTF1Valid,
384 unsigned int &nTF1_coeff,
385 const std::vector<RespFuncType> &SplineType) {
386// *****************************************
387 // Need to extract: the total number of events
388 // number of parameters
389 // maximum number of knots
390 MaxPoints = 0;
391 nEvents = 0;
392 numParams = 0;
393 nSplines = 0;
394 numKnots = 0;
395 NSplinesValid = 0;
396 nTF1Valid = 0;
397 nTF1_coeff = 0;
398
399 // Check the number of events
400 nEvents = int(MasterSpline.size());
401
402 // Maximum number of splines one event can have (scan through and find this number)
403 int nMaxSplines_PerEvent = 0;
404
405 //KS: We later check that each event has the same number of splines so this is fine
406 numParams = short(MasterSpline[0].size());
407 // Initialise
408 SplineInfoArray.resize(numParams);
409
410 // Loop over each parameter
411 for(unsigned int EventCounter = 0; EventCounter < MasterSpline.size(); ++EventCounter) {
412 // Check that each event has each spline saved
413 if (numParams > 0) {
414 int TempSize = int(MasterSpline[EventCounter].size());
415 if (TempSize != numParams) {
416 MACH3LOG_ERROR("Found {} parameters for event {}", TempSize, EventCounter);
417 MACH3LOG_ERROR("but was expecting {} since that's what I found for the previous event", numParams);
418 MACH3LOG_ERROR("Somehow this event has a different number of spline parameters... Please study further!");
419 throw MaCh3Exception(__FILE__ , __LINE__ );
420 }
421 }
422 numParams = short(MasterSpline[EventCounter].size());
423
424 int nSplines_SingleEvent = 0;
425 int nPoints = 0;
426 // Loop over each pointer
427 for(unsigned int ParamNumber = 0; ParamNumber < MasterSpline[EventCounter].size(); ++ParamNumber) {
428
429 if (MasterSpline[EventCounter][ParamNumber]) {
430 if(SplineType[ParamNumber] == kTSpline3_red)
431 {
432 TResponseFunction_red* TespFunc = MasterSpline[EventCounter][ParamNumber];
433 TSpline3_red* CurrSpline = dynamic_cast<TSpline3_red*>(TespFunc);
434 if(CurrSpline){
435 nPoints = CurrSpline->GetNp();
436 }
437
438 if (nPoints > MaxPoints) {
439 MaxPoints = static_cast<short int>(nPoints);
440 }
441 numKnots += nPoints;
442 nSplines_SingleEvent++;
443
444 // Fill the SplineInfoArray entries with information on each splinified parameter
445 if (SplineInfoArray[ParamNumber].xPts.size() == 0)
446 {
447 // Fill the number of points
448 SplineInfoArray[ParamNumber].nPts = CurrSpline->GetNp();
449
450 // Fill the x points
451 SplineInfoArray[ParamNumber].xPts.resize(SplineInfoArray[ParamNumber].nPts);
452 for (M3::int_t k = 0; k < SplineInfoArray[ParamNumber].nPts; ++k)
453 {
454 M3::float_t xtemp = M3::float_t(-999.99);
455 M3::float_t ytemp = M3::float_t(-999.99);
456 CurrSpline->GetKnot(k, xtemp, ytemp);
457 SplineInfoArray[ParamNumber].xPts[k] = xtemp;
458 }
459 }
460 NSplinesValid++;
461 }
462 else if (SplineType[ParamNumber] == kTF1_red)
463 {
464 TResponseFunction_red* TespFunc = MasterSpline[EventCounter][ParamNumber];
465 TF1_red* CurrSpline = dynamic_cast<TF1_red*>(TespFunc);
466 nPoints = CurrSpline->GetSize();
467 nTF1_coeff += nPoints;
468 nTF1Valid++;
469 }
470 } else {
471 // If NULL we don't have this spline for the event, so move to next spline
472 continue;
473 }
474 }
475 if (nSplines_SingleEvent > nMaxSplines_PerEvent) nMaxSplines_PerEvent = nSplines_SingleEvent;
476 }
477 nSplines = nMaxSplines_PerEvent;
478
479 int Counter = 0;
480 //KS: Sanity check that everything was set correctly
481 for (M3::int_t i = 0; i < numParams; ++i)
482 {
483 // KS: We don't find segment for TF1, so ignore this
484 if (SplineType[i] == kTF1_red) continue;
485
486 const M3::int_t nPoints = SplineInfoArray[i].nPts;
487 const std::vector<M3::float_t>& xArray = SplineInfoArray[i].xPts;
488 if (nPoints == -999 || xArray.size() == 0) {
489 Counter++;
490 if(Counter < 5) {
491 MACH3LOG_WARN("SplineInfoArray[{}] isn't set yet", i);
492 }
493 continue;
494 //throw MaCh3Exception(__FILE__ , __LINE__ );
495 }
496 }
497 MACH3LOG_WARN("In total SplineInfoArray for {} hasn't been initialised", Counter);
498}
499
500// *****************************************
501// Load SplineFile
502SMonolith::SMonolith(const std::string& FileName)
503 : SplineBase() {
504// *****************************************
505 Initialise();
506 MACH3LOG_INFO("-- GPUING WITH {X} and {Y,B,C,D} arrays and master spline containing TSpline3_red");
507 // Convert the TSpline3 pointers to the reduced form and call the reduced constructor
508 LoadSplineFile(FileName);
509}
510
511// *****************************************
512// Load SplineMonolith from ROOT file
513void SMonolith::LoadSplineFile(std::string FileName) {
514// *****************************************
515 #ifdef Weight_On_SplineBySpline_Basis
516 MACH3LOG_ERROR("Trying to load Monolith from file using weight by weight base, this is not supported right now, sorry");
517 throw MaCh3Exception(__FILE__ , __LINE__ );
518 #endif
519
520 if (std::getenv("MACH3") != nullptr) {
521 FileName.insert(0, std::string(std::getenv("MACH3"))+"/");
522 }
523
524 auto SplineFile = std::make_unique<TFile>(FileName.c_str(), "OPEN");
525 TTree *Settings = SplineFile->Get<TTree>("Settings");
526 TTree *Monolith_TF1 = SplineFile->Get<TTree>("Monolith_TF1");
527 TTree *EventInfo = SplineFile->Get<TTree>("EventInfo");
528 TTree *FastSplineInfoTree = SplineFile->Get<TTree>("FastSplineInfoTree");
529 TTree *SplineTree = SplineFile->Get<TTree>("SplineTree");
530
531 unsigned int NEvents_temp;
532 short int nParams_temp;
533 int _max_knots_temp;
534 unsigned int nKnots_temp;
535 unsigned int NSplines_valid_temp;
536 unsigned int nTF1Valid_temp;
537 unsigned int nTF1coeff_temp;
538
539 Settings->SetBranchAddress("NEvents", &NEvents_temp);
540 Settings->SetBranchAddress("nParams", &nParams_temp);
541 Settings->SetBranchAddress("_max_knots", &_max_knots_temp);
542 Settings->SetBranchAddress("nKnots", &nKnots_temp);
543 Settings->SetBranchAddress("NSplines_valid", &NSplines_valid_temp);
544 Settings->SetBranchAddress("NTF1_valid", &nTF1Valid_temp);
545 Settings->SetBranchAddress("nTF1coeff", &nTF1coeff_temp);
546
547 Settings->GetEntry(0);
548
549 NEvents = NEvents_temp;
550 nParams = nParams_temp;
551 _max_knots = static_cast<short int>(_max_knots_temp);
552 nKnots = nKnots_temp;
553 NSplines_valid = NSplines_valid_temp;
554 NTF1_valid = nTF1Valid_temp;
555 nTF1coeff = nTF1coeff_temp;
556
557 //KS: Since we are going to copy it each step use fancy CUDA memory allocation
558#ifdef MaCh3_CUDA
561#else
562 SplineSegments = new short int[nParams]();
563 ParamValues = new float[nParams]();
564#endif
565
566 cpu_nParamPerEvent.resize(2*NEvents);
569
570 //KS: This is tricky as this variable use both by CPU and GPU, however if use CUDA we use cudaMallocHost
571#ifndef MaCh3_CUDA
572 cpu_total_weights = new float[NEvents]();
574 cpu_weights_tf1_var = new float[NTF1_valid]();
575#endif
576
577 SplineTree->SetBranchAddress("SplineObject", &cpu_spline_handler);
578 SplineTree->GetEntry(0);
579
580 float coeff_tf1 = 0.;
581 Monolith_TF1->SetBranchAddress("cpu_coeff_TF1_many", &coeff_tf1);
582 for(unsigned int i = 0; i < nTF1coeff; i++)
583 {
584 Monolith_TF1->GetEntry(i);
585 cpu_coeff_TF1_many[i] = coeff_tf1;
586 }
587
588 unsigned int nParamPerEvent = 0;
589 unsigned int nParamPerEvent_tf1 = 0;
590
591 EventInfo->SetBranchAddress("cpu_nParamPerEvent", &nParamPerEvent);
592 EventInfo->SetBranchAddress("cpu_nParamPerEvent_tf1", &nParamPerEvent_tf1);
593 for(unsigned int i = 0; i < 2*NEvents; i++)
594 {
595 EventInfo->GetEntry(i);
596 cpu_nParamPerEvent[i] = nParamPerEvent;
597 cpu_nParamPerEvent_tf1[i] = nParamPerEvent_tf1;
598 }
599
600 M3::int_t nPoints = 0;
601 float xtemp[20];
602 FastSplineInfoTree->SetBranchAddress("nPts", &nPoints);
603 FastSplineInfoTree->SetBranchAddress("xPts", &xtemp);
604
605 SplineInfoArray.resize(nParams);
606 for (M3::int_t i = 0; i < nParams; ++i) {
607 FastSplineInfoTree->GetEntry(i);
608
609 // Fill the number of points
610 SplineInfoArray[i].nPts = nPoints;
611 if(nPoints == -999) continue;
612 SplineInfoArray[i].xPts.resize(SplineInfoArray[i].nPts);
613 for (M3::int_t k = 0; k < SplineInfoArray[i].nPts; ++k)
614 {
615 SplineInfoArray[i].xPts[k] = xtemp[k];
616 }
617 }
618 SplineFile->Close();
619
620 // Print some info; could probably make this to a separate function
622
623 MoveToGPU();
624}
625
626// *****************************************
627// Save SplineMonolith into ROOT file
629// *****************************************
630 std::string FileName = "SplineFile.root";
631 if (std::getenv("MACH3") != nullptr) {
632 FileName.insert(0, std::string(std::getenv("MACH3"))+"/");
633 }
634
635 auto SplineFile = std::make_unique<TFile>(FileName.c_str(), "recreate");
636 TTree *Settings = new TTree("Settings", "Settings");
637 TTree *Monolith_TF1 = new TTree("Monolith_TF1", "Monolith_TF1");
638 TTree *XKnots = new TTree("XKnots", "XKnots");
639 TTree *EventInfo = new TTree("EventInfo", "EventInfo");
640 TTree *FastSplineInfoTree = new TTree("FastSplineInfoTree", "FastSplineInfoTree");
641
642 unsigned int NEvents_temp = NEvents;
643 short int nParams_temp = nParams;
644 int _max_knots_temp = _max_knots;
645 unsigned int nKnots_temp = nKnots;
646 unsigned int NSplines_valid_temp = NSplines_valid;
647 unsigned int nTF1Valid_temp = NTF1_valid;
648 unsigned int nTF1coeff_temp = nTF1coeff;
649
650 Settings->Branch("NEvents", &NEvents_temp, "NEvents/i");
651 Settings->Branch("nParams", &nParams_temp, "nParams/S");
652 Settings->Branch("_max_knots", &_max_knots_temp, "_max_knots/I");
653 Settings->Branch("nKnots", &nKnots_temp, "nKnots/i");
654 Settings->Branch("NSplines_valid", &NSplines_valid_temp, "NSplines_valid/i");
655 Settings->Branch("NTF1_valid", &nTF1Valid_temp, "NTF1_valid/i");
656 Settings->Branch("nTF1coeff", &nTF1coeff_temp, "nTF1coeff/i");
657
658 Settings->Fill();
659
660 SplineFile->cd();
661 Settings->Write();
662
663 TTree *SplineTree = new TTree("SplineTree", "SplineTree");
664 // Create a branch for the SplineMonoStruct object
665 SplineTree->Branch("SplineObject", &cpu_spline_handler);
666 SplineTree->Fill();
667 SplineTree->Write();
668 delete SplineTree;
669
670 float coeff_tf1 = 0.;
671 Monolith_TF1->Branch("cpu_coeff_TF1_many", &coeff_tf1, "cpu_coeff_TF1_many/F");
672 for(unsigned int i = 0; i < nTF1coeff; i++)
673 {
674 coeff_tf1 = cpu_coeff_TF1_many[i];
675 Monolith_TF1->Fill();
676 }
677 SplineFile->cd();
678 Monolith_TF1->Write();
679
680 unsigned int nParamPerEvent = 0;
681 unsigned int nParamPerEvent_tf1 = 0;
682
683 EventInfo->Branch("cpu_nParamPerEvent", &nParamPerEvent, "cpu_nParamPerEvent/i");
684 EventInfo->Branch("cpu_nParamPerEvent_tf1", &nParamPerEvent_tf1, "cpu_nParamPerEvent_tf1/i");
685
686 for(unsigned int i = 0; i < 2*NEvents; i++)
687 {
688 nParamPerEvent = cpu_nParamPerEvent[i];
689 nParamPerEvent_tf1 = cpu_nParamPerEvent_tf1[i];
690 EventInfo->Fill();
691 }
692 SplineFile->cd();
693 EventInfo->Write();
694
695 M3::int_t nPoints = 0;
696 float xtemp[20];
697 FastSplineInfoTree->Branch("nPts", &nPoints, "nPts/I");
698 FastSplineInfoTree->Branch("xPts", xtemp, "xPts[nPts]/F");
699
700 for (M3::int_t i = 0; i < nParams; ++i)
701 {
702 nPoints = SplineInfoArray[i].nPts;
703
704 for (M3::int_t k = 0; k < SplineInfoArray[i].nPts; ++k)
705 {
706 xtemp[k] = float(SplineInfoArray[i].xPts[k]);
707 }
708 FastSplineInfoTree->Fill();
709 }
710
711 SplineFile->cd();
712 FastSplineInfoTree->Write();
713
714 delete Settings;
715 delete Monolith_TF1;
716 delete XKnots;
717 delete EventInfo;
718 delete FastSplineInfoTree;
719 SplineFile->Close();
720}
721
722// *****************************************
723// Destructor
724// Cleans up the allocated GPU memory
726// *****************************************
727 #ifdef MaCh3_CUDA
729 #ifndef Weight_On_SplineBySpline_Basis
731 #endif
732 );
733
734 //KS: Since we declared them using CUDA alloc we have to free memory using also cuda functions
736
737 delete gpu_spline_handler;
738 #else
739 if(SplineSegments != nullptr) delete[] SplineSegments;
740 if(ParamValues != nullptr) delete[] ParamValues;
741 if(cpu_total_weights != nullptr) delete[] cpu_total_weights;
742 #endif
743
744 if(cpu_weights != nullptr) delete[] cpu_weights;
745 if(cpu_weights_spline_var != nullptr) delete[] cpu_weights_spline_var;
746 if(cpu_weights_tf1_var != nullptr) delete[] cpu_weights_tf1_var;
747
748 if(cpu_spline_handler != nullptr) delete cpu_spline_handler;
749}
750
751// *****************************************
752// Get the spline coefficients from the TSpline3 so that we can load ONLY these onto the GPU, not the whole TSpline3 object
753// This loads up coefficients into two arrays: one x array and one yabcd array
754// This should maximize our cache hits!
755void SMonolith::getSplineCoeff_SepMany(TSpline3_red* &spl, int &nPoints, float *& xArray, float *& manyArray) {
756// *****************************************
757 // Initialise all arrays to 1.0
758 for (int i = 0; i < _max_knots; ++i) {
759 xArray[i] = 1.0;
760 for (int j = 0; j < _nCoeff_; j++) {
761 manyArray[i*_nCoeff_+j] = 1.0;
762 }
763 }
764 // Get number of points in spline
765 int Np = spl->GetNp();
766 // If spline is flat, set number of knots to 1.0,
767 // This is used later to expedite the calculations for flat splines
768 // tmpArray[0] is number of knots
769 nPoints = Np;
770 if (Np > _max_knots) {
771 MACH3LOG_ERROR("Error, number of points is greater than saved {}", _max_knots);
772 MACH3LOG_ERROR("This _WILL_ cause problems with GPU splines and _SHOULD_ be fixed!");
773 MACH3LOG_ERROR("nPoints = {}, _max_knots = {}", nPoints, _max_knots);
774 throw MaCh3Exception(__FILE__ , __LINE__ );
775 }
776
777 // The coefficients we're writing to
778 M3::float_t x, y, b, c, d;
779 // TSpline3 can only take doubles, not floats
780 // But our GPU is slow with doubles, so need to cast to float
781 for(int i = 0; i < Np; i++) {
782 // Get the coefficients from the TSpline3 object
783 spl->GetCoeff(i, x, y, b, c, d);
784 // Write the arrays
785 xArray[i] = float(x);
786 manyArray[i*_nCoeff_] = float(y); // 4 because manyArray stores y,b,c,d
787 manyArray[i*_nCoeff_+1] = float(b);
788 manyArray[i*_nCoeff_+2] = float(c);
789 manyArray[i*_nCoeff_+3] = float(d);
790 if((xArray[i] == -999) || (manyArray[i*_nCoeff_] == -999) || (manyArray[i*4+1] == -999) || (manyArray[i*_nCoeff_+2] == -999) || (manyArray[i*_nCoeff_+3] == -999)){
791 MACH3LOG_ERROR("*********** Bad params in getSplineCoeff_SepMany() ************");
792 MACH3LOG_ERROR("pre cast to float (x, y, b, c, d) = {:.2f}, {:.2f}, {:.2f}, {:.2f}, {:.2f}", x, y, b, c, d);
793 MACH3LOG_ERROR("pre cast to float (x, y, b, c, d) = {:.2f}, {:.2f}, {:.2f}, {:.2f}, {:.2f}", xArray[i], manyArray[i*4], manyArray[i*4+1], manyArray[i*4+2], manyArray[i*_nCoeff_+3]);
794 MACH3LOG_ERROR("This will cause problems when preparing for GPU");
795 MACH3LOG_ERROR("***************************************************************");
796 }
797 }
798}
799
800#ifdef MaCh3_CUDA
801// *****************************************
802// Tell the GPU to evaluate the weights
803// Load up the two x,{y,b,c,d} arrays into memory and have GPU read them with more coalescence instead of one monolithic array
804// This should be used when we're using separate x,y,a,b,c,d arrays
805// Also pass the segments for the parameter along with their parameter values
806// This avoids doing lots of binary searches on the GPU
807void SMonolith::Evaluate() {
808// *****************************************
809 // There's a parameter mapping that goes from spline parameter to a global parameter index
810 // Find the spline segments
812
813 // The main call to the GPU
815 #ifdef Weight_On_SplineBySpline_Basis
818 #else
820 #endif
824 NTF1_valid);
825
826 //KS: Normally it does nothing, in case you want to have weight for each spline it does the mapping, used mostly for debugging
828}
829#else
830//If CUDA is not enabled do the same on CPU
831// *****************************************
833// *****************************************
834 // There's a parameter mapping that goes from spline parameter to a global parameter index
835 // Find the spline segments
837
838 //KS: Huge MP loop over all valid splines
840
841 //KS: Huge MP loop over all events calculating total weight
843}
844#endif
845
846//*********************************************************
848//*********************************************************
849 #ifdef MULTITHREAD
850 //KS: Open parallel region
851 #pragma omp parallel
852 {
853 #endif
854 //KS: First we calculate
855 #ifdef MULTITHREAD
856 #pragma omp for simd nowait
857 #endif
858 for (unsigned int splineNum = 0; splineNum < NSplines_valid; ++splineNum)
859 {
860 //CW: Which Parameter we are accessing
861 const short int Param = cpu_spline_handler->paramNo_arr[splineNum];
862
863 //CW: Avoids doing costly binary search on GPU
864 const short int segment = SplineSegments[Param];
865
866 //KS: Segment for coeff_x is simply parameter*max knots + segment as each parameters has the same spacing
867 const short int segment_X = short(Param*_max_knots+segment);
868
869 //KS: Find knot position in out monolithical structure
870 const unsigned int CurrentKnotPos = cpu_spline_handler->nKnots_arr[splineNum]*_nCoeff_+segment*_nCoeff_;
871
872 // We've read the segment straight from CPU and is saved in segment_gpu
873 // polynomial parameters from the monolithic splineMonolith
874 const float fY = cpu_spline_handler->coeff_many[CurrentKnotPos];
875 const float fB = cpu_spline_handler->coeff_many[CurrentKnotPos + 1];
876 const float fC = cpu_spline_handler->coeff_many[CurrentKnotPos + 2];
877 const float fD = cpu_spline_handler->coeff_many[CurrentKnotPos + 3];
878 // The is the variation itself (needed to evaluate variation - stored spline point = dx)
879 const float dx = ParamValues[Param] - cpu_spline_handler->coeff_x[segment_X];
880
881 //CW: Wooow, let's use some fancy intrinsic and pull down the processing time by <1% from normal multiplication! HURRAY
882 cpu_weights_spline_var[splineNum] = fmaf(dx, fmaf(dx, fmaf(dx, fD, fC), fB), fY);
883 // Or for the more "easy to read" version:
884 //cpu_weights_spline_var[splineNum] = (fY+dx*(fB+dx*(fC+dx*fD)));
885 }
886
887 #ifdef MULTITHREAD
888 #pragma omp for simd
889 #endif
890 for (unsigned int tf1Num = 0; tf1Num < NTF1_valid; ++tf1Num)
891 {
892 // The is the variation itself (needed to evaluate variation - stored spline point = dx)
893 const float x = ParamValues[cpu_paramNo_TF1_arr[tf1Num]];
894
895 // Read the coefficients
896 const unsigned int TF1_Index = tf1Num * _nTF1Coeff_;
897 const float a = cpu_coeff_TF1_many[TF1_Index];
898 const float b = cpu_coeff_TF1_many[TF1_Index + 1];
899
900 cpu_weights_tf1_var[tf1Num] = fmaf(a, x, b);
901 // cpu_weights_tf1_var[tf1Num] = a*x + b;
902 //cpu_weights_tf1_var[splineNum] = 1 + a*x + b*x*x + c*x*x*x + d*x*x*x*x + e*x*x*x*x*x;
903 }
904 #ifdef MULTITHREAD
905 //KS: End parallel region
906 }
907 #endif
908}
909
910//*********************************************************
911//KS: Calc total event weight on CPU
913//*********************************************************
914#ifndef Weight_On_SplineBySpline_Basis
915 #ifdef MULTITHREAD
916 #pragma omp parallel for
917 #endif
918 for (unsigned int EventNum = 0; EventNum < NEvents; ++EventNum)
919 {
920 float totalWeight = 1.0f; // Initialize total weight for each event
921
922 const unsigned int Offset = 2 * EventNum;
923
924 // Extract the parameters for the current event
925 const unsigned int startIndex = cpu_nParamPerEvent[Offset + 1];
926 const unsigned int numParams = cpu_nParamPerEvent[Offset];
927
928 // Compute total weight for the current event
929 #ifdef MULTITHREAD
930 #pragma omp simd
931 #endif
932 for (unsigned int id = 0; id < numParams; ++id) {
933 totalWeight *= cpu_weights_spline_var[startIndex + id];
934 }
935 //Now TF1
936 // Extract the parameters for the current event
937 const unsigned int startIndex_tf1 = cpu_nParamPerEvent_tf1[Offset + 1];
938 const unsigned int numParams_tf1 = cpu_nParamPerEvent_tf1[Offset];
939
940 // Compute total weight for the current event
941 #ifdef MULTITHREAD
942 #pragma omp simd
943 #endif
944 for (unsigned int id = 0; id < numParams_tf1; ++id) {
945 totalWeight *= cpu_weights_tf1_var[startIndex_tf1 + id];
946 }
947
948 // Store the total weight for the current event
949 cpu_total_weights[EventNum] = totalWeight;
950 }
951#else
952 //KS: Name is confusing but what it does it make a nice mapping used for debugging
954#endif
955}
956
957//*********************************************************
958//KS: Normally it does nothing, in case you want to have weight for each spline it does the mapping, used mostly for debugging
960//*********************************************************
961#ifdef Weight_On_SplineBySpline_Basis
962 // Multi-thread here because _numIndex is really quite large!
963 #ifdef MULTITHREAD
964 #pragma omp parallel for
965 #endif
966 for (unsigned int i = 0; i < NSplines_total_large; ++i) {
967 if (index_spline_cpu[i] >= 0) {
969 } else if (index_TF1_cpu[i] >= 0) {
971 } else {
972 cpu_weights[i] = 1.;
973 }
974 }
975#endif
976}
977
978//*********************************************************
979//KS: Print info about how much knots etc has been initialised
981//*********************************************************
982 unsigned int event_size_max = _max_knots * nParams;
983
984 MACH3LOG_INFO("--- INITIALISED Spline Monolith ---");
985 MACH3LOG_INFO("{} events with {} splines", NEvents, NSplines_valid);
986 MACH3LOG_INFO("On average {:.2f} splines per event ({}/{})", float(NSplines_valid)/float(NEvents), NSplines_valid, NEvents);
987 MACH3LOG_INFO("Size of x array = {:.4f} MB", double(sizeof(float)*event_size_max)/1.E6);
988 MACH3LOG_INFO("Size of coefficient (y,b,c,d) array = {:.2f} MB", double(sizeof(float)*nKnots*_nCoeff_)/1.E6);
989 MACH3LOG_INFO("Size of parameter # array = {:.2f} MB", double(sizeof(short int)*NSplines_valid)/1.E6);
990
991 MACH3LOG_INFO("On average {:.2f} TF1 per event ({}/{})", float(NTF1_valid)/float(NEvents), NTF1_valid, NEvents);
992 MACH3LOG_INFO("Size of TF1 coefficient (a,b,c,d,e) array = {:.2f} MB", double(sizeof(float)*NTF1_valid*_nTF1Coeff_)/1.E6);
993}
994
995//*********************************************************
996//KS: After calculations are done on GPU we copy memory to CPU. This operation is asynchronous meaning while memory is being copied some operations are being carried. Memory must be copied before actual reweight. This function make sure all has been copied.
998//*********************************************************
999 #ifdef MaCh3_CUDA
1001 #endif
1002}
int size
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:25
#define MACH3LOG_INFO
Definition: MaCh3Logger.h:23
#define MACH3LOG_WARN
Definition: MaCh3Logger.h:24
@ kTF1_red
Uses TF1_red for interpolation.
@ kTSpline3_red
Uses TSpline3_red for interpolation.
void CleanVector(std::vector< T > &vec)
Generic cleanup function.
#define _nCoeff_
KS: We store coefficients {y,b,c,d} in one array one by one, this is only to define it once rather th...
Definition: SplineCommon.h:13
#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
Custom exception class for MaCh3 errors.
Class responsible for calculating spline weight on GPU.
__host__ void CleanupGPU_SplineMonolith(float *cpu_total_weights)
This function deallocates the resources allocated for the separate {x} and {ybcd} arrays in the and T...
__host__ void InitGPU_Vals(float **vals)
Allocate memory for spline segments.
__host__ void InitGPU_SplineMonolith(float **cpu_total_weights, int n_events, unsigned int total_nknots, unsigned int n_splines, unsigned int n_tf1, int Eve_size)
Allocate memory on gpu for spline monolith.
__host__ void RunGPU_SplineMonolith(float *cpu_total_weights, float *vals, short int *segment, const unsigned int h_n_splines, const unsigned int h_n_tf1)
Run the GPU code for the separate many arrays. As in separate {x}, {y,b,c,d} arrays Pass the segment ...
__host__ void InitGPU_Segments(short int **segment)
Allocate memory for spline segments.
__host__ void CopyToGPU_SplineMonolith(SplineMonoStruct *cpu_spline_handler, std::vector< float > cpu_many_array_TF1, std::vector< short int > cpu_paramNo_arr_TF1, int n_events, std::vector< unsigned int > cpu_nParamPerEvent, std::vector< unsigned int > cpu_nParamPerEvent_TF1, int n_params, unsigned int n_splines, short int spline_size, unsigned int total_nknots, unsigned int n_tf1)
Copies data from CPU to GPU for the spline monolith.
__host__ void CleanupGPU_Segments(short int *segment, float *vals)
Clean up pinned variables at CPU.
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.
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.
SMonolith(std::vector< std::vector< TResponseFunction_red * > > &MasterSpline, const std::vector< RespFuncType > &SplineType, const bool SaveFlatTree=false)
Constructor.
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.
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!)
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
void FindSplineSegment()
CW:Code used in step by step reweighting, Find Spline Segment for each param.
Definition: SplineBase.cpp:24
short int * SplineSegments
Definition: SplineBase.h:60
std::vector< FastSplineInfo > SplineInfoArray
Definition: SplineBase.h:57
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
void getTF1Coeff(TF1_red *&spl, int &nPoints, float *&coeffs)
CW: Gets the polynomial coefficients for TF1.
Definition: SplineBase.cpp:106
CW: A reduced TF1 class only. Only saves parameters for each TF1 and how many parameters each paramet...
int GetSize()
Get the size.
KS: A reduced ResponseFunction Generic function used for evaluating weight.
CW: Reduced TSpline3 class.
void GetKnot(int i, M3::float_t &xtmp, M3::float_t &ytmp)
M3::int_t GetNp() override
CW: Get the number of points.
void GetCoeff(int segment, M3::float_t &x, M3::float_t &y, M3::float_t &b, M3::float_t &c, M3::float_t &d)
CW: Get the coefficient of a given segment.
__host__ void SynchroniseSplines()
Make sure all Cuda threads finished execution.
MaCh3 event-by-event cross-section spline code.
double float_t
Definition: Core.h:28
int int_t
Definition: Core.h:29
KS: Struct storing information for spline monolith.
Definition: SplineCommon.h:30
std::vector< unsigned int > nKnots_arr
KS: CPU Number of knots per spline.
Definition: SplineCommon.h:41
std::vector< float > coeff_x
KS: CPU arrays to hold X coefficient.
Definition: SplineCommon.h:35
std::vector< float > coeff_many
CPU arrays to hold other coefficients.
Definition: SplineCommon.h:38
std::vector< short int > paramNo_arr
CW: CPU array with the number of points per spline (not per spline point!)
Definition: SplineCommon.h:44