MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
ParameterStructs.h
Go to the documentation of this file.
1#pragma once
2
3// C++ includes
4#include <set>
5#include <list>
6#include <unordered_map>
7
8// MaCh3 includes
10#include "Manager/MaCh3Logger.h"
11#include "Manager/Core.h"
12
14// ROOT include
15#include "TSpline.h"
16#include "TObjString.h"
17#include "TFile.h"
18#include "TF1.h"
19#include "TH2Poly.h"
20#include "TH1.h"
22
30
31// *******************
33template< typename T, size_t N >
34std::vector<T> MakeVector( const T (&data)[N] ) {
35// *******************
36 return std::vector<T>(data, data+N);
37}
38
39// *******************
41template <typename T>
42void CleanVector(std::vector<T>& vec) {
43// *******************
44 vec.clear();
45 vec.shrink_to_fit();
46}
47
48// *******************
50template <typename T>
51void CleanVector(std::vector<std::vector<T>>& vec) {
52// *******************
53 for (auto& innerVec : vec) {
54 innerVec.clear();
55 innerVec.shrink_to_fit();
56 }
57 vec.clear();
58 vec.shrink_to_fit();
59}
60
61// *******************
63template <typename T>
64void CleanVector(std::vector<std::vector<std::vector<T>>>& vec) {
65// *******************
66 for (auto& inner2DVec : vec) {
67 for (auto& innerVec : inner2DVec) {
68 innerVec.clear();
69 innerVec.shrink_to_fit();
70 }
71 inner2DVec.clear();
72 inner2DVec.shrink_to_fit();
73 }
74 vec.clear();
75 vec.shrink_to_fit();
76}
77
78// *******************
81template <typename T>
82void CleanVector(std::vector<std::vector<std::vector<std::vector<std::vector<std::vector<std::vector<T>>>>>>> &vec) {
83// *******************
84 for (auto& v6 : vec) {
85 for (auto& v5 : v6) {
86 for (auto& v4 : v5) {
87 for (auto& v3 : v4) {
88 for (auto& v2 : v3) {
89 for (auto& v1 : v2) {
90 v1.clear();
91 v1.shrink_to_fit();
92 }
93 v2.clear();
94 v2.shrink_to_fit();
95 }
96 v3.clear();
97 v3.shrink_to_fit();
98 }
99 v4.clear();
100 v4.shrink_to_fit();
101 }
102 v5.clear();
103 v5.shrink_to_fit();
104 }
105 v6.clear();
106 v6.shrink_to_fit();
107 }
108 vec.clear();
109 vec.shrink_to_fit();
110}
111
112// *******************
114template <typename T>
115void CleanContainer(std::vector<T*>& container) {
116// *******************
117 for (T* ptr : container) {
118 if(ptr) delete ptr;
119 ptr = nullptr;
120 }
121 container.clear();
122 container.shrink_to_fit();
123}
124
125// *******************
127template <typename T>
128void CleanContainer(std::vector< std::vector< std::vector<T*> > >& container) {
129// *******************
130 for (auto& container2D : container) {
131 for (auto& container1D : container2D) {
132 for (T* ptr : container1D) {
133 if(ptr) delete ptr;
134 }
135 container1D.clear();
136 container1D.shrink_to_fit();
137 }
138 container2D.clear();
139 container2D.shrink_to_fit();
140 }
141 container.clear();
142 container.shrink_to_fit();
143}
144
145// *******************
148// *******************
150 std::string name;
151
154};
155
156// *******************
161 std::vector<int> modes;
163 std::vector<int> horncurrents;
165 std::vector<int> pdgs;
167 std::vector<int> preoscpdgs;
169 std::vector<int> targets;
174 std::vector< std::vector< std::vector<double> > > Selection;
175
180 std::vector< std::string > KinematicVarStr;
181};
182
183// HH - a shorthand type for funcpar functions
184using FuncParFuncType = std::function<void (const double*, std::size_t)>;
185// *******************
189// *******************
191 std::vector<int> modes;
193 std::vector<int> horncurrents;
195 std::vector<int> pdgs;
197 std::vector<int> preoscpdgs;
199 std::vector<int> targets;
204 std::vector< std::vector< std::vector<double> > > Selection;
205
210 std::vector< std::string > KinematicVarStr;
211
213 const double* valuePtr;
214
217};
218
225
235
236// **************************************************
239inline std::string GetTF1(const SplineInterpolation i) {
240 // **************************************************
241 std::string Func = "";
242 switch(i) {
244 Func = "([1]+[0]*x)";
245 break;
251 MACH3LOG_ERROR("Interpolation type {} not supported for TF1!", static_cast<int>(i));
252 throw MaCh3Exception(__FILE__, __LINE__);
253 default:
254 MACH3LOG_ERROR("UNKNOWN SPLINE INTERPOLATION SPECIFIED!");
255 MACH3LOG_ERROR("You gave {}", static_cast<int>(i));
256 throw MaCh3Exception(__FILE__ , __LINE__ );
257 }
258 return Func;
259}
260
261// **************************************************
265// **************************************************
267 switch(i) {
273 break;
276 break;
278 MACH3LOG_ERROR("kSplineInterpolations is not a valid interpolation type!");
279 throw MaCh3Exception(__FILE__, __LINE__);
280 default:
281 MACH3LOG_ERROR("UNKNOWN SPLINE INTERPOLATION SPECIFIED!");
282 MACH3LOG_ERROR("You gave {}", static_cast<int>(i));
283 throw MaCh3Exception(__FILE__, __LINE__);
284 }
285 return Type;
286}
287
288// **************************************************
291// **************************************************
292 std::string name = "";
293 switch(i) {
294 // TSpline3 (third order spline in ROOT)
296 name = "TSpline3";
297 break;
299 name = "Linear";
300 break;
302 name = "Monotonic";
303 break;
304 // (Experimental) Akima_Spline (crd order spline which is allowed to be discontinuous in 2nd deriv)
306 name = "Akima";
307 break;
309 name = "LinearFunc";
310 break;
312 MACH3LOG_ERROR("kSplineInterpolations is not a valid interpolation type!");
313 throw MaCh3Exception(__FILE__, __LINE__);
314 default:
315 MACH3LOG_ERROR("UNKNOWN SPLINE INTERPOLATION SPECIFIED!");
316 MACH3LOG_ERROR("You gave {}", static_cast<int>(i));
317 throw MaCh3Exception(__FILE__ , __LINE__ );
318 }
319 return name;
320}
321
331
332// *******************
335// *******************
338
340 std::vector<int> _fSplineModes;
341
346};
347
348// **************************************************
350inline std::string SystType_ToString(const SystType i) {
351// **************************************************
352 std::string name = "";
353 switch(i) {
354 case SystType::kNorm:
355 name = "Norm";
356 break;
358 name = "Spline";
359 break;
360 case SystType::kFunc:
361 name = "Functional";
362 break;
363 case SystType::kOsc:
364 name = "Oscillation";
365 break;
367 MACH3LOG_ERROR("kSystTypes is not a valid SystType!");
368 throw MaCh3Exception(__FILE__, __LINE__);
369 default:
370 MACH3LOG_ERROR("UNKNOWN SYST TYPE SPECIFIED!");
371 MACH3LOG_ERROR("You gave {}", static_cast<int>(i));
372 throw MaCh3Exception(__FILE__ , __LINE__ );
373 }
374 return name;
375}
376
377// *******************
381// *******************
382
383};
#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
KS: Based on this https://github.com/gabime/spdlog/blob/a2b4262090fd3f005c2315dcb5be2f0f1774a005/incl...
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:25
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.
@ kMonotonic
EM: DOES NOT make the entire spline monotonic, only the segments.
@ kSplineInterpolations
This only enumerates.
@ kLinear
Linear interpolation between knots.
@ kLinearFunc
Liner interpolation using TF1 not spline.
@ kAkima
EM: Akima spline iis allowed to be discontinuous in 2nd derivative and coefficients in any segment.
RespFuncType
Make an enum of the spline interpolation type.
@ kTF1_red
Uses TF1_red for interpolation.
@ kTSpline3_red
Uses TSpline3_red for interpolation.
@ kRespFuncTypes
This only enumerates.
std::vector< T > MakeVector(const T(&data)[N])
Template to make vector out of an array of any length.
void CleanVector(std::vector< T > &vec)
Generic cleanup function.
std::string GetTF1(const SplineInterpolation i)
Get function for TF1_red.
std::function< void(const double *, std::size_t)> FuncParFuncType
std::string SystType_ToString(const SystType i)
Convert a Syst type type to a string.
void CleanContainer(std::vector< T * > &container)
Generic cleanup function.
RespFuncType SplineInterpolation_ToRespFuncType(const SplineInterpolation i)
Convert a RespFuncType type to a SplineInterpolation.
SystType
@ kNorm
For normalisation parameters.
@ kSpline
For splined parameters (1D)
@ kSystTypes
This only enumerates.
@ kOsc
For oscillation parameters.
@ kFunc
For functional parameters.
Custom exception class for MaCh3 errors.
static constexpr const int _BAD_INT_
Default value used for int initialisation.
Definition: Core.h:45
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.
bool hasKinBounds
Does this parameter have kinematic bounds.
std::vector< int > targets
Targets which parameter applies to.
std::vector< int > horncurrents
Horn currents which parameter applies to.
std::vector< int > pdgs
PDG which parameter applies to.
FuncParFuncType * funcPtr
Function pointer.
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 > horncurrents
Horn currents 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.