MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
ParameterHandlerGeneric.cpp
Go to the documentation of this file.
2
3// ********************************************
4// ETA - YAML constructor
5// this will replace the root file constructor but let's keep it in
6// to do some validations
7ParameterHandlerGeneric::ParameterHandlerGeneric(const std::vector<std::string>& YAMLFile, std::string name, double threshold, int FirstPCA, int LastPCA)
8 : ParameterHandlerBase(YAMLFile, name, threshold, FirstPCA, LastPCA){
9// ********************************************
11
12 //ETA - again this really doesn't need to be hear...
13 for (int i = 0; i < _fNumPar; i++)
14 {
15 // Sort out the print length
16 if(int(_fNames[i].length()) > PrintLength) PrintLength = int(_fNames[i].length());
17 } // end the for loop
18
19 MACH3LOG_DEBUG("Constructing instance of ParameterHandler");
20 InitParams();
21 // Print
22 Print();
23}
24
25// ********************************************
27// ********************************************
28 _fSystToGlobalSystIndexMap.resize(SystType::kSystTypes);
29
30 _fParamType = std::vector<SystType>(_fNumPar);
31 _ParameterGroup = std::vector<std::string>(_fNumPar);
32
33 //KS: We know at most how params we expect so reserve memory for max possible params. Later we will shrink to size to not waste memory. Reserving means slightly faster loading and possible less memory fragmentation.
34 NormParams.reserve(_fNumPar);
35 SplineParams.reserve(_fNumPar);
36 FuncParams.reserve(_fNumPar);
37 OscParams.reserve(_fNumPar);
38
39 int i = 0;
40 unsigned int ParamCounter[SystType::kSystTypes] = {0};
41 //ETA - read in the systematics. Would be good to add in some checks to make sure
42 //that there are the correct number of entries i.e. are the _fNumPars for Names,
43 //PreFitValues etc etc.
44 for (auto const &param : _fYAMLDoc["Systematics"])
45 {
46 _ParameterGroup[i] = Get<std::string>(param["Systematic"]["ParameterGroup"], __FILE__ , __LINE__);
47
48 //Fill the map to get the correlations later as well
49 auto ParamType = Get<std::string>(param["Systematic"]["Type"], __FILE__ , __LINE__);
50 //Now load in variables for spline systematics only
51 if (ParamType.find(SystType_ToString(SystType::kSpline)) != std::string::npos)
52 {
53 //Set param type
54 _fParamType[i] = SystType::kSpline;
55 // Fill Spline info
56 SplineParams.push_back(GetSplineParameter(param["Systematic"], i));
57
58 if (param["Systematic"]["SplineInformation"]["SplineName"]) {
59 _fSplineNames.push_back(param["Systematic"]["SplineInformation"]["SplineName"].as<std::string>());
60 }
61
62 //Insert the mapping from the spline index i.e. the length of _fSplineNames etc
63 //to the Systematic index i.e. the counter for things like _fSampleID
64 _fSystToGlobalSystIndexMap[SystType::kSpline].insert(std::make_pair(ParamCounter[SystType::kSpline], i));
65 ParamCounter[SystType::kSpline]++;
66 } else if(param["Systematic"]["Type"].as<std::string>() == SystType_ToString(SystType::kNorm)) {
67 _fParamType[i] = SystType::kNorm;
68 NormParams.push_back(GetNormParameter(param["Systematic"], i));
69 _fSystToGlobalSystIndexMap[SystType::kNorm].insert(std::make_pair(ParamCounter[SystType::kNorm], i));
70 ParamCounter[SystType::kNorm]++;
71 } else if(param["Systematic"]["Type"].as<std::string>() == SystType_ToString(SystType::kFunc)){
72 _fParamType[i] = SystType::kFunc;
73 FuncParams.push_back(GetFunctionalParameters(param["Systematic"], i));
74 _fSystToGlobalSystIndexMap[SystType::kFunc].insert(std::make_pair(ParamCounter[SystType::kFunc], i));
75 ParamCounter[SystType::kFunc]++;
76 } else if(param["Systematic"]["Type"].as<std::string>() == SystType_ToString(SystType::kOsc)){
77 _fParamType[i] = SystType::kOsc;
78 OscParams.push_back(GetOscillationParameters(param["Systematic"], i));
79 _fSystToGlobalSystIndexMap[SystType::kOsc].insert(std::make_pair(ParamCounter[SystType::kOsc], i));
80 ParamCounter[SystType::kOsc]++;
81 } else{
82 MACH3LOG_ERROR("Given unrecognised systematic type: {}", param["Systematic"]["Type"].as<std::string>());
83 std::string expectedTypes = "Expecting ";
84 for (int s = 0; s < SystType::kSystTypes; ++s) {
85 if (s > 0) expectedTypes += ", ";
86 expectedTypes += SystType_ToString(static_cast<SystType>(s)) + "\"";
87 }
88 expectedTypes += ".";
89 MACH3LOG_ERROR(expectedTypes);
90 throw MaCh3Exception(__FILE__, __LINE__);
91 }
92 i++;
93 } //end loop over params
94
95 //Add a sanity check,
96 if(_fSplineNames.size() != ParamCounter[SystType::kSpline]){
97 MACH3LOG_ERROR("_fSplineNames is of size {} but found {} spline parameters", _fSplineNames.size(), ParamCounter[SystType::kSpline]);
98 throw MaCh3Exception(__FILE__, __LINE__);
99 }
100 //KS We resized them above to all params to fight memory fragmentation, now let's resize to fit only allocated memory to save RAM
101 NormParams.shrink_to_fit();
102 SplineParams.shrink_to_fit();
103 FuncParams.shrink_to_fit();
104 OscParams.shrink_to_fit();
105}
106
107// ********************************************
109// ********************************************
110 MACH3LOG_DEBUG("Deleting ParameterHandler");
111}
112
113// ********************************************
114// DB Grab the Spline Names for the relevant SampleName
115const std::vector<std::string> ParameterHandlerGeneric::GetSplineParsNamesFromSampleName(const std::string& SampleName) {
116// ********************************************
117 std::vector<std::string> returnVec;
118 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kSpline]) {
119 auto &SplineIndex = pair.first;
120 auto &SystIndex = pair.second;
121 if (AppliesToSample(SystIndex, SampleName)) { //If parameter applies to required Sample
122 returnVec.push_back(_fSplineNames.at(SplineIndex));
123 }
124 }
125 return returnVec;
126}
127
128// ********************************************
129const std::vector<SplineInterpolation> ParameterHandlerGeneric::GetSplineInterpolationFromSampleName(const std::string& SampleName) {
130// ********************************************
131 std::vector<SplineInterpolation> returnVec;
132 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kSpline]) {
133 auto &SplineIndex = pair.first;
134 auto &SystIndex = pair.second;
135
136 if (AppliesToSample(SystIndex, SampleName)) { //If parameter applies to required SampleID
137 returnVec.push_back(SplineParams.at(SplineIndex)._SplineInterpolationType);
138 }
139 }
140 return returnVec;
141}
142
143// ********************************************
144// DB Grab the Spline Modes for the relevant SampleName
145const std::vector< std::vector<int> > ParameterHandlerGeneric::GetSplineModeVecFromSampleName(const std::string& SampleName) {
146// ********************************************
147 std::vector< std::vector<int> > returnVec;
148 //Need a counter or something to correctly get the index in _fSplineModes since it's not of length nPars
149 //Should probably just make a std::map<std::string, int> for param name to FD spline index
150 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kSpline]) {
151 auto &SplineIndex = pair.first;
152 auto &SystIndex = pair.second;
153 if (AppliesToSample(SystIndex, SampleName)) { //If parameter applies to required SampleID
154 returnVec.push_back(SplineParams.at(SplineIndex)._fSplineModes);
155 }
156 }
157 return returnVec;
158}
159
160// ********************************************
161// Get Norm params
162NormParameter ParameterHandlerGeneric::GetNormParameter(const YAML::Node& param, const int Index) {
163// ********************************************
164 NormParameter norm;
165
166 GetBaseParameter(param, Index, norm);
167
170 norm.modes = GetFromManager<std::vector<int>>(param["Mode"], {}, __FILE__ , __LINE__);
171 norm.pdgs = GetFromManager<std::vector<int>>(param["NeutrinoFlavour"], {}, __FILE__ , __LINE__);
172 norm.preoscpdgs = GetFromManager<std::vector<int>>(param["NeutrinoFlavourUnosc"], {}, __FILE__ , __LINE__);
173 norm.targets = GetFromManager<std::vector<int>>(param["TargetNuclei"], {}, __FILE__ , __LINE__);
174
175 if(_fLowBound[Index] < 0.) {
176 MACH3LOG_ERROR("Normalisation Parameter {} ({}), has lower parameters bound which can go below 0 and is equal {}",
177 GetParFancyName(Index), Index, _fLowBound[Index]);
178 MACH3LOG_ERROR("Normalisation parameters can't go bellow 0 as this is unphysical");
179 throw MaCh3Exception(__FILE__, __LINE__);
180 }
181 int NumKinematicCuts = 0;
182 if(param["KinematicCuts"]) {
183 NumKinematicCuts = int(param["KinematicCuts"].size());
184
185 std::vector<std::string> TempKinematicStrings;
186 std::vector<std::vector<std::vector<double>>> TempKinematicBounds;
187 //First element of TempKinematicBounds is always -999, and size is then 3
188 for(int KinVar_i = 0 ; KinVar_i < NumKinematicCuts ; ++KinVar_i) {
189 //ETA: This is a bit messy, Kinematic cuts is a list of maps
190 for (YAML::const_iterator it = param["KinematicCuts"][KinVar_i].begin();it!=param["KinematicCuts"][KinVar_i].end();++it) {
191 TempKinematicStrings.push_back(it->first.as<std::string>());
192 TempKinematicBounds.push_back(Get2DBounds(it->second));
193 }
194 if(TempKinematicStrings.size() == 0) {
195 MACH3LOG_ERROR("Received a KinematicCuts node but couldn't read the contents (it's a list of single-element dictionaries (python) = map of pairs (C++))");
196 MACH3LOG_ERROR("For Param {}", norm.name);
197 throw MaCh3Exception(__FILE__, __LINE__);
198 }
199 }//KinVar_i
200 norm.KinematicVarStr = TempKinematicStrings;
201 norm.Selection = TempKinematicBounds;
202 }
203
204 //Next ones are kinematic bounds on where normalisation parameter should apply
205 //We set a bool to see if any bounds exist so we can short-circuit checking all of them every step
206 bool HasKinBounds = false;
207
208 if(norm.KinematicVarStr.size() > 0) HasKinBounds = true;
209
210 norm.hasKinBounds = HasKinBounds;
211 //End of kinematic bound checking
212
213 return norm;
214}
215
216// ********************************************
217// Get Base Param
218void ParameterHandlerGeneric::GetBaseParameter(const YAML::Node& param, const int Index, TypeParameterBase& Parameter) {
219// ********************************************
220 // KS: For now we don't use so avoid compilation error
221 (void) param;
222
223 Parameter.name = GetParFancyName(Index);
224
225 // Set the global parameter index of the normalisation parameter
226 Parameter.index = Index;
227}
228
229
230// ********************************************
231// Grab the global syst index for the relevant SampleName
232// i.e. get a vector of size nSplines where each entry is filled with the global syst number
233const std::vector<int> ParameterHandlerGeneric::GetGlobalSystIndexFromSampleName(const std::string& SampleName, const SystType Type) {
234// ********************************************
235 std::vector<int> returnVec;
236 for (auto &pair : _fSystToGlobalSystIndexMap[Type]) {
237 auto &SystIndex = pair.second;
238 if (AppliesToSample(SystIndex, SampleName)) { //If parameter applies to required SampleID
239 returnVec.push_back(SystIndex);
240 }
241 }
242 return returnVec;
243}
244
245// ********************************************
246// Grab the global syst index for the relevant SampleName
247// i.e. get a vector of size nSplines where each entry is filled with the global syst number
248const std::vector<int> ParameterHandlerGeneric::GetSystIndexFromSampleName(const std::string& SampleName, const SystType Type) const {
249// ********************************************
250 std::vector<int> returnVec;
251 for (auto &pair : _fSystToGlobalSystIndexMap[Type]) {
252 auto &SplineIndex = pair.first;
253 auto &systIndex = pair.second;
254 if (AppliesToSample(systIndex, SampleName)) { //If parameter applies to required SampleID
255 returnVec.push_back(SplineIndex);
256 }
257 }
258 return returnVec;
259}
260
261// ********************************************
262// Get Norm params
263SplineParameter ParameterHandlerGeneric::GetSplineParameter(const YAML::Node& param, const int Index) {
264// ********************************************
265 SplineParameter Spline;
266
267 GetBaseParameter(param, Index, Spline);
268 //Now get the Spline interpolation type
269 if (param["SplineInformation"]["InterpolationType"]){
270 for(int InterpType = 0; InterpType < kSplineInterpolations ; ++InterpType){
271 if(param["SplineInformation"]["InterpolationType"].as<std::string>() == SplineInterpolation_ToString(SplineInterpolation(InterpType)))
273 }
274 } else { //KS: By default use TSpline3
276 }
277 Spline._SplineKnotUpBound = GetFromManager<double>(param["SplineInformation"]["SplineKnotUpBound"], M3::DefSplineKnotUpBound, __FILE__ , __LINE__);
278 Spline._SplineKnotLowBound = GetFromManager<double>(param["SplineInformation"]["SplineKnotLowBound"], M3::DefSplineKnotLowBound, __FILE__ , __LINE__);
279
281 MACH3LOG_WARN("Spline knot capping enabled with bounds [{}, {}]. For reliable fits, consider modifying the input generation instead.",
283 }
284 //If there is no mode information given then this will be an empty vector
285 Spline._fSplineModes = GetFromManager(param["SplineInformation"]["Mode"], std::vector<int>(), __FILE__ , __LINE__);
286
287 return Spline;
288}
289
290// ********************************************
291// Get Func params
293// ********************************************
295 GetBaseParameter(param, Index, func);
296
297 func.pdgs = GetFromManager<std::vector<int>>(param["NeutrinoFlavour"], std::vector<int>(), __FILE__ , __LINE__);
298 func.targets = GetFromManager<std::vector<int>>(param["TargetNuclei"], std::vector<int>(), __FILE__ , __LINE__);
299 func.modes = GetFromManager<std::vector<int>>(param["Mode"], std::vector<int>(), __FILE__ , __LINE__);
300 func.preoscpdgs = GetFromManager<std::vector<int>>(param["NeutrinoFlavourUnosc"], std::vector<int>(), __FILE__ , __LINE__);
301
302 // HH - Copied from GetXsecNorm
303 int NumKinematicCuts = 0;
304 if(param["KinematicCuts"]){
305
306 NumKinematicCuts = int(param["KinematicCuts"].size());
307
308 std::vector<std::string> TempKinematicStrings;
309 std::vector<std::vector<std::vector<double>>> TempKinematicBounds;
310 //First element of TempKinematicBounds is always -999, and size is then 3
311 for(int KinVar_i = 0 ; KinVar_i < NumKinematicCuts ; ++KinVar_i){
312 //ETA: This is a bit messy, Kinematic cuts is a list of maps
313 for (YAML::const_iterator it = param["KinematicCuts"][KinVar_i].begin();it!=param["KinematicCuts"][KinVar_i].end();++it) {
314 TempKinematicStrings.push_back(it->first.as<std::string>());
315 TempKinematicBounds.push_back(Get2DBounds(it->second));
316 }
317 if(TempKinematicStrings.size() == 0) {
318 MACH3LOG_ERROR("Received a KinematicCuts node but couldn't read the contents (it's a list of single-element dictionaries (python) = map of pairs (C++))");
319 MACH3LOG_ERROR("For Param {}", func.name);
320 throw MaCh3Exception(__FILE__, __LINE__);
321 }
322 }//KinVar_i
323 func.KinematicVarStr = TempKinematicStrings;
324 func.Selection = TempKinematicBounds;
325 }
326 func.valuePtr = RetPointer(Index);
327 return func;
328}
329
330// ********************************************
331// Get Osc params
333// ********************************************
334 OscillationParameter OscParamInfo;
335 GetBaseParameter(param, Index, OscParamInfo);
336
337 return OscParamInfo;
338}
339
340// ********************************************
341// HH: Grab the Functional parameters for the relevant SampleName
342const std::vector<FunctionalParameter> ParameterHandlerGeneric::GetFunctionalParametersFromSampleName(const std::string& SampleName) const {
343// ********************************************
344 return GetTypeParamsFromSampleName(_fSystToGlobalSystIndexMap[SystType::kFunc], FuncParams, SampleName);
345}
346
347// ********************************************
348// DB Grab the Normalisation parameters for the relevant SampleName
349const std::vector<NormParameter> ParameterHandlerGeneric::GetNormParsFromSampleName(const std::string& SampleName) const {
350// ********************************************
351 return GetTypeParamsFromSampleName(_fSystToGlobalSystIndexMap[SystType::kNorm], NormParams, SampleName);
352}
353
354// ********************************************
355// KS Grab the Spline parameters for the relevant SampleName
356const std::vector<SplineParameter> ParameterHandlerGeneric::GetSplineParsFromSampleName(const std::string& SampleName) const {
357// ********************************************
358 return GetTypeParamsFromSampleName(_fSystToGlobalSystIndexMap[SystType::kSpline], SplineParams, SampleName);
359}
360
361// ********************************************
362template<typename ParamT>
363std::vector<ParamT> ParameterHandlerGeneric::GetTypeParamsFromSampleName(const std::map<int, int>& indexMap, const std::vector<ParamT>& params, const std::string& SampleName) const {
364// ********************************************
365 std::vector<ParamT> returnVec;
366 for (const auto& pair : indexMap) {
367 const auto& localIndex = pair.first;
368 const auto& globalIndex = pair.second;
369 if (AppliesToSample(globalIndex, SampleName)) {
370 returnVec.push_back(params[localIndex]);
371 }
372 }
373 return returnVec;
374}
375
376// ********************************************
377// DB Grab the number of parameters for the relevant SampleName
378int ParameterHandlerGeneric::GetNumParamsFromSampleName(const std::string& SampleName, const SystType Type) {
379// ********************************************
380 int returnVal = 0;
381 IterateOverParams(SampleName,
382 [&](int i) { return GetParamType(i) == Type; }, // Filter condition
383 [&](int) { returnVal += 1; } // Action to perform if filter passes
384 );
385 return returnVal;
386}
387
388// ********************************************
389// DB Grab the parameter names for the relevant SampleName
390const std::vector<std::string> ParameterHandlerGeneric::GetParsNamesFromSampleName(const std::string& SampleName, const SystType Type) {
391// ********************************************
392 std::vector<std::string> returnVec;
393 IterateOverParams(SampleName,
394 [&](int i) { return GetParamType(i) == Type; }, // Filter condition
395 [&](int i) { returnVec.push_back(GetParFancyName(i)); } // Action to perform if filter passes
396 );
397 return returnVec;
398}
399
400// ********************************************
401// DB DB Grab the parameter indices for the relevant SampleName
402const std::vector<int> ParameterHandlerGeneric::GetParsIndexFromSampleName(const std::string& SampleName, const SystType Type) {
403// ********************************************
404 std::vector<int> returnVec;
405 IterateOverParams(SampleName,
406 [&](int i) { return GetParamType(i) == Type; }, // Filter condition
407 [&](int i) { returnVec.push_back(i); } // Action to perform if filter passes
408 );
409 return returnVec;
410}
411
412// ********************************************
413template <typename FilterFunc, typename ActionFunc>
414void ParameterHandlerGeneric::IterateOverParams(const std::string& SampleName, FilterFunc filter, ActionFunc action) {
415// ********************************************
416 for (int i = 0; i < _fNumPar; ++i) {
417 if ((AppliesToSample(i, SampleName)) && filter(i)) { // Common filter logic
418 action(i); // Specific action for each function
419 }
420 }
421}
422
423// ********************************************
425// ********************************************
426 for (int i = 0; i < _fNumPar; ++i) {
427 //ETA - set the name to be xsec_% as this is what ProcessorMCMC expects
428 _fNames[i] = "xsec_"+std::to_string(i);
429
430 // KS: Plenty
431 if(_fParamType[i] == kOsc){
432 _fNames[i] = _fFancyNames[i];
433
434 if(_ParameterGroup[i] != "Osc"){
435 MACH3LOG_ERROR("Parameter {}, is of type Oscillation but doesn't belong to Osc group", _fFancyNames[i]);
436 MACH3LOG_ERROR("It belongs to {} group", _ParameterGroup[i]);
437 throw MaCh3Exception(__FILE__ , __LINE__ );
438 }
439 }
440 // Set ParameterHandler parameters (Curr = current, Prop = proposed, Sigma = step)
441 _fCurrVal[i] = _fPreFitValue[i];
442 _fPropVal[i] = _fCurrVal[i];
443 }
444 Randomize();
445 //KS: Transfer the starting parameters to the PCA basis, you don't want to start with zero..
446 if (pca) {
447 PCAObj->SetInitialParameters(_fIndivStepScale);
448 }
449}
450
451// ********************************************
452// Print everything we know about the inputs we're Getting
454// ********************************************
455 MACH3LOG_INFO("#################################################");
456 MACH3LOG_INFO("Printing ParameterHandlerGeneric:");
457
459
461
463
465
467
469
470 MACH3LOG_INFO("Finished");
471 MACH3LOG_INFO("#################################################");
472
474} // End
475
476// ********************************************
478// ********************************************
479 MACH3LOG_INFO("============================================================================================================================================================");
480 MACH3LOG_INFO("{:<5} {:2} {:<40} {:2} {:<10} {:2} {:<10} {:2} {:<10} {:2} {:<10} {:2} {:<10} {:2} {:<20} {:2} {:<10}", "#", "|", "Name", "|", "Prior", "|", "Error", "|", "Lower", "|", "Upper", "|", "StepScale", "|", "SampleNames", "|", "Type");
481 MACH3LOG_INFO("------------------------------------------------------------------------------------------------------------------------------------------------------------");
482 for (int i = 0; i < GetNumParams(); i++) {
483 std::string ErrString = fmt::format("{:.2f}", _fError[i]);
484 std::string SampleNameString = "";
485 for (const auto& SampleName : _fSampleNames[i]) {
486 if (!SampleNameString.empty()) {
487 SampleNameString += ", ";
488 }
489 SampleNameString += SampleName;
490 }
491 MACH3LOG_INFO("{:<5} {:2} {:<40} {:2} {:<10} {:2} {:<10} {:2} {:<10} {:2} {:<10} {:2} {:<10} {:2} {:<20} {:2} {:<10}", i, "|", GetParFancyName(i), "|", _fPreFitValue[i], "|", "+/- " + ErrString, "|", _fLowBound[i], "|", _fUpBound[i], "|", _fIndivStepScale[i], "|", SampleNameString, "|", SystType_ToString(_fParamType[i]));
492 }
493 MACH3LOG_INFO("============================================================================================================================================================");
494}
495
496// ********************************************
498// ********************************************
499 // Output the normalisation parameters as a sanity check!
500 MACH3LOG_INFO("Normalisation parameters: {}", NormParams.size());
501 if(_fSystToGlobalSystIndexMap[SystType::kNorm].size() == 0) return;
502
503 bool have_parameter_with_kin_bounds = false;
504
505 //KS: Consider making some class producing table..
506 MACH3LOG_INFO("┌────┬──────────┬────────────────────────────────────────┬────────────────────┬────────────────────┬────────────────────┐");
507 MACH3LOG_INFO("│{0:4}│{1:10}│{2:40}│{3:20}│{4:20}│{5:20}│", "#", "Global #", "Name", "Int. mode", "Target", "pdg");
508 MACH3LOG_INFO("├────┼──────────┼────────────────────────────────────────┼────────────────────┼────────────────────┼────────────────────┤");
509
510 for (unsigned int i = 0; i < NormParams.size(); ++i)
511 {
512 std::string intModeString;
513 for (unsigned int j = 0; j < NormParams[i].modes.size(); j++) {
514 intModeString += std::to_string(NormParams[i].modes[j]);
515 intModeString += " ";
516 }
517 if (NormParams[i].modes.empty()) intModeString += "all";
518
519 std::string targetString;
520 for (unsigned int j = 0; j < NormParams[i].targets.size(); j++) {
521 targetString += std::to_string(NormParams[i].targets[j]);
522 targetString += " ";
523 }
524 if (NormParams[i].targets.empty()) targetString += "all";
525
526 std::string pdgString;
527 for (unsigned int j = 0; j < NormParams[i].pdgs.size(); j++) {
528 pdgString += std::to_string(NormParams[i].pdgs[j]);
529 pdgString += " ";
530 }
531 if (NormParams[i].pdgs.empty()) pdgString += "all";
532
533 MACH3LOG_INFO("│{: <4}│{: <10}│{: <40}│{: <20}│{: <20}│{: <20}│", i, NormParams[i].index, NormParams[i].name, intModeString, targetString, pdgString);
534
535 if(NormParams[i].hasKinBounds) have_parameter_with_kin_bounds = true;
536 }
537 MACH3LOG_INFO("└────┴──────────┴────────────────────────────────────────┴────────────────────┴────────────────────┴────────────────────┘");
538
539 if(have_parameter_with_kin_bounds) {
540 MACH3LOG_INFO("Normalisation parameters KinematicCuts information");
541 MACH3LOG_INFO("┌────┬──────────┬────────────────────────────────────────┬────────────────────┬────────────────────────────────────────┐");
542 MACH3LOG_INFO("│{0:4}│{1:10}│{2:40}│{3:20}│{4:40}│", "#", "Global #", "Name", "KinematicCut", "Value");
543 MACH3LOG_INFO("├────┼──────────┼────────────────────────────────────────┼────────────────────┼────────────────────────────────────────┤");
544 for (unsigned int i = 0; i < NormParams.size(); ++i)
545 {
546 //skip parameters with no KinematicCuts
547 if(!NormParams[i].hasKinBounds) continue;
548
549 const long unsigned int ncuts = NormParams[i].KinematicVarStr.size();
550 for(long unsigned int icut = 0; icut < ncuts; icut++) {
551 std::string kinematicCutValueString;
552 for(const auto & value : NormParams[i].Selection[icut]) {
553 for (const auto& v : value) {
554 kinematicCutValueString += fmt::format("{:.2f} ", v);
555 }
556 }
557 if(icut == 0)
558 MACH3LOG_INFO("│{: <4}│{: <10}│{: <40}│{: <20}│{: <40}│", i, NormParams[i].index, NormParams[i].name, NormParams[i].KinematicVarStr[icut], kinematicCutValueString);
559 else
560 MACH3LOG_INFO("│{: <4}│{: <10}│{: <40}│{: <20}│{: <40}│", "", "", "", NormParams[i].KinematicVarStr[icut], kinematicCutValueString);
561 }//icut
562 }//i
563 MACH3LOG_INFO("└────┴──────────┴────────────────────────────────────────┴────────────────────┴────────────────────────────────────────┘");
564 }
565 else
566 MACH3LOG_INFO("No normalisation parameters have KinematicCuts defined");
567}
568
569// ********************************************
571// ********************************************
572 MACH3LOG_INFO("Spline parameters: {}", _fSystToGlobalSystIndexMap[SystType::kSpline].size());
573 if(_fSystToGlobalSystIndexMap[SystType::kSpline].size() == 0) return;
574 MACH3LOG_INFO("=====================================================================================================================================================================");
575 MACH3LOG_INFO("{:<4} {:<2} {:<40} {:<2} {:<40} {:<2} {:<20} {:<2} {:<20} {:<2} {:<20} {:<2}", "#", "|", "Name", "|", "Spline Name", "|", "Spline Interpolation", "|", "Low Knot Bound", "|", "Up Knot Bound", "|");
576 MACH3LOG_INFO("---------------------------------------------------------------------------------------------------------------------------------------------------------------------");
577 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kSpline]) {
578 auto &SplineIndex = pair.first;
579 auto &GlobalIndex = pair.second;
580
581 MACH3LOG_INFO("{:<4} {:<2} {:<40} {:<2} {:<40} {:<2} {:<20} {:<2} {:<20} {:<2} {:<20} {:<2}",
582 SplineIndex, "|", GetParFancyName(GlobalIndex), "|",
583 _fSplineNames[SplineIndex], "|",
585 GetParSplineKnotLowerBound(SplineIndex), "|",
586 GetParSplineKnotUpperBound(SplineIndex), "|");
587 }
588 MACH3LOG_INFO("=====================================================================================================================================================================");
589}
590
591// ********************************************
593// ********************************************
594 MACH3LOG_INFO("Functional parameters: {}", _fSystToGlobalSystIndexMap[SystType::kFunc].size());
595 if(_fSystToGlobalSystIndexMap[SystType::kFunc].size() == 0) return;
596 MACH3LOG_INFO("┌────┬──────────┬────────────────────────────────────────┐");
597 MACH3LOG_INFO("│{0:4}│{1:10}│{2:40}│", "#", "Global #", "Name");
598 MACH3LOG_INFO("├────┼──────────┼────────────────────────────────────────┤");
599 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kFunc]) {
600 auto &FuncIndex = pair.first;
601 auto &GlobalIndex = pair.second;
602 MACH3LOG_INFO("│{0:4}│{1:<10}│{2:40}│", std::to_string(FuncIndex), GlobalIndex, GetParFancyName(GlobalIndex));
603 }
604 MACH3LOG_INFO("└────┴──────────┴────────────────────────────────────────┘");
605}
606
607// ********************************************
609// ********************************************
610 MACH3LOG_INFO("Oscillation parameters: {}", _fSystToGlobalSystIndexMap[SystType::kOsc].size());
611 if(_fSystToGlobalSystIndexMap[SystType::kOsc].size() == 0) return;
612 MACH3LOG_INFO("┌────┬──────────┬────────────────────────────────────────┐");
613 MACH3LOG_INFO("│{0:4}│{1:10}│{2:40}│", "#", "Global #", "Name");
614 MACH3LOG_INFO("├────┼──────────┼────────────────────────────────────────┤");
615 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kOsc]) {
616 auto &OscIndex = pair.first;
617 auto &GlobalIndex = pair.second;
618 MACH3LOG_INFO("│{0:4}│{1:<10}│{2:40}│", std::to_string(OscIndex), GlobalIndex, GetParFancyName(GlobalIndex));
619 }
620 MACH3LOG_INFO("└────┴──────────┴────────────────────────────────────────┘");
621}
622
623// ********************************************
625// ********************************************
626 // KS: Create a map to store the counts of unique strings, in principle this could be in header file
627 std::unordered_map<std::string, int> paramCounts;
628
629 std::for_each(_ParameterGroup.begin(), _ParameterGroup.end(),
630 [&paramCounts](const std::string& param) {
631 paramCounts[param]++;
632 });
633
634 MACH3LOG_INFO("Printing parameter groups");
635 // Output the counts
636 for (const auto& pair : paramCounts) {
637 MACH3LOG_INFO("Found {}: {} params", pair.second, pair.first);
638 }
639}
640
641// ********************************************
642// KS: Check if matrix is correctly initialised
644// ********************************************
645 // KS: Lambda Function which simply checks if there are no duplicates in std::vector
646 auto CheckForDuplicates = [](const std::vector<std::string>& names, const std::string& nameType) {
647 std::unordered_map<std::string, size_t> seenStrings;
648 for (size_t i = 0; i < names.size(); ++i) {
649 const auto& name = names[i];
650 if (seenStrings.find(name) != seenStrings.end()) {
651 size_t firstIndex = seenStrings[name];
652 MACH3LOG_CRITICAL("There are two systematics with the same {} '{}', first at index {}, and again at index {}", nameType, name, firstIndex, i);
653 throw MaCh3Exception(__FILE__, __LINE__);
654 }
655 seenStrings[name] = i;
656 }
657 };
658
659 // KS: Checks if there are no duplicates in fancy names etc, this can happen if we merge configs etc
660 CheckForDuplicates(_fFancyNames, "_fFancyNames");
661 CheckForDuplicates(_fSplineNames, "_fSplineNames");
662}
663
664// ********************************************
665// Function to set to prior parameters of a given group
666void ParameterHandlerGeneric::SetGroupOnlyParameters(const std::vector< std::string>& Groups) {
667// ********************************************
668 for(size_t i = 0; i < Groups.size(); i++){
669 SetGroupOnlyParameters(Groups[i]);
670 }
671}
672
673// ********************************************
674// Function to set to prior parameters of a given group
675void ParameterHandlerGeneric::SetGroupOnlyParameters(const std::string& Group, const std::vector<double>& Pars) {
676// ********************************************
677 // If empty, set the proposed to prior
678 if (Pars.empty()) {
679 for (int i = 0; i < _fNumPar; i++) {
680 if(IsParFromGroup(i, Group)) _fPropVal[i] = _fPreFitValue[i];
681 }
682 } else{
683 const size_t ExpectedSize = static_cast<size_t>(GetNumParFromGroup(Group));
684 if (Pars.size() != ExpectedSize) {
685 MACH3LOG_ERROR("Number of param in group {} is {}, while you passed {}", Group, ExpectedSize, Pars.size());
686 throw MaCh3Exception(__FILE__, __LINE__);
687 }
688 int Counter = 0;
689 for (int i = 0; i < _fNumPar; i++) {
690 // If belongs to group set value from parsed vector, otherwise use propose value
691 if(IsParFromGroup(i, Group)){
692 _fPropVal[i] = Pars[Counter];
693 Counter++;
694 }
695 }
696 }
697 // And if pca make the transfer
698 if (pca) {
699 PCAObj->TransferToPCA();
700 PCAObj->TransferToParam();
701 }
702}
703
704// ********************************************
705// Checks if parameter belongs to a given group
706bool ParameterHandlerGeneric::IsParFromGroup(const int i, const std::string& Group) const {
707// ********************************************
708 std::string groupLower = Group;
709 std::string paramGroupLower = _ParameterGroup[i];
710
711 // KS: Convert both strings to lowercase, this way comparison will be case insensitive
712 std::transform(groupLower.begin(), groupLower.end(), groupLower.begin(), ::tolower);
713 std::transform(paramGroupLower.begin(), paramGroupLower.end(), paramGroupLower.begin(), ::tolower);
714
715 return groupLower == paramGroupLower;
716}
717
718// ********************************************
719int ParameterHandlerGeneric::GetNumParFromGroup(const std::string& Group) const {
720// ********************************************
721 int Counter = 0;
722 for (int i = 0; i < _fNumPar; i++) {
723 if(IsParFromGroup(i, Group)) Counter++;
724 }
725 return Counter;
726}
727
728// ********************************************
729// DB Grab the Normalisation parameters for the relevant sample name
730std::vector<const double*> ParameterHandlerGeneric::GetOscParsFromSampleName(const std::string& SampleName) {
731// ********************************************
732 std::vector<const double*> returnVec;
733 for (const auto& pair : _fSystToGlobalSystIndexMap[SystType::kOsc]) {
734 const auto& globalIndex = pair.second;
735 if (AppliesToSample(globalIndex, SampleName)) {
736 returnVec.push_back(RetPointer(globalIndex));
737 }
738 }
739 return returnVec;
740}
741
742// ********************************************
743// Dump Matrix to ROOT file, useful when we need to pass matrix info to another fitting group
744void ParameterHandlerGeneric::DumpMatrixToFile(const std::string& Name) {
745// ********************************************
746 TFile* outputFile = new TFile(Name.c_str(), "RECREATE");
747
748 TObjArray* xsec_param_names = new TObjArray();
749 TObjArray* xsec_spline_interpolation = new TObjArray();
750 TObjArray* xsec_spline_names = new TObjArray();
751
752 TVectorD* xsec_param_prior = new TVectorD(_fNumPar);
753 TVectorD* xsec_flat_prior = new TVectorD(_fNumPar);
754 TVectorD* xsec_stepscale = new TVectorD(_fNumPar);
755 TVectorD* xsec_param_lb = new TVectorD(_fNumPar);
756 TVectorD* xsec_param_ub = new TVectorD(_fNumPar);
757
758 TVectorD* xsec_param_knot_weight_lb = new TVectorD(_fNumPar);
759 TVectorD* xsec_param_knot_weight_ub = new TVectorD(_fNumPar);
760 TVectorD* xsec_error = new TVectorD(_fNumPar);
761
762 for(int i = 0; i < _fNumPar; ++i)
763 {
764 TObjString* nameObj = new TObjString(_fFancyNames[i].c_str());
765 xsec_param_names->AddLast(nameObj);
766
767 TObjString* splineType = new TObjString("TSpline3");
768 xsec_spline_interpolation->AddLast(splineType);
769
770 TObjString* splineName = new TObjString("");
771 xsec_spline_names->AddLast(splineName);
772
773 (*xsec_param_prior)[i] = _fPreFitValue[i];
774 (*xsec_flat_prior)[i] = _fFlatPrior[i];
775 (*xsec_stepscale)[i] = _fIndivStepScale[i];
776 (*xsec_error)[i] = _fError[i];
777
778 (*xsec_param_lb)[i] = _fLowBound[i];
779 (*xsec_param_ub)[i] = _fUpBound[i];
780
781 //Default values
782 (*xsec_param_knot_weight_lb)[i] = -9999;
783 (*xsec_param_knot_weight_ub)[i] = +9999;
784 }
785
786 for (auto &pair : _fSystToGlobalSystIndexMap[SystType::kSpline]) {
787 auto &SplineIndex = pair.first;
788 auto &SystIndex = pair.second;
789
790 (*xsec_param_knot_weight_lb)[SystIndex] = SplineParams.at(SplineIndex)._SplineKnotLowBound;
791 (*xsec_param_knot_weight_ub)[SystIndex] = SplineParams.at(SplineIndex)._SplineKnotUpBound;
792
793 TObjString* splineType = new TObjString(SplineInterpolation_ToString(SplineParams.at(SplineIndex)._SplineInterpolationType).c_str());
794 xsec_spline_interpolation->AddAt(splineType, SystIndex);
795
796 TObjString* splineName = new TObjString(_fSplineNames[SplineIndex].c_str());
797 xsec_spline_names->AddAt(splineName, SystIndex);
798 }
799 xsec_param_names->Write("xsec_param_names", TObject::kSingleKey);
800 delete xsec_param_names;
801 xsec_spline_interpolation->Write("xsec_spline_interpolation", TObject::kSingleKey);
802 delete xsec_spline_interpolation;
803 xsec_spline_names->Write("xsec_spline_names", TObject::kSingleKey);
804 delete xsec_spline_names;
805
806 xsec_param_prior->Write("xsec_param_prior");
807 delete xsec_param_prior;
808 xsec_flat_prior->Write("xsec_flat_prior");
809 delete xsec_flat_prior;
810 xsec_stepscale->Write("xsec_stepscale");
811 delete xsec_stepscale;
812 xsec_param_lb->Write("xsec_param_lb");
813 delete xsec_param_lb;
814 xsec_param_ub->Write("xsec_param_ub");
815 delete xsec_param_ub;
816
817 xsec_param_knot_weight_lb->Write("xsec_param_knot_weight_lb");
818 delete xsec_param_knot_weight_lb;
819 xsec_param_knot_weight_ub->Write("xsec_param_knot_weight_ub");
820 delete xsec_param_knot_weight_ub;
821 xsec_error->Write("xsec_error");
822 delete xsec_error;
823
824 covMatrix->Write("xsec_cov");
825 TH2D* CorrMatrix = GetCorrelationMatrix();
826 CorrMatrix->Write("hcov");
827 delete CorrMatrix;
828
829 outputFile->Close();
830 delete outputFile;
831
832 MACH3LOG_INFO("Finished dumping ParameterHandler object");
833}
int size
#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
#define MACH3LOG_WARN
Definition: MaCh3Logger.h:24
std::string SplineInterpolation_ToString(const SplineInterpolation i)
Convert a LLH type to a string.
SplineInterpolation
Make an enum of the spline interpolation type.
@ kTSpline3
Default TSpline3 interpolation.
@ kSplineInterpolations
This only enumerates.
std::string SystType_ToString(const SystType i)
Convert a Syst type type to a string.
SystType
@ kOsc
For oscillation parameters.
Type GetFromManager(const YAML::Node &node, Type defval, const std::string File="", const int Line=1)
Get content of config file if node is not found take default value specified.
Definition: YamlHelper.h:298
#define Get2DBounds(filename)
Definition: YamlHelper.h:562
Custom exception class for MaCh3 errors.
Base class responsible for handling of systematic error parameters. Capable of using PCA or using ada...
const double * RetPointer(const int iParam)
DB Pointer return to param position.
bool AppliesToSample(const int SystIndex, const std::string &SampleName) const
Check if parameter is affecting given sample name.
std::vector< std::string > _fFancyNames
Fancy name for example rather than xsec_0 it is MAQE, useful for human reading.
std::vector< bool > _fFlatPrior
Whether to apply flat prior or not.
std::vector< double > _fError
Prior error on the parameter.
YAML::Node _fYAMLDoc
Stores config describing systematics.
std::unique_ptr< PCAHandler > PCAObj
Struct containing information about PCA.
int _fNumPar
Number of systematic parameters.
std::vector< double > _fLowBound
Lowest physical bound, parameter will not be able to go beyond it.
std::vector< double > _fCurrVal
Current value of the parameter.
TMatrixDSym * covMatrix
The covariance matrix.
std::vector< double > _fPreFitValue
Parameter value dictated by the prior model. Based on it penalty term is calculated.
void Randomize() _noexcept_
"Randomize" the parameters in the covariance class for the proposed step. Used the proposal kernel an...
std::vector< std::string > _fNames
ETA _fNames is set automatically in the covariance class to be something like xsec_i,...
std::vector< double > _fUpBound
Upper physical bound, parameter will not be able to go beyond it.
bool pca
perform PCA or not
std::vector< double > _fIndivStepScale
Individual step scale used by MCMC algorithm.
int PrintLength
KS: This is used when printing parameters, sometimes we have super long parameters name,...
std::vector< double > _fPropVal
Proposed value of the parameter.
std::vector< std::vector< std::string > > _fSampleNames
Tells to which samples object param should be applied.
void PrintNormParams()
Prints normalization parameters.
FunctionalParameter GetFunctionalParameters(const YAML::Node &param, const int Index)
Get Func params.
std::vector< SplineParameter > SplineParams
Vector containing info for normalisation systematics.
void PrintParameterGroups()
Prints groups of parameters.
SplineParameter GetSplineParameter(const YAML::Node &param, const int Index)
Get Spline params.
void Print()
Print information about the whole object once it is set.
ParameterHandlerGeneric(const std::vector< std::string > &FileNames, std::string name="xsec_cov", double threshold=-1, int FirstPCAdpar=-999, int LastPCAdpar=-999)
Constructor.
std::vector< std::string > _fSplineNames
Name of spline in TTree (TBranch),.
void PrintOscillationParams()
Prints oscillation parameters.
void PrintSplineParams()
Prints spline parameters.
void InitParams()
Initializes the systematic parameters from the configuration file. This function loads parameters lik...
std::vector< NormParameter > NormParams
Vector containing info for normalisation systematics.
std::vector< std::map< int, int > > _fSystToGlobalSystIndexMap
Map between number of given parameter type with global parameter numbering. For example 2nd norm para...
std::vector< SystType > _fParamType
Type of parameter like norm, spline etc.
void PrintFunctionalParams()
Prints functional parameters.
void CheckCorrectInitialisation()
KS: Check if matrix is correctly initialised.
void PrintGlobablInfo()
Prints general information about the ParameterHandler object.
NormParameter GetNormParameter(const YAML::Node &param, const int Index)
Get Norm params.
void InitParametersTypeFromConfig()
Parses the YAML configuration to set up cross-section parameters. The YAML file defines the types of ...
std::vector< ParamT > GetTypeParamsFromSampleName(const std::map< int, int > &indexMap, const std::vector< ParamT > &params, const std::string &SampleName) const
Retrieve parameters that apply to a given sample name.
std::vector< FunctionalParameter > FuncParams
Vector containing info for functional systematics.
std::vector< OscillationParameter > OscParams
Vector containing info for functional systematics.
void DumpMatrixToFile(const std::string &Name)
Dump Matrix to ROOT file, useful when we need to pass matrix info to another fitting group.
SplineInterpolation GetParSplineInterpolation(const int i) const
Get interpolation type for a given parameter.
OscillationParameter GetOscillationParameters(const YAML::Node &param, const int Index)
Get Osc params.
void IterateOverParams(const std::string &SampleName, FilterFunc filter, ActionFunc action)
Iterates over parameters and applies a filter and action function.
void GetBaseParameter(const YAML::Node &param, const int Index, TypeParameterBase &Parameter)
Fill base parameters.
std::vector< std::string > _ParameterGroup
KS: Allow to group parameters for example to affect only cross-section or only flux etc.
int GetNumParams() const
Get total number of parameters.
const std::vector< NormParameter > GetNormParsFromSampleName(const std::string &SampleName) const
DB Get norm/func parameters depending on given SampleName.
SystType GetParamType(const int i) const
Returns enum describing our param type.
const std::vector< int > GetSystIndexFromSampleName(const std::string &SampleName, const SystType Type) const
Grab the index of the syst relative to global numbering.
const std::vector< int > GetGlobalSystIndexFromSampleName(const std::string &SampleName, const SystType Type)
DB Get spline parameters depending on given SampleName.
double GetParSplineKnotUpperBound(const int i) const
EM: value at which we cap spline knot weight.
std::string GetParFancyName(const int i) const
Get fancy name of the Parameter.
TH2D * GetCorrelationMatrix()
KS: Convert covariance matrix to correlation matrix and return TH2D which can be used for fancy plott...
const std::vector< int > GetParsIndexFromSampleName(const std::string &SampleName, const SystType Type)
DB Grab the parameter indices for the relevant SampleName.
const std::vector< SplineParameter > GetSplineParsFromSampleName(const std::string &SampleName) const
KS: Grab the Spline parameters for the relevant SampleName.
const std::vector< FunctionalParameter > GetFunctionalParametersFromSampleName(const std::string &SampleName) const
HH Get functional parameters for the relevant 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 > GetParsNamesFromSampleName(const std::string &SampleName, const SystType Type)
DB Grab the parameter names for the relevant SampleName.
const std::vector< std::string > GetSplineParsNamesFromSampleName(const std::string &SampleName)
DB Get spline parameters depending on given SampleName.
bool IsParFromGroup(const int i, const std::string &Group) const
Checks if parameter belongs to a given group.
const std::vector< SplineInterpolation > GetSplineInterpolationFromSampleName(const std::string &SampleName)
Get the interpolation types for splines affecting a particular SampleName.
double GetParSplineKnotLowerBound(const int i) const
EM: value at which we cap spline knot weight.
std::vector< const double * > GetOscParsFromSampleName(const std::string &SampleName)
Get pointers to Osc params from Sample name.
const std::vector< std::vector< int > > GetSplineModeVecFromSampleName(const std::string &SampleName)
DB Grab the Spline Modes for the relevant SampleName.
int GetNumParFromGroup(const std::string &Group) const
KS: Check how many parameters are associated with given group.
void SetGroupOnlyParameters(const std::string &Group, const std::vector< double > &Pars={})
KS Function to set to prior parameters of a given group or values from vector.
static constexpr const double DefSplineKnotUpBound
Default value for spline knot capping, default mean not capping is being applied.
Definition: Core.h:76
static constexpr const double DefSplineKnotLowBound
Default value for spline knot capping, default mean not capping is being applied.
Definition: Core.h:78
HH - Functional parameters Carrier for whether you want to apply a systematic to an event or not.
std::vector< int > modes
Mode which parameter applies to.
std::vector< std::string > KinematicVarStr
const double * valuePtr
Parameter value pointer.
std::vector< int > targets
Targets which parameter applies to.
std::vector< int > pdgs
PDG which parameter applies to.
std::vector< std::vector< std::vector< double > > > Selection
std::vector< int > preoscpdgs
Preosc PDG which parameter applies to.
ETA - Normalisations for cross-section parameters Carrier for whether you want to apply a systematic ...
std::vector< int > preoscpdgs
Preosc PDG which parameter applies to.
bool hasKinBounds
Does this parameter have kinematic bounds.
std::vector< std::vector< std::vector< double > > > Selection
std::vector< int > modes
Mode which parameter applies to.
std::vector< int > pdgs
PDG which parameter applies to.
std::vector< std::string > KinematicVarStr
std::vector< int > targets
Targets which parameter applies to.
KS: Struct holding info about oscillation Systematics.
KS: Struct holding info about Spline Systematics.
SplineInterpolation _SplineInterpolationType
Spline interpolation vector.
double _SplineKnotUpBound
EM: Cap spline knot higher value.
double _SplineKnotLowBound
EM: Cap spline knot lower value.
std::vector< int > _fSplineModes
Modes to which spline applies (valid only for binned splines)
Base class storing info for parameters types, helping unify codebase.
int index
Parameter number of this normalisation in current systematic model.
std::string name
Name of parameters.