-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathJetCorrectorParameters.h
206 lines (175 loc) · 8.27 KB
/
JetCorrectorParameters.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
//
// Original Author: Fedor Ratnikov Nov 9, 2007
// $Id: JetCorrectorParameters.h,v 1.14 2011/01/27 12:14:13 kkousour Exp $
//
// Generic parameters for Jet corrections
//
#ifndef JetCorrectorParameters_h
#define JetCorrectorParameters_h
#include <string>
#include <vector>
#include <algorithm>
#include <iostream>
//#include "FWCore/Utilities/interface/Exception.h"
#include "Utilities.cc"
class JetCorrectorParameters
{
//---------------- JetCorrectorParameters class ----------------
//-- Encapsulates all the information of the parametrization ---
public:
//---------------- Definitions class ---------------------------
//-- Global iformation about the parametrization is kept here --
class Definitions
{
public:
//-------- Constructors --------------
Definitions() {}
Definitions(const std::vector<std::string>& fBinVar, const std::vector<std::string>& fParVar, const std::string& fFormula, bool fIsResponse);
Definitions(const std::string& fLine);
//-------- Member functions ----------
unsigned nBinVar() const {return mBinVar.size(); }
unsigned nParVar() const {return mParVar.size(); }
std::vector<std::string> parVar() const {return mParVar; }
std::vector<std::string> binVar() const {return mBinVar; }
std::string parVar(unsigned fIndex) const {return mParVar[fIndex];}
std::string binVar(unsigned fIndex) const {return mBinVar[fIndex];}
std::string formula() const {return mFormula; }
std::string level() const {return mLevel; }
bool isResponse() const {return mIsResponse; }
private:
//-------- Member variables ----------
bool mIsResponse;
std::string mLevel;
std::string mFormula;
std::vector<std::string> mParVar;
std::vector<std::string> mBinVar;
};
//---------------- Record class --------------------------------
//-- Each Record holds the properties of a bin -----------------
class Record
{
public:
//-------- Constructors --------------
Record() : mNvar(0),mMin(0),mMax(0) {}
Record(unsigned fNvar, const std::vector<float>& fXMin, const std::vector<float>& fXMax, const std::vector<float>& fParameters) : mNvar(fNvar),mMin(fXMin),mMax(fXMax),mParameters(fParameters) {}
Record(const std::string& fLine, unsigned fNvar);
//-------- Member functions ----------
float xMin(unsigned fVar) const {return mMin[fVar]; }
float xMax(unsigned fVar) const {return mMax[fVar]; }
float xMiddle(unsigned fVar) const {return 0.5*(xMin(fVar)+xMax(fVar));}
float parameter(unsigned fIndex) const {return mParameters[fIndex]; }
std::vector<float> parameters() const {return mParameters; }
unsigned nParameters() const {return mParameters.size(); }
int operator< (const Record& other) const {return xMin(0) < other.xMin(0); }
private:
//-------- Member variables ----------
unsigned mNvar;
std::vector<float> mMin;
std::vector<float> mMax;
std::vector<float> mParameters;
};
//-------- Constructors --------------
JetCorrectorParameters() { valid_ = false;}
JetCorrectorParameters(const std::string& fFile, const std::string& fSection = "");
JetCorrectorParameters(const JetCorrectorParameters::Definitions& fDefinitions,
const std::vector<JetCorrectorParameters::Record>& fRecords)
: mDefinitions(fDefinitions),mRecords(fRecords) { valid_ = true;}
//-------- Member functions ----------
const Record& record(unsigned fBin) const {return mRecords[fBin]; }
const Definitions& definitions() const {return mDefinitions; }
unsigned size() const {return mRecords.size();}
unsigned size(unsigned fVar) const;
int binIndex(const std::vector<float>& fX) const;
int neighbourBin(unsigned fIndex, unsigned fVar, bool fNext) const;
std::vector<float> binCenters(unsigned fVar) const;
void printScreen() const;
void printFile(const std::string& fFileName) const;
bool isValid() const { return valid_; }
private:
//-------- Member variables ----------
JetCorrectorParameters::Definitions mDefinitions;
std::vector<JetCorrectorParameters::Record> mRecords;
bool valid_; /// is this a valid set?
};
class JetCorrectorParametersCollection {
//---------------- JetCorrectorParametersCollection class ----------------
//-- Adds several JetCorrectorParameters together by algorithm type ---
//-- to reduce the number of payloads in the Database ---
public:
enum Level_t { L1Offset=0,
L1JPTOffset=7,
L1FastJet = 10,
L2Relative=1,
L3Absolute=2,
L2L3Residual=8,
L4EMF=3,
L5Flavor=4,
L6UE=5,
L7Parton=6,
Uncertainty=9,
N_LEVELS=11
};
enum L5_Species_t {L5_bJ=0,L5_cJ,L5_qJ,L5_gJ,L5_bT,L5_cT,L5_qT,L5_gT,N_L5_SPECIES};
enum L7_Species_t {L7_gJ=0,L7_qJ,L7_cJ,L7_bJ,L7_jJ,L7_qT,L7_cT,L7_bT,L7_jT,N_L7_SPECIES};
typedef int key_type;
typedef std::string label_type;
typedef JetCorrectorParameters value_type;
typedef std::pair<key_type,value_type> pair_type;
typedef std::vector<pair_type> collection_type;
// Constructor... initialize all three vectors to zero
JetCorrectorParametersCollection() { corrections_.clear(); correctionsL5_.clear(); correctionsL7_.clear(); }
// Add a JetCorrectorParameter object, possibly with flavor.
void push_back( key_type i, value_type const & j, label_type const & flav = "" );
// Access the JetCorrectorParameter via the key k.
// key_type is hashed to deal with the three collections
JetCorrectorParameters const & operator[]( key_type k ) const;
// Access the JetCorrectorParameter via a string.
// Will find the hashed value for the label, and call via that
// operator.
JetCorrectorParameters const & operator[]( std::string const & label ) const {
return operator[]( findKey(label) );
}
// Get a list of valid keys. These will contain hashed keys
// that are aware of all three collections.
void validKeys(std::vector<key_type> & keys ) const;
// Helper method to find all of the sections in a given
// parameters file
static void getSections( std::string inputFile,
std::vector<std::string> & outputs );
// Find the L5 bin for hashing
static key_type getL5Bin( std::string const & flav );
// Find the L7 bin for hashing
static key_type getL7Bin( std::string const & flav );
// Check if this is an L5 hashed value
static bool isL5( key_type k );
// Check if this is an L7 hashed value
static bool isL7( key_type k );
static std::string findLabel( key_type k ){
if ( isL5(k) ) return findL5Flavor(k);
else if ( isL7(k) ) return findL7Parton(k);
else return labels_[k];
}
static std::string findL5Flavor( key_type k ){
if ( k == L5Flavor ) return labels_[L5Flavor];
else
return l5Flavors_[k / 100 - 1];
}
static std::string findL7Parton( key_type k ){
if ( k == L7Parton ) return labels_[L7Parton];
else
return l7Partons_[k / 1000 - 1];
}
protected:
// Find the key corresponding to each label
key_type findKey( std::string const & label ) const;
collection_type corrections_;
collection_type correctionsL5_;
collection_type correctionsL7_;
static const char * labelsArray_[N_LEVELS];
static std::vector<std::string> labels_;
static const char * l5FlavorArray_[N_L5_SPECIES];
static std::vector<std::string> l5Flavors_;
static const char * l7PartonArray_[N_L7_SPECIES];
static std::vector<std::string> l7Partons_;
};
#endif