MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
BinnedSplineHandler.cpp
Go to the documentation of this file.
2#include <memory>
3
4#pragma GCC diagnostic ignored "-Wuseless-cast"
5#pragma GCC diagnostic ignored "-Wfloat-conversion"
6
8#include "TROOT.h"
9#include "TKey.h"
11
12//****************************************
14//****************************************
15 if (!xsec_) {
16 MACH3LOG_ERROR("Trying to create BinnedSplineHandler with uninitialized covariance object");
17 throw MaCh3Exception(__FILE__, __LINE__);
18 }
19 xsec = xsec_;
20
21 if (!Modes_) {
22 MACH3LOG_ERROR("Trying to create BinnedSplineHandler with uninitialized MaCh3Modes object");
23 throw MaCh3Exception(__FILE__, __LINE__);
24 }
25 Modes = Modes_;
26
27 // Keep these in class scope, important for using 1 monolith/sample!
28 MonolithIndex = 0; //Keeps track of the monolith index we're on when filling arrays (declared here so we can have multiple FillSampleArray calls)
29 CoeffIndex = 0; //Keeps track of our indexing the coefficient arrays [x, ybcd]
30}
31//****************************************
33//****************************************
34 if(manycoeff_arr != nullptr) delete[] manycoeff_arr;
35 if(xcoeff_arr != nullptr) delete[] xcoeff_arr;
36 if(SplineSegments != nullptr) delete[] SplineSegments;
37 if(ParamValues != nullptr) delete[] ParamValues;
38}
39//****************************************
41//****************************************
42 //Call once everything's been allocated in SampleHandlerFDBase, cleans up junk from memory!
43 //Not a huge saving but it's better than leaving everything up to the compiler
44 MACH3LOG_INFO("Cleaning up spline memory");
59 if(isflatarray) delete [] isflatarray;
60}
61
62//****************************************
63void BinnedSplineHandler::AddSample(const std::string& SampleName,
64 const std::vector<std::string>& OscChanFileNames,
65 const std::vector<std::string>& SplineVarNames)
66//Adds samples to the large array
67//****************************************
68{
69 SampleNames.push_back(SampleName);
70 Dimensions.push_back(int(SplineVarNames.size()));
71 DimensionLabels.push_back(SplineVarNames);
72
73 int nSplineParam = xsec->GetNumParamsFromSampleName(SampleName, SystType::kSpline);
74 nSplineParams.push_back(nSplineParam);
75
76 //This holds the global index of the spline i.e. 0 -> _fNumPar
77 std::vector<int> GlobalSystIndex_Sample = xsec->GetGlobalSystIndexFromSampleName(SampleName, SystType::kSpline);
78 //Keep track of this for all the samples
79 GlobalSystIndex.push_back(GlobalSystIndex_Sample);
80
81 std::vector<SplineInterpolation> SplineInterpolation_Sample = xsec->GetSplineInterpolationFromSampleName(SampleName);
82 // Keep track of this for all samples
83 SplineInterpolationTypes.push_back(SplineInterpolation_Sample);
84
85 std::vector<std::string> SplineFileParPrefixNames_Sample = xsec->GetSplineParsNamesFromSampleName(SampleName);
86 SplineFileParPrefixNames.push_back(SplineFileParPrefixNames_Sample);
87
88 MACH3LOG_INFO("Create SplineModeVecs_Sample");
89 std::vector<std::vector<int>> SplineModeVecs_Sample = StripDuplicatedModes(xsec->GetSplineModeVecFromSampleName(SampleName));
90 MACH3LOG_INFO("SplineModeVecs_Sample is of size {}", SplineModeVecs_Sample.size());
91 SplineModeVecs.push_back(SplineModeVecs_Sample);
92
93 MACH3LOG_INFO("SplineModeVecs is of size {}", SplineModeVecs.size());
94
95 int nOscChan = int(OscChanFileNames.size());
96 nOscChans.push_back(nOscChan);
97
98 PrintSampleDetails(SampleName);
99
100 std::vector<std::vector<TAxis *>> SampleBinning(nOscChan);
101 for (int iOscChan = 0; iOscChan < nOscChan; iOscChan++)
102 {
103 SampleBinning[iOscChan] = FindSplineBinning(OscChanFileNames[iOscChan], SampleName);
104 }
105 MACH3LOG_INFO("#----------------------------------------------------------------------------------------------------------------------------------#");
106 SplineBinning.push_back(SampleBinning);
107
108 BuildSampleIndexingArray(SampleName);
109 PrintArrayDetails(SampleName);
110 MACH3LOG_INFO("#----------------------------------------------------------------------------------------------------------------------------------#");
111
112 FillSampleArray(SampleName, OscChanFileNames);
113 MACH3LOG_INFO("#----------------------------------------------------------------------------------------------------------------------------------#");
114}
115
116//****************************************
118//****************************************
119{
122
124 MACH3LOG_ERROR("Something's gone wrong when we tried to get the size of your monolith");
125 MACH3LOG_ERROR("MonolithSize is {}", MonolithSize);
126 MACH3LOG_ERROR("MonolithIndex is {}", MonolithIndex);
127 throw MaCh3Exception(__FILE__ , __LINE__ );
128 }
129
130 MACH3LOG_INFO("Now transferring splines to a monolith if size {}", MonolithSize);
131
134 isflatarray = new bool[MonolithSize];
135
138
139 for (unsigned int iSample = 0; iSample < indexvec.size(); iSample++)
140 { // Loop over sample
141 for (unsigned int iOscChan = 0; iOscChan < indexvec[iSample].size(); iOscChan++)
142 { // Loop over oscillation channels
143 for (unsigned int iSyst = 0; iSyst < indexvec[iSample][iOscChan].size(); iSyst++)
144 { // Loop over systematics
145 for (unsigned int iMode = 0; iMode < indexvec[iSample][iOscChan][iSyst].size(); iMode++)
146 { // Loop over modes
147 for (unsigned int iVar1 = 0; iVar1 < indexvec[iSample][iOscChan][iSyst][iMode].size(); iVar1++)
148 { // Loop over first dimension
149 for (unsigned int iVar2 = 0; iVar2 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1].size(); iVar2++)
150 { // Loop over second dimension
151 for (unsigned int iVar3 = 0; iVar3 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2].size(); iVar3++)
152 {
153 int splineindex=indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2][iVar3];
154 weightvec_Monolith[splineindex] = 1.0;
155
156 bool foundUniqueSpline = false;
157 for (int iUniqueSyst = 0; iUniqueSyst < nParams; iUniqueSyst++)
158 {
159 if (SplineFileParPrefixNames[iSample][iSyst] == UniqueSystNames[iUniqueSyst])
160 {
161 uniquesplinevec_Monolith[splineindex] = iUniqueSyst;
162 foundUniqueSpline = true;
163 }
164 }//unique syst loop end
165
166 if (!foundUniqueSpline)
167 {
168 MACH3LOG_ERROR("Unique spline index not found");
169 MACH3LOG_ERROR("For Spline {}", SplineFileParPrefixNames[iSample][iSyst]);
170 MACH3LOG_ERROR("Couldn't match {} with any of the following {} systs:", SplineFileParPrefixNames[iSample][iSyst], nParams);
171 for (int iUniqueSyst = 0; iUniqueSyst < nParams; iUniqueSyst++)
172 {
173 MACH3LOG_ERROR("{},", UniqueSystNames.at(iUniqueSyst));
174 }//unique syst loop end
175 throw MaCh3Exception(__FILE__ , __LINE__ );
176 }
177
178 int splineKnots;
179 if(splinevec_Monolith[splineindex]){
180 isflatarray[splineindex]=false;
181 splineKnots=splinevec_Monolith[splineindex]->GetNp();
182
183 //Now to fill up our coefficient arrayss
184 M3::float_t* tmpXCoeffArr = new M3::float_t[splineKnots];
185 M3::float_t* tmpManyCoeffArr = new M3::float_t[splineKnots*4];
186
187 int iCoeff=coeffindexvec[splineindex];
188 getSplineCoeff_SepMany(splineindex, tmpXCoeffArr, tmpManyCoeffArr);
189
190 for(int i = 0; i < splineKnots; i++){
191 xcoeff_arr[iCoeff+i]=tmpXCoeffArr[i];
192
193 for(int j=0; j<4; j++){
194 manycoeff_arr[(iCoeff+i)*4+j]=tmpManyCoeffArr[i*4+j];
195 }
196 }
197 delete[] tmpXCoeffArr;
198 delete[] tmpManyCoeffArr;
199 }
200 else {
201 isflatarray[splineindex]=true;
202 }
203 }//3d loop end
204 }//2d loop end
205 }//1d loop end
206 }//mode loop end
207 }//syst2 loop end
208 }//osc loop end
209 }//syst1 loop end
210}
211
212// *****************************************
214// *****************************************
215 // There's a parameter mapping that goes from spline parameter to a global parameter index
216 // Find the spline segments
218
219 //KS: Huge MP loop over all valid splines
221
222 //KS: Huge MP loop over all events calculating total weight
224}
225
226//****************************************
228//****************************************
229{
230 #ifdef MULTITHREAD
231 #pragma omp parallel for simd
232 #endif
233 for (size_t iCoeff = 0; iCoeff < uniquecoeffindices.size(); ++iCoeff)
234 {
235 const int iSpline = uniquecoeffindices[iCoeff];
236 const short int uniqueIndex = short(uniquesplinevec_Monolith[iSpline]);
237 const short int currentsegment = short(SplineSegments[uniqueIndex]);
238
239 const int segCoeff = coeffindexvec[iSpline]+currentsegment;
240 const int coeffOffset = segCoeff * 4;
241 // These are what we can extract from the TSpline3
242 const M3::float_t y = manycoeff_arr[coeffOffset+kCoeffY];
243 const M3::float_t b = manycoeff_arr[coeffOffset+kCoeffB];
244 const M3::float_t c = manycoeff_arr[coeffOffset+kCoeffC];
245 const M3::float_t d = manycoeff_arr[coeffOffset+kCoeffD];
246
247 // Get the variation for this reconfigure for the ith parameter
249 const M3::float_t xvar = (*SplineInfoArray[uniqueIndex].splineParsPointer);
250 // The Delta(x) = xvar - x
251 const M3::float_t dx = xvar - xcoeff_arr[segCoeff];
252
253 //Speedy 1% time boost https://en.cppreference.com/w/c/numeric/math/fma (see ND code!)
254 M3::float_t weight = M3::fmaf_t(dx, M3::fmaf_t(dx, M3::fmaf_t(dx, d, c), b), y);
255 //This is the speedy version of writing dx^3+b*dx^2+c*dx+d
256
257 //ETA - do we need this? We check later for negative weights and I wonder if this is even
258 //possible with the fmaf line above?
259 if(weight < 0){weight = 0.;} //Stops is getting negative weights
260
261 weightvec_Monolith[iSpline] = weight;
262 }
263}
264
265//****************************************
266//Creates an array to be filled with monolith indexes for each sample (allows for indexing between 7D binning and 1D Vector)
267//Only need 1 indexing array everything else interfaces with this to get binning properties
268void BinnedSplineHandler::BuildSampleIndexingArray(const std::string& SampleName)
269//****************************************
270{
271 int iSample = getSampleIndex(SampleName);
272 int nSplineSysts = nSplineParams[iSample];
273 int nOscChannels = nOscChans[iSample];
274
275 // Resize the main indexing structure
276 indexvec.emplace_back(nOscChannels);
277
278 for (int iOscChan = 0; iOscChan < nOscChannels; ++iOscChan)
279 {
280 indexvec.back()[iOscChan].resize(nSplineSysts);
281 for (int iSyst = 0; iSyst < nSplineSysts; ++iSyst)
282 {
283 int nModesInSyst = static_cast<int>(SplineModeVecs[iSample][iSyst].size());
284 indexvec.back()[iOscChan][iSyst].resize(nModesInSyst);
285
286 for (int iMode = 0; iMode < nModesInSyst; ++iMode)
287 {
288 const int nBins1 = SplineBinning[iSample][iOscChan][0]->GetNbins();
289 const int nBins2 = SplineBinning[iSample][iOscChan][1]->GetNbins();
290 const int nBins3 = SplineBinning[iSample][iOscChan][2]->GetNbins();
291
292 indexvec.back()[iOscChan][iSyst][iMode]
293 .resize(nBins1,std::vector<std::vector<int>>(nBins2, std::vector<int>(nBins3, 0)));
294 }
295 } // end of iSyst loop
296 } // end of iOscChan loop
297}
298
299//****************************************
300std::vector<TAxis *> BinnedSplineHandler::FindSplineBinning(const std::string& FileName, const std::string& SampleName)
301//****************************************
302{
303 std::vector<TAxis *> ReturnVec;
304 int iSample=getSampleIndex(SampleName);
305
306 //Try declaring these outside of TFile so they aren't owned by File
307 int nDummyBins = 1;
308 const double DummyEdges[2] = {-1e15, 1e15};
309 TAxis* DummyAxis = new TAxis(nDummyBins, DummyEdges);
310 TH2F* Hist2D = nullptr;
311 TH3F* Hist3D = nullptr;
312
313 auto File = std::unique_ptr<TFile>(TFile::Open(FileName.c_str(), "READ"));
314 if (!File || File->IsZombie())
315 {
316 MACH3LOG_ERROR("File {} not found", FileName);
317 MACH3LOG_ERROR("This is caused by something here! {} : {}", __FILE__, __LINE__);
318 throw MaCh3Exception(__FILE__ , __LINE__ );
319 }
320
321 MACH3LOG_INFO("Finding binning for:");
322 MACH3LOG_INFO("{}", FileName);
323
324 std::string TemplateName = "dev_tmp_0_0";
325 TObject *Obj = File->Get(TemplateName.c_str());
326 //If you can't find dev_tmp_0_0 then this will cause a problem
327 if (!Obj)
328 {
329 TemplateName = "dev_tmp.0.0";
330 Obj = File->Get(TemplateName.c_str());
331 if (!Obj)
332 {
333 MACH3LOG_ERROR("Could not find dev_tmp_0_0 in spline file. Spline binning cannot be set!");
334 MACH3LOG_ERROR("FileName: {}", FileName);
335 throw MaCh3Exception(__FILE__ , __LINE__ );
336 }
337 }
338
339 //Now check if dev_tmp_0_0 is a TH2 i.e. specifying the dimensions of the splines is 2D
340 bool isHist2D = Obj->IsA() == TH2F::Class();
341 //For T2K annoyingly all objects are TH3Fs
342 bool isHist3D = Obj->IsA() == TH3F::Class();
343 if (!isHist2D && !isHist3D)
344 {
345 MACH3LOG_ERROR("Object doesn't inherit from either TH2D and TH3D - Odd A");
346 throw MaCh3Exception(__FILE__ , __LINE__ );
347 }
348
349 if (isHist2D)
350 {
351 if (Dimensions[iSample] != 2)
352 {
353 MACH3LOG_ERROR("Trying to load a 2D spline template when nDim={}", Dimensions[iSample]);
354 throw MaCh3Exception(__FILE__, __LINE__);
355 }
356 Hist2D = File->Get<TH2F>(TemplateName.c_str());
357 }
358
359 if (isHist3D)
360 {
361 Hist3D = File->Get<TH3F>((TemplateName.c_str()));
362
363 if (Dimensions[iSample] != 3 && Hist3D->GetZaxis()->GetNbins() != 1)
364 {
365 MACH3LOG_ERROR("Trying to load a 3D spline template when nDim={}", Dimensions[iSample]);
366 throw MaCh3Exception(__FILE__ , __LINE__ );
367 }
368 Hist3D = File->Get<TH3F>(TemplateName.c_str());
369 }
370
371 if (Dimensions[iSample] == 2) {
372 if(isHist2D){
373 ReturnVec.push_back(static_cast<TAxis*>(Hist2D->GetXaxis()->Clone()));
374 ReturnVec.push_back(static_cast<TAxis*>(Hist2D->GetYaxis()->Clone()));
375 ReturnVec.push_back(static_cast<TAxis*>(DummyAxis->Clone()));
376 } else if (isHist3D) {
377 ReturnVec.push_back(static_cast<TAxis*>(Hist3D->GetXaxis()->Clone()));
378 ReturnVec.push_back(static_cast<TAxis*>(Hist3D->GetYaxis()->Clone()));
379 ReturnVec.push_back(static_cast<TAxis*>(DummyAxis->Clone()));
380 }
381 } else if (Dimensions[iSample] == 3) {
382 ReturnVec.push_back(static_cast<TAxis*>(Hist3D->GetXaxis()->Clone()));
383 ReturnVec.push_back(static_cast<TAxis*>(Hist3D->GetYaxis()->Clone()));
384 ReturnVec.push_back(static_cast<TAxis*>(Hist3D->GetZaxis()->Clone()));
385 } else {
386 MACH3LOG_ERROR("Number of dimensions not valid! Given: {}", Dimensions[iSample]);
387 throw MaCh3Exception(__FILE__ , __LINE__ );
388 }
389
390 for (unsigned int iAxis = 0; iAxis < ReturnVec.size(); ++iAxis) {
391 PrintBinning(ReturnVec[iAxis]);
392 }
393
394 MACH3LOG_INFO("Left PrintBinning now tidying up");
395 delete DummyAxis;
396
397 return ReturnVec;
398}
399
400//****************************************
402//****************************************
403{
404 int SampleCounter_NonFlat = 0;
405 int SampleCounter_All = 0;
406 int FullCounter_NonFlat = 0;
407 int FullCounter_All = 0;
408
409 for (unsigned int iSample = 0; iSample < indexvec.size(); iSample++)
410 { // Loop over systematics
411 SampleCounter_NonFlat = 0;
412 SampleCounter_All = 0;
413 std::string SampleName = SampleNames[iSample];
414 for (unsigned int iOscChan = 0; iOscChan < indexvec[iSample].size(); iOscChan++)
415 { // Loop over oscillation channels
416 for (unsigned int iSyst = 0; iSyst < indexvec[iSample][iOscChan].size(); iSyst++)
417 { // Loop over systematics
418 for (unsigned int iMode = 0; iMode < indexvec[iSample][iOscChan][iSyst].size(); iMode++)
419 { // Loop over modes
420 for (unsigned int iVar1 = 0; iVar1 < indexvec[iSample][iOscChan][iSyst][iMode].size(); iVar1++)
421 { // Loop over first dimension
422 for (unsigned int iVar2 = 0; iVar2 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1].size(); iVar2++)
423 { // Loop over second dimension
424 for (unsigned int iVar3 = 0; iVar3 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2].size(); iVar3++)
425 { // Loop over third dimension
426 if (isValidSplineIndex(SampleName, iOscChan, iSyst, iMode, iVar1, iVar2, iVar3))
427 {
428 int splineindex = indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2][iVar3];
429 if (splinevec_Monolith[splineindex])
430 {
431 SampleCounter_NonFlat += 1;
432 }
433 SampleCounter_All += 1;
434 }
435 }
436 }
437 }
438 }
439 }
440 }
441 MACH3LOG_DEBUG("{:<10} has {:<10} splines, of which {:<10} are not flat", SampleNames[iSample], SampleCounter_All, SampleCounter_NonFlat);
442
443 FullCounter_NonFlat += SampleCounter_NonFlat;
444 FullCounter_All += SampleCounter_All;
445 }
446
447 if (Verbosity > 0)
448 {
449 MACH3LOG_INFO("Total number of splines loaded: {}", FullCounter_All);
450 MACH3LOG_INFO("Total number of non-flat splines loaded: {}", FullCounter_NonFlat);
451 }
452
453 if (NonFlat) {
454 return FullCounter_NonFlat;
455 } else {
456 return FullCounter_All;
457 }
458}
459
460//****************************************
462//****************************************
463 std::vector<TSpline3_red*> UniqueSystSplines;
464
465 // DB Find all the Unique systs across each sample and oscillation channel
466 // This assumes that each occurence of the same systematic spline has the same knot spacing
467 // Which is a reasonable assumption for the current implementation of spline evaluations
468 for (unsigned int iSample = 0; iSample < indexvec.size(); iSample++)
469 { // Loop over systematics
470 for (unsigned int iSyst = 0; iSyst < indexvec[iSample][0].size(); iSyst++)
471 { // Loop over systematics
472 std::string SystName = SplineFileParPrefixNames[iSample][iSyst];
473 bool FoundSyst = false;
474
475 //ETA - this always seems to be empty to begin with??
476 //so this loop never gets used?
477 for (unsigned int iFoundSyst = 0; iFoundSyst < UniqueSystNames.size(); iFoundSyst++)
478 {
479 if (SystName == UniqueSystNames[iFoundSyst])
480 {
481 FoundSyst = true;
482 }
483 }
484 if (!FoundSyst)
485 {
486 bool FoundNonFlatSpline = false;
487 //This is all basically to check that you have some resposne from a spline
488 //systematic somewhere
489 for (unsigned int iOscChan = 0; iOscChan < indexvec[iSample].size(); iOscChan++)
490 { // Loop over oscillation channels
491 for (unsigned int iMode = 0; iMode < indexvec[iSample][iOscChan][iSyst].size(); iMode++)
492 { // Loop over modes
493 for (unsigned int iVar1 = 0; iVar1 < indexvec[iSample][iOscChan][iSyst][iMode].size(); iVar1++)
494 { // Loop over first dimension
495 for (unsigned int iVar2 = 0; iVar2 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1].size(); iVar2++)
496 { // Loop over second dimension
497 for (unsigned int iVar3 = 0; iVar3 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2].size(); iVar3++)
498 { // Loop over third dimension
499 int splineindex=indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2][iVar3];
500 if (splinevec_Monolith[splineindex])
501 {
502 UniqueSystSplines.push_back(splinevec_Monolith[splineindex]);
503 UniqueSystIndices.push_back(GlobalSystIndex[iSample][iSyst]);
504 FoundNonFlatSpline = true;
505 }
506 if (FoundNonFlatSpline) { break;}
507 }//3D loop end
508 if (FoundNonFlatSpline) { break;}
509 }//2D loop end
510 if (FoundNonFlatSpline){ break; }
511 }//1D loop end
512 if (FoundNonFlatSpline){ break; }
513 }//mode loop end
514 if (FoundNonFlatSpline) { break; }
515 }//osc loop end
516 //ETA - only push back unique name if a non-flat response has been found
517 if(FoundNonFlatSpline){
518 UniqueSystNames.push_back(SystName);
519 }
520
521 if (!FoundNonFlatSpline)
522 {
523 MACH3LOG_INFO("{} syst has no response in sample {}", SystName, iSample);
524 MACH3LOG_INFO("Whilst this isn't necessarily a problem, it seems odd");
525 continue;
526 }
527 }
528 }//Syst loop end
529 }
530
531 nParams = static_cast<short int>(UniqueSystSplines.size());
532
533 // DB Find the number of splines knots which assumes each instance of the syst has the same number of knots
534 SplineSegments = new short int[nParams]();
535 ParamValues = new float[nParams]();
536 SplineInfoArray.resize(nParams);
537 for (int iSpline = 0; iSpline < nParams; iSpline++)
538 {
539 SplineInfoArray[iSpline].nPts = static_cast<M3::int_t>(UniqueSystSplines[iSpline]->GetNp());
540 SplineInfoArray[iSpline].xPts.resize(SplineInfoArray[iSpline].nPts);
541 SplineInfoArray[iSpline].splineParsPointer = xsec->RetPointer(UniqueSystIndices[iSpline]);
542 for (int iKnot = 0; iKnot < SplineInfoArray[iSpline].nPts; iKnot++)
543 {
544 M3::float_t xPoint;
545 M3::float_t yPoint;
546 UniqueSystSplines[iSpline]->GetKnot(iKnot, xPoint, yPoint);
547 SplineInfoArray[iSpline].xPts[iKnot] = xPoint;
548 }
549 //ETA - let this just be set as the first segment by default
550 SplineSegments[iSpline] = 0;
551 ParamValues[iSpline] = 0.;
552 }
553
554 MACH3LOG_INFO("nUniqueSysts: {}", nParams);
555 MACH3LOG_INFO("{:<15} | {:<20} | {:<6}", "Spline Index", "Syst Name", "nKnots");
556 for (int iUniqueSyst = 0; iUniqueSyst < nParams; iUniqueSyst++)
557 {
558 MACH3LOG_INFO("{:<15} | {:<20} | {:<6}", iUniqueSyst, UniqueSystNames[iUniqueSyst], SplineInfoArray[iUniqueSyst].nPts);
559 }
560
561 //ETA
562 //Isn't this just doing what CountNumberOfLoadedSplines() does?
563 int nCombinations_FlatSplines = 0;
564 int nCombinations_All = 0;
565 // DB Now actually loop over splines to determine which are all null i.e. flat
566 for (unsigned int iSample = 0; iSample < indexvec.size(); iSample++)
567 { // Loop over systematics
568 for (unsigned int iOscChan = 0; iOscChan < indexvec[iSample].size(); iOscChan++)
569 { // Loop over oscillation channels
570 for (unsigned int iSyst = 0; iSyst < indexvec[iSample][iOscChan].size(); iSyst++)
571 { // Loop over systematics
572 for (unsigned int iMode = 0; iMode < indexvec[iSample][iOscChan][iSyst].size(); iMode++)
573 { // Loop over modes
574 //bool isFlat = false;
575 for (unsigned int iVar1 = 0; iVar1 < indexvec[iSample][iOscChan][iSyst][iMode].size(); iVar1++)
576 { // Loop over first dimension
577 for (unsigned int iVar2 = 0; iVar2 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1].size(); iVar2++)
578 { // Loop over second dimension
579 for (unsigned int iVar3 = 0; iVar3 < indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2].size(); iVar3++)
580 { // Loop over third dimension
581 int splineindex=indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2][iVar3];
582 if (splinevec_Monolith[splineindex])
583 {
584 nCombinations_All += 1;
585 } else{
586 nCombinations_FlatSplines += 1;
587 nCombinations_All += 1;
588 }
589 }
590 }
591 }
592 }
593 }
594 }
595 }
596 // We need to grab the maximum number of knots
597 MACH3LOG_INFO("Number of combinations of Sample, OscChan, Syst and Mode which have entirely flat response: {} / {}", nCombinations_FlatSplines, nCombinations_All);
598}
599
600//****************************************
601// Rather work with spline coefficients in the splines, let's copy ND and use coefficient arrays
602void BinnedSplineHandler::getSplineCoeff_SepMany(int splineindex, M3::float_t* &xArray, M3::float_t* &manyArray){
603//****************************************
604 //No point evaluating a flat spline
605 int nPoints = splinevec_Monolith[splineindex]->GetNp();
606
607 for (int i = 0; i < nPoints; i++) {
608 xArray[i] = 1.0;
609 for (int j = 0; j < 4; j++) {
610 manyArray[i*4+j] = 1.0;
611 }
612 }
613
614 for(int i=0; i<nPoints; i++) {
615 M3::float_t x = M3::float_t(-999.99);
616 M3::float_t y = M3::float_t(-999.99);
617 M3::float_t b = M3::float_t(-999.99);
618 M3::float_t c = M3::float_t(-999.99);
619 M3::float_t d = M3::float_t(-999.99);
620 splinevec_Monolith[splineindex]->GetCoeff(i, x, y, b, c, d);
621
622 // Store the coefficients for each knot contiguously in memory
623 // 4 because manyArray stores y,b,c,d
624 xArray[i] = x;
625 manyArray[i*4] = y;
626 manyArray[i*4+1] = b;
627 manyArray[i*4+2] = c;
628 manyArray[i*4+3] = d;
629 }
630
631 //We now clean up the splines!
632 delete splinevec_Monolith[splineindex];
633 splinevec_Monolith[splineindex] = nullptr;
634}
635
636//****************************************
637//Equally though could just use KinematicVariable to map back
638std::string BinnedSplineHandler::getDimLabel(const int iSample, const unsigned int Axis) const
639//****************************************
640{
641 if(Axis > DimensionLabels[iSample].size()){
642 MACH3LOG_ERROR("The spline Axis you are trying to get the label of is larger than the number of dimensions");
643 MACH3LOG_ERROR("You are trying to get axis {} but have only got {}", Axis, Dimensions[iSample]);
644 throw MaCh3Exception(__FILE__ , __LINE__ );
645 }
646 return DimensionLabels.at(iSample).at(Axis);
647}
648
649//****************************************
650//Returns sample index in
651int BinnedSplineHandler::getSampleIndex(const std::string& SampleName) const{
652//****************************************
653 for (size_t iSample = 0; iSample < SampleNames.size(); ++iSample) {
654 if (SampleName == SampleNames[iSample]) {
655 return static_cast<int>(iSample);
656 }
657 }
658 MACH3LOG_ERROR("Sample name not found: {}", SampleName);
659 throw MaCh3Exception(__FILE__, __LINE__);
660}
661
662//****************************************
663void BinnedSplineHandler::PrintSampleDetails(const std::string& SampleName) const
664//****************************************
665{
666 const int iSample = getSampleIndex(SampleName);
667
668 MACH3LOG_INFO("Details about sample: {:<20}", SampleNames[iSample]);
669 MACH3LOG_INFO("\t Dimension: {:<35}", Dimensions[iSample]);
670 MACH3LOG_INFO("\t nSplineParam: {:<35}", nSplineParams[iSample]);
671 MACH3LOG_INFO("\t nOscChan: {:<35}", nOscChans[iSample]);
672}
673
674//****************************************
675void BinnedSplineHandler::PrintArrayDetails(const std::string& SampleName) const
676//****************************************
677{
678 int iSample = getSampleIndex(SampleName);
679 int nOscChannels = int(indexvec[iSample].size());
680 MACH3LOG_INFO("Sample {} has {} oscillation channels", SampleName, nOscChannels);
681
682 for (int iOscChan = 0; iOscChan < nOscChannels; iOscChan++)
683 {
684 int nSysts = int(indexvec[iSample][iOscChan].size());
685 MACH3LOG_INFO("Oscillation channel {} has {} systematics", iOscChan, nSysts);
686 }
687
688 MACH3LOG_INFO("#----------------------------------------------------------------------------------------------------------------------------------#");
689 MACH3LOG_INFO("Printing no. of modes affected by each systematic for each oscillation channel");
690 for (unsigned int iOscChan = 0; iOscChan < indexvec[iSample].size(); iOscChan++) {
691 std::string modes = fmt::format("OscChan: {}\t", iOscChan);
692 for (unsigned int iSyst = 0; iSyst < indexvec[iSample][iOscChan].size(); iSyst++) {
693 modes += fmt::format("{} ", indexvec[iSample][iOscChan][iSyst].size());
694 }
695 MACH3LOG_INFO("{}", modes);
696 }
697 MACH3LOG_INFO("#----------------------------------------------------------------------------------------------------------------------------------#");
698}
699
700//****************************************
701bool BinnedSplineHandler::isValidSplineIndex(const std::string& SampleName, int iOscChan, int iSyst, int iMode, int iVar1, int iVar2, int iVar3)
702//****************************************
703{
704 int iSample = getSampleIndex(SampleName);
705 bool isValid = true;
706
707 // Lambda to check if an index is valid for a specific dimension
708 auto checkIndex = [&isValid](int index, size_t size, const std::string& name) {
709 if (index < 0 || index >= int(size)) {
710 MACH3LOG_ERROR("{} index is invalid! 0 <= Index < {} ", name, size);
711 isValid = false;
712 }
713 };
714
715 checkIndex(iSample, indexvec.size(), "Sample");
716 if (isValid) checkIndex(iOscChan, indexvec[iSample].size(), "OscChan");
717 if (isValid) checkIndex(iSyst, indexvec[iSample][iOscChan].size(), "Syst");
718 if (isValid) checkIndex(iMode, indexvec[iSample][iOscChan][iSyst].size(), "Mode");
719 if (isValid) checkIndex(iVar1, indexvec[iSample][iOscChan][iSyst][iMode].size(), "Var1");
720 if (isValid) checkIndex(iVar2, indexvec[iSample][iOscChan][iSyst][iMode][iVar1].size(), "Var2");
721 if (isValid) checkIndex(iVar3, indexvec[iSample][iOscChan][iSyst][iMode][iVar1][iVar2].size(), "Var3");
722
723 if (!isValid)
724 {
725 MACH3LOG_ERROR("Given iSample: {}", iSample);
726 MACH3LOG_ERROR("Given iOscChan: {}", iOscChan);
727 MACH3LOG_ERROR("Given iSyst: {}", iSyst);
728 MACH3LOG_ERROR("Given iMode: {}", iMode);
729 MACH3LOG_ERROR("Given iVar1: {}", iVar1);
730 MACH3LOG_ERROR("Given iVar2: {}", iVar2);
731 MACH3LOG_ERROR("Given iVar3: {}", iVar3);
732 MACH3LOG_ERROR("Come visit me at : {} : {}", __FILE__, __LINE__);
733 throw MaCh3Exception(__FILE__, __LINE__);
734 }
735
736 return true;
737}
738
739//****************************************
741//****************************************
742{
743 const int NBins = Axis->GetNbins();
744 std::string text = "";
745 for (int iBin = 0; iBin <= NBins; iBin++) {
746 text += fmt::format("{} ", Axis->GetXbins()->GetAt(iBin));
747 }
748 MACH3LOG_INFO("{}", text);
749}
750
751//****************************************
752std::vector< std::vector<int> > BinnedSplineHandler::GetEventSplines(const std::string& SampleName, int iOscChan, int EventMode, double Var1Val, double Var2Val, double Var3Val)
753//****************************************
754{
755 std::vector<std::vector<int>> ReturnVec;
756 int SampleIndex = -1;
757 for (unsigned int iSample = 0; iSample < SampleNames.size(); iSample++) {
758 if (SampleName == SampleNames[iSample]) {
759 SampleIndex = iSample;
760 }
761 }
762
763 if (SampleIndex == -1) {
764 MACH3LOG_ERROR("Sample not found: {}", SampleName);
765 throw MaCh3Exception(__FILE__, __LINE__);
766 }
767
768 int nSplineSysts = static_cast<int>(indexvec[SampleIndex][iOscChan].size());
769
770
771 int Mode = -1;
772 std::string SuffixForEventMode = Modes->GetSplineSuffixFromMaCh3Mode(EventMode);
773 for (int iMode=0;iMode<Modes->GetNModes();iMode++) {
774 if (SuffixForEventMode == Modes->GetSplineSuffixFromMaCh3Mode(iMode)) {
775 Mode = iMode;
776 break;
777 }
778 }
779 if (Mode == -1) {
780 return ReturnVec;
781 }
782
783 int Var1Bin = SplineBinning[SampleIndex][iOscChan][0]->FindBin(Var1Val)-1;
784 if (Var1Bin < 0 || Var1Bin >= SplineBinning[SampleIndex][iOscChan][0]->GetNbins()) {
785 return ReturnVec;
786 }
787
788 int Var2Bin = SplineBinning[SampleIndex][iOscChan][1]->FindBin(Var2Val)-1;
789 if (Var2Bin < 0 || Var2Bin >= SplineBinning[SampleIndex][iOscChan][1]->GetNbins()) {
790 return ReturnVec;
791 }
792
793 int Var3Bin = SplineBinning[SampleIndex][iOscChan][2]->FindBin(Var3Val)-1;
794 if (Var3Bin < 0 || Var3Bin >= SplineBinning[SampleIndex][iOscChan][2]->GetNbins()){
795 return ReturnVec;
796 }
797
798 for(int iSyst=0; iSyst<nSplineSysts; iSyst++){
799 std::vector<int> spline_modes = SplineModeVecs[SampleIndex][iSyst];
800 int nSampleModes = static_cast<int>(spline_modes.size());
801
802 //ETA - look here at the length of spline_modes and what you're actually comparing against
803 for(int iMode = 0; iMode<nSampleModes ; iMode++){
804 //Only consider if the event mode (Mode) matches ones of the spline modes
805 if (Mode == spline_modes[iMode]) {
806 int splineID=indexvec[SampleIndex][iOscChan][iSyst][iMode][Var1Bin][Var2Bin][Var3Bin];
807 //Also check that the spline isn't flat
808 if(!isflatarray[splineID]){
809 ReturnVec.push_back({SampleIndex, iOscChan, iSyst, iMode, Var1Bin, Var2Bin, Var3Bin});
810 }
811 }
812 }
813 }
814
815 return ReturnVec;
816}
817
818// checks if there are multiple modes with the same SplineSuffix
819// (for example if CCRES and CCCoherent are treated as one spline mode)
820std::vector< std::vector<int> > BinnedSplineHandler::StripDuplicatedModes(const std::vector< std::vector<int> >& InputVector) {
821
822 //ETA - this is of size nPars from the xsec model
823 size_t InputVectorSize = InputVector.size();
824 std::vector< std::vector<int> > ReturnVec(InputVectorSize);
825
826 //ETA - loop over all systematics
827 for (size_t iSyst=0;iSyst<InputVectorSize;iSyst++) {
828 std::vector<int> TmpVec;
829 std::vector<std::string> TestVec;
830
831 //Loop over the modes that we've listed in xsec cov
832 for (unsigned int iMode = 0 ; iMode < InputVector[iSyst].size() ; iMode++) {
833 int Mode = InputVector[iSyst][iMode];
834 std::string ModeName = Modes->GetSplineSuffixFromMaCh3Mode(Mode);
835
836 bool IncludeMode = true;
837 for (auto TestString : TestVec) {
838 if (ModeName == TestString) {
839 IncludeMode = false;
840 break;
841 }
842 }
843
844 if (IncludeMode) {
845 TmpVec.push_back(Mode);
846 TestVec.push_back(ModeName);
847 }
848 }
849
850 ReturnVec[iSyst] = TmpVec;
851 }
852 return ReturnVec;
853}
854
855void BinnedSplineHandler::FillSampleArray(std::string SampleName, std::vector<std::string> OscChanFileNames)
856{
857 int iSample = getSampleIndex(SampleName);
858 int nOscChannels = nOscChans[iSample];
859
860 for (int iOscChan = 0; iOscChan < nOscChannels; iOscChan++) {
861 MACH3LOG_INFO("Processing: {}", OscChanFileNames[iOscChan]);
862 TSpline3* mySpline = nullptr;
863 TSpline3_red* Spline = nullptr;
864 TString Syst, Mode;
865 int nKnots, SystNum, ModeNum, Var1Bin, Var2Bin, Var3Bin = M3::_BAD_INT_;
866 double x,y = M3::_BAD_DOUBLE_;
867 bool isFlat = true;
868
869 std::set<std::string> SplineFileNames;
870
871 auto File = std::unique_ptr<TFile>(TFile::Open(OscChanFileNames[iOscChan].c_str()));
872
873 if (!File || File->IsZombie()) {
874 MACH3LOG_ERROR("File {} not found", OscChanFileNames[iOscChan]);
875 throw MaCh3Exception(__FILE__, __LINE__);
876 }
877
878 //This is the MC specific part of the code
879 //i.e. we always assume that the splines are just store in single TDirectory and they're all in there as single objects
880 for (auto k : *File->GetListOfKeys()) {
881 auto Key = static_cast<TKey*>(k);
882 TClass *Class = gROOT->GetClass(Key->GetClassName(), false);
883 if(!Class->InheritsFrom("TSpline3")) {
884 continue;
885 }
886
887 std::string FullSplineName = std::string(Key->GetName());
888
889 if (SplineFileNames.count(FullSplineName) > 0) {
890 MACH3LOG_CRITICAL("Skipping spline - Found a spline whose name has already been encountered before: {}", FullSplineName);
891 continue;
892 }
893 SplineFileNames.insert(FullSplineName);
894
895 std::vector<std::string> Tokens = GetTokensFromSplineName(FullSplineName);
896
897 if (Tokens.size() != kNTokens) {
898 MACH3LOG_ERROR("Invalid tokens from spline name - Expected {} tokens. Check implementation in GetTokensFromSplineName()", static_cast<int>(kNTokens));
899 throw MaCh3Exception(__FILE__, __LINE__);
900 }
901
902 Syst = Tokens[kSystToken];
903 Mode = Tokens[kModeToken];
904 Var1Bin = std::stoi(Tokens[kVar1BinToken]);
905 Var2Bin = std::stoi(Tokens[kVar2BinToken]);
906 Var3Bin = std::stoi(Tokens[kVar3BinToken]);
907
908 SystNum = -1;
909 for (unsigned iSyst = 0; iSyst < SplineFileParPrefixNames[iSample].size(); iSyst++) {
910 if (Syst == SplineFileParPrefixNames[iSample][iSyst]) {
911 SystNum = iSyst;
912 break;
913 }
914 }
915
916 // If the syst doesn't match any of the spline names then skip it
917 if (SystNum == -1){
918 MACH3LOG_DEBUG("Couldn't match!!");
919 MACH3LOG_DEBUG("Couldn't Match any systematic name in ParameterHandler with spline name: {}" , FullSplineName);
920 continue;
921 }
922
923 ModeNum = -1;
924 for (unsigned int iMode = 0; iMode < SplineModeVecs[iSample][SystNum].size(); iMode++) {
925 if (Mode == Modes->GetSplineSuffixFromMaCh3Mode(SplineModeVecs[iSample][SystNum][iMode])) {
926 ModeNum = iMode;
927 break;
928 }
929 }
930
931 if (ModeNum == -1) {
932 //DB - If you have splines in the root file that you don't want to use (e.g. removing a mode from a syst), this will cause a throw
933 // Therefore include as debug warning and continue instead
934 MACH3LOG_DEBUG("Couldn't find mode for {} in {}. Problem Spline is : {} ", Mode, Syst, FullSplineName);
935 continue;
936 }
937
938 mySpline = Key->ReadObject<TSpline3>();
939
940 if (isValidSplineIndex(SampleName, iOscChan, SystNum, ModeNum, Var1Bin, Var2Bin, Var3Bin)) { // loop over all the spline knots and check their value
941 MACH3LOG_DEBUG("Pushed back monolith for spline {}", FullSplineName);
942 // if the value is 1 then set the flat bool to false
943 nKnots = mySpline->GetNp();
944 isFlat = true;
945 for (int iKnot = 0; iKnot < nKnots; iKnot++) {
946 mySpline->GetKnot(iKnot, x, y);
947 if (y < 0.99999 || y > 1.00001)
948 {
949 isFlat = false;
950 break;
951 }
952 }
953
954 //Rather than keeping a mega vector of splines then converting, this should just keep everything nice in memory!
955 indexvec[iSample][iOscChan][SystNum][ModeNum][Var1Bin][Var2Bin][Var3Bin]=MonolithIndex;
956 coeffindexvec.push_back(CoeffIndex);
957 // Should save memory rather saving [x_i_0 ,... x_i_maxknots] for every spline!
958 if (isFlat) {
959 splinevec_Monolith.push_back(nullptr);
960 delete mySpline;
961 } else {
962 ApplyKnotWeightCapTSpline3(mySpline, SystNum, xsec);
963 Spline = new TSpline3_red(mySpline, SplineInterpolationTypes[iSample][SystNum]);
964 if(mySpline) delete mySpline;
965
966 splinevec_Monolith.push_back(Spline);
967 uniquecoeffindices.push_back(MonolithIndex); //So we can get the unique coefficients and skip flat splines later on!
968 CoeffIndex+=nKnots;
969 }
970 //Incrementing MonolithIndex to keep track of number of valid spline indices
971 MonolithIndex+=1;
972 } else {
973 //Potentially you are not a valid spline index
974 delete mySpline;
975 }
976 }//End of loop over all TKeys in file
977
978 //A bit of clean up
979 File->Delete("*");
980 File->Close();
981 } //End of oscillation channel loop
982}
#define _MaCh3_Safe_Include_Start_
KS: Avoiding warning checking for headers.
Definition: Core.h:106
#define _MaCh3_Safe_Include_End_
KS: Restore warning checking after including external headers.
Definition: Core.h:117
int size
std::vector< bool > isFlat
#define MACH3LOG_CRITICAL
Definition: MaCh3Logger.h:26
#define MACH3LOG_DEBUG
Definition: MaCh3Logger.h:22
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:25
#define MACH3LOG_INFO
Definition: MaCh3Logger.h:23
void CleanVector(std::vector< T > &vec)
Generic cleanup function.
void CleanContainer(std::vector< T * > &container)
Generic cleanup function.
@ kCoeffB
Coefficient B.
Definition: SplineCommon.h:21
@ kCoeffD
Coefficient D.
Definition: SplineCommon.h:23
@ kCoeffY
Coefficient Y.
Definition: SplineCommon.h:20
@ kCoeffC
Coefficient C.
Definition: SplineCommon.h:22
void ApplyKnotWeightCapTSpline3(TSpline3 *&Spline, const int splineParsIndex, ParameterHandlerGeneric *XsecCov)
EM: Apply capping to knot weight for specified spline parameter. param graph needs to have been set i...
void getSplineCoeff_SepMany(int splineindex, M3::float_t *&xArray, M3::float_t *&manyArray)
std::string getDimLabel(const int BinningOpt, const unsigned int Axis) const
void Evaluate()
CW: This Eval should be used when using two separate x,{y,a,b,c,d} arrays to store the weights; proba...
bool * isflatarray
Need to keep track of which splines are flat and which aren't.
int CountNumberOfLoadedSplines(bool NonFlat=false, int Verbosity=0)
std::vector< std::vector< std::vector< std::vector< std::vector< std::vector< std::vector< int > > > > > > > indexvec
Variables related to determined which modes have splines and which piggy-back of other modes.
std::vector< std::vector< std::string > > SplineFileParPrefixNames
void CalcSplineWeights() override
CPU based code which eval weight for each spline.
std::vector< std::vector< int > > GetEventSplines(const std::string &SampleName, int iOscChan, int EventMode, double Var1Val, double Var2Val, double Var3Val)
Return the splines which affect a given event.
std::vector< int > coeffindexvec
void cleanUpMemory()
Remove setup variables not needed for spline evaluations.
void PrintBinning(TAxis *Axis) const
ParameterHandlerGeneric * xsec
Pointer to covariance from which we get information about spline params.
virtual std::vector< std::string > GetTokensFromSplineName(std::string FullSplineName)=0
std::vector< TAxis * > FindSplineBinning(const std::string &FileName, const std::string &SampleName)
Grab histograms with spline binning.
std::vector< int > nSplineParams
void ModifyWeights() override
Calc total event weight, not used by Bin-by-bin splines.
std::vector< std::vector< std::vector< int > > > SplineModeVecs
std::vector< int > nOscChans
M3::float_t * xcoeff_arr
x coefficients for each spline
std::vector< std::string > UniqueSystNames
name of each spline parameter
int getSampleIndex(const std::string &SampleName) const
Get index of sample based on name.
void AddSample(const std::string &SampleName, const std::vector< std::string > &OscChanFileNames, const std::vector< std::string > &SplineVarNames)
add oscillation channel to spline monolith
std::vector< std::vector< SplineInterpolation > > SplineInterpolationTypes
spline interpolation types for each sample. These vectors are from a call to GetSplineInterpolationFr...
std::vector< int > Dimensions
std::vector< int > uniquesplinevec_Monolith
Maps single spline object with single parameter.
void PrintArrayDetails(const std::string &SampleName) const
std::vector< std::vector< int > > StripDuplicatedModes(const std::vector< std::vector< int > > &InputVector)
Check if there are any repeated modes. This is used to reduce the number of modes in case many intera...
std::vector< M3::float_t > weightvec_Monolith
Stores weight from spline evaluation for each single spline.
void BuildSampleIndexingArray(const std::string &SampleName)
std::vector< std::vector< std::string > > DimensionLabels
BinnedSplineHandler(ParameterHandlerGeneric *xsec_, MaCh3Modes *Modes_)
Constructor.
std::vector< std::vector< std::vector< TAxis * > > > SplineBinning
virtual void FillSampleArray(std::string SampleName, std::vector< std::string > OscChanFileNames)
Loads and processes splines from ROOT files for a given sample.
std::vector< TSpline3_red * > splinevec_Monolith
holds each spline object before stripping into coefficient monolith
std::vector< int > uniquecoeffindices
Unique coefficient indices.
MaCh3Modes * Modes
pointer to MaCh3 Mode from which we get spline suffix
std::vector< std::vector< int > > GlobalSystIndex
This holds the global spline index and is used to grab the current parameter value to evaluate spline...
void PrintSampleDetails(const std::string &SampleName) const
Print info like Sample ID of spline params etc.
virtual ~BinnedSplineHandler()
Destructor.
std::vector< std::string > SampleNames
void TransferToMonolith()
flatten multidimensional spline array into proper monolith
std::vector< int > UniqueSystIndices
Global index of each spline param, it allows us to match spline ordering with global.
bool isValidSplineIndex(const std::string &SampleName, int iSyst, int iOscChan, int iMode, int iVar1, int iVar2, int iVar3)
Ensure we have spline for a given bin.
M3::float_t * manycoeff_arr
ybcd coefficients for each spline
Custom exception class for MaCh3 errors.
KS: Class describing MaCh3 modes used in the analysis, it is being initialised from config.
Definition: MaCh3Modes.h:41
int GetNModes() const
KS: Get number of modes, keep in mind actual number is +1 greater due to unknown category.
Definition: MaCh3Modes.h:52
std::string GetSplineSuffixFromMaCh3Mode(const int Index)
DB: Get binned spline mode suffic from MaCh3 Mode.
Definition: MaCh3Modes.cpp:226
const double * RetPointer(const int iParam)
DB Pointer return to param position.
Class responsible for handling of systematic error parameters with different types defined in the con...
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
CW: Reduced TSpline3 class.
const std::vector< int > GetGlobalSystIndexFromSampleName(const std::string &SampleName, const SystType Type)
DB Get spline parameters depending on given SampleName.
int GetNumParamsFromSampleName(const std::string &SampleName, const SystType Type)
DB Grab the number of parameters for the relevant SampleName.
const std::vector< std::string > GetSplineParsNamesFromSampleName(const std::string &SampleName)
DB Get spline parameters depending on given SampleName.
const std::vector< SplineInterpolation > GetSplineInterpolationFromSampleName(const std::string &SampleName)
Get the interpolation types for splines affecting a particular SampleName.
const std::vector< std::vector< int > > GetSplineModeVecFromSampleName(const std::string &SampleName)
DB Grab the Spline Modes for the relevant SampleName.
double float_t
Definition: Core.h:28
static constexpr const double _BAD_DOUBLE_
Default value used for double initialisation.
Definition: Core.h:43
constexpr T fmaf_t(T x, T y, T z)
Function template for fused multiply-add.
Definition: Core.h:35
int int_t
Definition: Core.h:29
static constexpr const int _BAD_INT_
Default value used for int initialisation.
Definition: Core.h:45