CMS 3D CMS Logo
Functions.h
Go to the documentation of this file.
1 #ifndef FUNCTIONS_H
2 #define FUNCTIONS_H
3 
4 #include <iostream>
5 #include <vector>
6 #include <cmath>
7 #include "TMath.h"
8 #include "TString.h"
9 #include "TF1.h"
10 #include "TRandom.h"
11 #include "MuonAnalysis/MomentumScaleCalibration/interface/SigmaPtDiff.h"
12 
16 struct ParameterSet
17 {
18  ParameterSet() {}
19  ParameterSet(const TString & inputName, const double & inputStep, const double & inputMini, const double & inputMaxi) :
20  step(inputStep),
21  mini(inputMini),
22  maxi(inputMaxi)
23  {
24  std::cout << "setting name = " << inputName << std::endl;
25  name = inputName;
26  }
27  TString name;
28  double step, mini, maxi;
29 };
30 
31 // ----------------------- //
32 // Bias and scale functors //
33 // ----------------------- //
43 template <class T>
44 class scaleFunctionBase {
45  public:
46  virtual double scale(const double & pt, const double & eta, const double & phi, const int chg, const T & parScale) const = 0;
47  virtual ~scaleFunctionBase() = 0;
49  virtual void resetParameters(std::vector<double> * scaleVec) const {
50  std::cout << "ERROR: the resetParameters method must be defined in each scale function" << std::endl;
51  std::cout << "Please add it to the scaleFunction you are using" << std::endl;
52  exit(1);
53  }
55  virtual void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const T & parScale, const std::vector<int> & parScaleOrder, const int muonType) = 0;
56  virtual void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
57  const T & parResol, const std::vector<int> & parResolOrder,
58  const std::vector<double> & parStep,
59  const std::vector<double> & parMin,
60  const std::vector<double> & parMax,
61  const int muonType)
62  {
63  std::cout << "The method setParameters must be implemented for this scale function" << std::endl;
64  exit(1);
65  }
66  virtual int parNum() const { return parNum_; }
67  protected:
68  int parNum_;
70  virtual void setPar(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
71  TString* parname, const T & parScale, const std::vector<int> & parScaleOrder,
72  double* thisStep, double* thisMini, double* thisMaxi, TString* thisParName ) {
73  for( int iPar=0; iPar<this->parNum_; ++iPar ) {
74  Start[iPar] = parScale[iPar];
75  Step[iPar] = thisStep[iPar];
76  Mini[iPar] = thisMini[iPar];
77  Maxi[iPar] = thisMaxi[iPar];
78  ind[iPar] = parScaleOrder[iPar];
79  parname[iPar] = thisParName[iPar];
80  }
81  }
82  virtual void setPar(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
83  TString* parname, const T & parResol, const std::vector<int> & parResolOrder, const std::vector<ParameterSet> & parSet ) {
84  if( int(parSet.size()) != this->parNum_ ) {
85  std::cout << "Error: wrong number of parameter initializations = " << parSet.size() << ". Number of parameters is " << this->parNum_ << std::endl;
86  exit(1);
87  }
88  for( int iPar=0; iPar<this->parNum_; ++iPar ) {
89  Start[iPar] = parResol[iPar];
90  Step[iPar] = parSet[iPar].step;
91  Mini[iPar] = parSet[iPar].mini;
92  Maxi[iPar] = parSet[iPar].maxi;
93  ind[iPar] = parResolOrder[iPar];
94  parname[iPar] = parSet[iPar].name;
95  }
96  }
97 };
98 
99 // SLIMMED VERSION
100 // The full set of scaleFunction developed in the past can be found at
101 // https://github.com/scasasso/cmssw/blob/test_binned_function/MuonAnalysis/MomentumScaleCalibration/interface/Functions.h#L100-L5177
102 
103 template <class T> inline scaleFunctionBase<T>::~scaleFunctionBase() { } // defined even though it's pure virtual; should be faster this way.
104 
105 // No scale
106 // --------
107 template <class T>
108 class scaleFunctionType0 : public scaleFunctionBase<T> {
109 public:
110  scaleFunctionType0() {
111  // One of the two is required. This follows from when templates are used by the compiler and the names lookup rules in c++.
112  // scaleFunctionBase<T>::parNum_ = 0;
113  this->parNum_ = 0;
114  }
115  double scale(const double & pt, const double & eta, const double & phi, const int chg, const T & parScale) const override { return pt; }
116  void resetParameters(std::vector<double> * scaleVec) const override {}
117  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const T & parScale, const std::vector<int> & parScaleOrder, const int muonType) override {}
118 };
119 
120 //
121 // Curvature: (linear eta + sinusoidal in phi (both in 5 eta bins)) * global scale
122 // ------------------------------------------------------------
123 template <class T>
124 class scaleFunctionType50 : public scaleFunctionBase<T> {
125 public:
126  scaleFunctionType50() { this->parNum_ = 27; }
127  double scale(const double & pt, const double & eta, const double & phi, const int chg, const T & parScale) const override {
128  double ampl(0), phase(0), twist(0), ampl2(0), freq2(0), phase2(0);
129 
130 // very bwd bin
131  if ( eta < parScale[4] ) {
132  ampl = parScale[1]; phase = parScale[2]; ampl2 = parScale[21]; freq2 = parScale[22]; phase2 = parScale[23];
133  twist = parScale[3]*(eta-parScale[4])+parScale[7]*(parScale[4]-parScale[8])+parScale[11]*parScale[8];
134 // bwd bin
135  } else if ( parScale[4] <= eta && eta < parScale[8] ) {
136  ampl = parScale[5]; phase = parScale[6];
137  twist = parScale[7]*(eta-parScale[8])+parScale[11]*parScale[8] ;
138 // barrel bin
139  } else if ( parScale[8] <= eta && eta < parScale[12] ) {
140  ampl = parScale[9]; phase = parScale[10];
141  twist = parScale[11]*eta;
142 // fwd bin
143  } else if ( parScale[12] <= eta && eta < parScale[16] ) {
144  ampl = parScale[13]; phase = parScale[14];
145  twist = parScale[15]*(eta-parScale[12])+parScale[11]*parScale[12];
146 // very fwd bin
147  } else if ( parScale[16] < eta ) {
148  ampl = parScale[17]; phase = parScale[18]; ampl2 = parScale[24]; freq2 = parScale[25]; phase2 = parScale[26];
149  twist = parScale[19]*(eta-parScale[16])+parScale[15]*(parScale[16]-parScale[12])+parScale[11]*parScale[12];
150  }
151 
152 // apply the correction
153  double curv = (1.+parScale[0])*((double)chg/pt
154  -twist
155  -ampl*sin(phi+phase)
156  -ampl2*sin((int)freq2*phi+phase2)
157  -0.5*parScale[20]);
158  return 1./((double)chg*curv);
159  }
160  // Fill the scaleVec with neutral parameters
161  void resetParameters(std::vector<double> * scaleVec) const override {
162  // scaleVec->push_back(1);
163  for( int i=0; i<this->parNum_; ++i ) {
164  scaleVec->push_back(0);
165  }
166  }
167  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
168  TString* parname, const T & parScale, const std::vector<int> & parScaleOrder, const int muonType) override {
169 
170  double thisStep[] = {0.000001,
171  0.000001, 0.01, 0.000001, 0.01,
172  0.000001, 0.01, 0.000001, 0.01,
173  0.000001, 0.01, 0.000001, 0.01,
174  0.000001, 0.01, 0.000001, 0.01,
175  0.000001, 0.01, 0.000001,
176  0.000001,
177  0.000001, 0.01, 0.01,
178  0.000001, 0.01, 0.01};
179 
180  TString thisParName[] = { "Curv global scale"
181  , "Phi ampl eta vbwd" , "Phi phase eta vbwd" , "Twist eta vbwd", "vbwd/bwd bndry"
182  , "Phi ampl eta bwd" , "Phi phase eta bwd" , "Twist eta bwd", "bwd/bar bndry"
183  , "Phi ampl eta bar" , "Phi phase eta bar" , "Twist eta bar", "bar/fwd bndry"
184  , "Phi ampl eta fwd" , "Phi phase eta fwd" , "Twist eta fwd", "fwd/vfwd bndry"
185  , "Phi ampl eta vfwd" , "Phi phase eta vfwd" , "Twist eta vfwd"
186  , "Charge depend bias"
187  , "Phi ampl eta vbwd (2nd harmon.)", "Phi freq. eta vbwd (2nd harmon.)", "Phi phase eta vbwd (2nd harmon.)"
188  , "Phi ampl eta vfwd (2nd harmon.)", "Phi freq. eta vfwd (2nd harmon.)", "Phi phase eta vfwd (2nd harmon.)"};
189 
190  if( muonType == 1 ) {
191  double thisMini[] = {-0.1, -0.3, -3.1416, -0.3, -2.6, -0.3, -3.1416, -0.3, -2.6, -0.3, -3.1416, -0.3, 0., -0.3, -3.1416, -0.3, 0. , -0.3, -3.1416, -0.3, -0.1, -0.3, 1., -3.1416, -0.3, 1., -3.1416};
192  double thisMaxi[] = { 0.1, 0.3, 3.1416, 0.3, 0. , 0.3, 3.1416, 0.3, 0. , 0.3, 3.1416, 0.3, 2.6 , 0.3, 3.1416, 0.3, 2.6 , 0.3, 3.1416, 0.3, 0.1, 0.3, 5., 3.1416, 0.3, 5., 3.1416};
193  this->setPar( Start, Step, Mini, Maxi, ind, parname, parScale, parScaleOrder, thisStep, thisMini, thisMaxi, thisParName );
194  } else {
195  double thisMini[] = {-0.1, -0.3, -3.1416, -0.3, -2.6, -0.3, -3.1416, -0.3, -2.6, -0.3, -3.1416, -0.3, 0., -0.3, -3.1416, -0.3, 0. , -0.3, -3.1416, -0.3, -0.1, -0.3, 1., -3.1416, -0.3, 1., -3.1416};
196  double thisMaxi[] = { 0.1, 0.3, 3.1416, 0.3, 0. , 0.3, 3.1416, 0.3, 0. , 0.3, 3.1416, 0.3, 2.6 , 0.3, 3.1416, 0.3, 2.6 , 0.3, 3.1416, 0.3, 0.1, 0.3, 5., 3.1416, 0.3, 5., 3.1416};
197  this->setPar( Start, Step, Mini, Maxi, ind, parname, parScale, parScaleOrder, thisStep, thisMini, thisMaxi, thisParName );
198  }
199  }
200  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
201  const T & parScale, const std::vector<int> & parScaleOrder,
202  const std::vector<double> & parStep,
203  const std::vector<double> & parMin,
204  const std::vector<double> & parMax,
205  const int muonType) override
206  {
207  if( (int(parStep.size()) != this->parNum_) || (int(parMin.size()) != this->parNum_) || (int(parMax.size()) != this->parNum_) ) {
208  std::cout << "Error: par step or min or max do not match with number of parameters" << std::endl;
209  std::cout << "parNum = " << this->parNum_ << std::endl;
210  std::cout << "parStep.size() = " << parStep.size() << std::endl;
211  std::cout << "parMin.size() = " << parMin.size() << std::endl;
212  std::cout << "parMax.size() = " << parMax.size() << std::endl;
213  exit(1);
214  }
215  std::vector<ParameterSet> parSet(this->parNum_);
216  // name, step, mini, maxi
217  parSet[0] = ParameterSet( "Curv global scale", parStep[0], parMin[0], parMax[0] );
218  parSet[1] = ParameterSet( "Phi ampl vbwd", parStep[1], parMin[1], parMax[1] );
219  parSet[2] = ParameterSet( "Phi phase vbwd", parStep[2], parMin[2], parMax[2] );
220  parSet[3] = ParameterSet( "Twist vbwd", parStep[3], parMin[3], parMax[3] );
221  parSet[4] = ParameterSet( "vbwd/bwd bndry", parStep[4], parMin[4], parMax[4] );
222  parSet[5] = ParameterSet( "Phi ampl bwd", parStep[5], parMin[5], parMax[5] );
223  parSet[6] = ParameterSet( "Phi phase bwd", parStep[6], parMin[6], parMax[6] );
224  parSet[7] = ParameterSet( "Twist bwd", parStep[7], parMin[7], parMax[7] );
225  parSet[8] = ParameterSet( "bwd/bar bndry", parStep[8], parMin[8], parMax[8] );
226  parSet[9] = ParameterSet( "Phi ampl bar", parStep[9], parMin[9], parMax[9] );
227  parSet[10] = ParameterSet( "Phi phase bar", parStep[10],parMin[10],parMax[10] );
228  parSet[11] = ParameterSet( "Twist bar", parStep[11],parMin[11],parMax[11] );
229  parSet[12] = ParameterSet( "bar/fwd bndry", parStep[12],parMin[12],parMax[12] );
230  parSet[13] = ParameterSet( "Phi ampl fwd", parStep[13],parMin[13],parMax[13] );
231  parSet[14] = ParameterSet( "Phi phase fwd", parStep[14],parMin[14],parMax[14] );
232  parSet[15] = ParameterSet( "Twist fwd", parStep[15],parMin[15],parMax[15] );
233  parSet[16] = ParameterSet( "fwd/vfwd bndry", parStep[16],parMin[16],parMax[16] );
234  parSet[17] = ParameterSet( "Phi ampl vfwd", parStep[17],parMin[17],parMax[17] );
235  parSet[18] = ParameterSet( "Phi phase vfwd", parStep[18],parMin[18],parMax[18] );
236  parSet[19] = ParameterSet( "Twist vfwd", parStep[19],parMin[19],parMax[19] );
237  parSet[20] = ParameterSet( "Charge depend bias", parStep[20],parMin[20],parMax[20] );
238  parSet[21] = ParameterSet( "Phi ampl vbwd (2nd)", parStep[21],parMin[21],parMax[21] );
239  parSet[22] = ParameterSet( "Phi freq vbwd (2nd)", parStep[22],parMin[22],parMax[22] );
240  parSet[23] = ParameterSet( "Phi phase vbwd (2nd)",parStep[23],parMin[23],parMax[23] );
241  parSet[24] = ParameterSet( "Phi ampl vfwd (2nd)", parStep[24],parMin[24],parMax[24] );
242  parSet[25] = ParameterSet( "Phi freq vfwd (2nd)", parStep[25],parMin[25],parMax[25] );
243  parSet[26] = ParameterSet( "Phi phase vfwd (2nd)",parStep[26],parMin[26],parMax[26] );
244 
245 
246  std::cout << "setting parameters" << std::endl;
247  for( int i=0; i<this->parNum_; ++i ) {
248  std::cout << "parStep["<<i<<"] = " << parStep[i]
249  << ", parMin["<<i<<"] = " << parMin[i]
250  << ", parMax["<<i<<"] = " << parMin[i] << std::endl;
251  }
252  this->setPar( Start, Step, Mini, Maxi, ind, parname, parScale, parScaleOrder, parSet );
253  }
254 
255 };
256 
257 //
258 // Curvature: binned function w/ constraints on eta rings
259 // good results on Legacy 2011 MC
260 // phi bins = [8,8,8,8,8,8,8] in 7 eta bins
261 // ------------------------------------------------------------
262 template <class T>
263 class scaleFunctionType64 : public scaleFunctionBase<T> {
264 public:
265  scaleFunctionType64() {
266  this->parNum_ = 50;
267  }
268  double scale(const double & pt, const double & eta, const double & phi, const int chg, const T & parScale) const override {
269  double deltaK(0);
270  double p0 = parScale[0];
271 
272  if ( eta<-2.1 && eta>-2.4 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[1];
273  else if ( eta<-2.1 && eta>-2.4 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[2];
274  else if ( eta<-2.1 && eta>-2.4 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[3];
275  else if ( eta<-2.1 && eta>-2.4 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[4];
276  else if ( eta<-2.1 && eta>-2.4 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[5];
277  else if ( eta<-2.1 && eta>-2.4 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[6];
278  else if ( eta<-2.1 && eta>-2.4 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[7];
279  else if ( eta<-2.1 && eta>-2.4 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[1]+parScale[2]+parScale[3]+parScale[4]+parScale[5]+parScale[6]+parScale[7]);
280 
281  else if ( eta<-1.5 && eta>-2.1 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[8];
282  else if ( eta<-1.5 && eta>-2.1 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[9];
283  else if ( eta<-1.5 && eta>-2.1 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[10];
284  else if ( eta<-1.5 && eta>-2.1 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[11];
285  else if ( eta<-1.5 && eta>-2.1 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[12];
286  else if ( eta<-1.5 && eta>-2.1 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[13];
287  else if ( eta<-1.5 && eta>-2.1 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[14];
288  else if ( eta<-1.5 && eta>-2.1 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[8]+parScale[9]+parScale[10]+parScale[11]+parScale[12]+parScale[13]+parScale[14]);
289 
290  else if ( eta<-0.9 && eta>-1.5 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[15];
291  else if ( eta<-0.9 && eta>-1.5 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[16];
292  else if ( eta<-0.9 && eta>-1.5 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[17];
293  else if ( eta<-0.9 && eta>-1.5 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[18];
294  else if ( eta<-0.9 && eta>-1.5 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[19];
295  else if ( eta<-0.9 && eta>-1.5 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[20];
296  else if ( eta<-0.9 && eta>-1.5 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[21];
297  else if ( eta<-0.9 && eta>-1.5 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[15]+parScale[16]+parScale[17]+parScale[18]+parScale[19]+parScale[20]+parScale[21]);
298 
299  else if ( eta<0.9 && eta>-0.9 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[22];
300  else if ( eta<0.9 && eta>-0.9 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[23];
301  else if ( eta<0.9 && eta>-0.9 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[24];
302  else if ( eta<0.9 && eta>-0.9 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[25];
303  else if ( eta<0.9 && eta>-0.9 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[26];
304  else if ( eta<0.9 && eta>-0.9 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[27];
305  else if ( eta<0.9 && eta>-0.9 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[28];
306  else if ( eta<0.9 && eta>-0.9 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[22]+parScale[23]+parScale[24]+parScale[25]+parScale[26]+parScale[27]+parScale[28]);
307 
308  else if ( eta<1.5 && eta>0.9 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[29];
309  else if ( eta<1.5 && eta>0.9 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[30];
310  else if ( eta<1.5 && eta>0.9 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[31];
311  else if ( eta<1.5 && eta>0.9 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[32];
312  else if ( eta<1.5 && eta>0.9 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[33];
313  else if ( eta<1.5 && eta>0.9 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[34];
314  else if ( eta<1.5 && eta>0.9 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[35];
315  else if ( eta<1.5 && eta>0.9 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[29]+parScale[30]+parScale[31]+parScale[32]+parScale[33]+parScale[34]+parScale[35]);
316 
317  else if ( eta<2.1 && eta>1.5 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[36];
318  else if ( eta<2.1 && eta>1.5 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[37];
319  else if ( eta<2.1 && eta>1.5 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[38];
320  else if ( eta<2.1 && eta>1.5 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[39];
321  else if ( eta<2.1 && eta>1.5 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[40];
322  else if ( eta<2.1 && eta>1.5 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[41];
323  else if ( eta<2.1 && eta>1.5 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[42];
324  else if ( eta<2.1 && eta>1.5 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[36]+parScale[37]+parScale[38]+parScale[39]+parScale[40]+parScale[41]+parScale[42]);
325 
326  else if ( eta<2.4 && eta>2.1 && phi<-2.35619449019 && phi>-3.14159265359 ) deltaK = parScale[43];
327  else if ( eta<2.4 && eta>2.1 && phi<-1.57079632679 && phi>-2.35619449019 ) deltaK = parScale[44];
328  else if ( eta<2.4 && eta>2.1 && phi<-0.785398163397 && phi>-1.57079632679 ) deltaK = parScale[45];
329  else if ( eta<2.4 && eta>2.1 && phi<0.0 && phi>-0.785398163397 ) deltaK = parScale[46];
330  else if ( eta<2.4 && eta>2.1 && phi<0.785398163397 && phi>0.0 ) deltaK = parScale[47];
331  else if ( eta<2.4 && eta>2.1 && phi<1.57079632679 && phi>0.785398163397 ) deltaK = parScale[48];
332  else if ( eta<2.4 && eta>2.1 && phi<2.35619449019 && phi>1.57079632679 ) deltaK = parScale[49];
333  else if ( eta<2.4 && eta>2.1 && phi<3.14159265359 && phi>2.35619449019 ) deltaK = -(parScale[43]+parScale[44]+parScale[45]+parScale[46]+parScale[47]+parScale[48]+parScale[49]);
334  else {
335  std::cout << "This should really not happen, this muon has eta = " << eta << "and phi = " << phi << std::endl;
336  exit(1);
337  }
338 
339  // apply the correction
340  double curv = (double)chg/pt;
341  return 1./((double)chg*(1+p0)*(curv+deltaK));
342  }
343  // Fill the scaleVec with neutral parameters
344  void resetParameters(std::vector<double> * scaleVec) const override {
345  // scaleVec->push_back(1);
346  for( int i=0; i<this->parNum_; ++i ) {
347  scaleVec->push_back(0);
348  }
349  }
350  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
351  TString* parname, const T & parScale, const std::vector<int> & parScaleOrder, const int muonType) override {
352 
353  double thisStep[] = {
354  0.000001,
355  0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001
356  };
357 
358  TString thisParName[] = {
359  "Curv global scale",
360  "deltaK bin1", "deltaK bin2", "deltaK bin3", "deltaK bin4", "deltaK bin5", "deltaK bin6", "deltaK bin7", "deltaK bin8", "deltaK bin9", "deltaK bin10", "deltaK bin11", "deltaK bin12", "deltaK bin13", "deltaK bin14", "deltaK bin15", "deltaK bin16", "deltaK bin17", "deltaK bin18", "deltaK bin19", "deltaK bin20", "deltaK bin21", "deltaK bin22", "deltaK bin23", "deltaK bin24", "deltaK bin25", "deltaK bin26", "deltaK bin27", "deltaK bin28", "deltaK bin29", "deltaK bin30", "deltaK bin31", "deltaK bin32", "deltaK bin33", "deltaK bin34", "deltaK bin35", "deltaK bin36", "deltaK bin37", "deltaK bin38", "deltaK bin39", "deltaK bin40", "deltaK bin41", "deltaK bin42", "deltaK bin43", "deltaK bin44", "deltaK bin45", "deltaK bin46", "deltaK bin47", "deltaK bin48", "deltaK bin49"
361  };
362 
363  if( muonType == 1 ) {
364  double thisMini[] = {
365  -0.1,
366  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
367  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
368  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
369  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
370  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
371  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
372  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005
373  };
374  double thisMaxi[] = {
375  0.1,
376  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
377  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
378  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
379  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
380  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
381  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
382  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005
383  };
384  this->setPar( Start, Step, Mini, Maxi, ind, parname, parScale, parScaleOrder, thisStep, thisMini, thisMaxi, thisParName );
385  } else {
386  double thisMini[] = {
387  -0.1,
388  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
389  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
390  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
391  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
392  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
393  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005,
394  -0.005, -0.005, -0.005, -0.005, -0.005, -0.005, -0.005
395  };
396  double thisMaxi[] = {
397  0.1,
398  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
399  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
400  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
401  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
402  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
403  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005,
404  0.005, 0.005, 0.005, 0.005, 0.005, 0.005, 0.005
405  };
406  this->setPar( Start, Step, Mini, Maxi, ind, parname, parScale, parScaleOrder, thisStep, thisMini, thisMaxi, thisParName );
407  }
408  }
409  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
410  const T & parScale, const std::vector<int> & parScaleOrder,
411  const std::vector<double> & parStep,
412  const std::vector<double> & parMin,
413  const std::vector<double> & parMax,
414  const int muonType) override
415  {
416  if( (int(parStep.size()) != this->parNum_) || (int(parMin.size()) != this->parNum_) || (int(parMax.size()) != this->parNum_) ) {
417  std::cout << "Error: par step or min or max do not match with number of parameters" << std::endl;
418  std::cout << "parNum = " << this->parNum_ << std::endl;
419  std::cout << "parStep.size() = " << parStep.size() << std::endl;
420  std::cout << "parMin.size() = " << parMin.size() << std::endl;
421  std::cout << "parMax.size() = " << parMax.size() << std::endl;
422  exit(1);
423  }
424  std::vector<ParameterSet> parSet(this->parNum_);
425  // name, step, mini, maxi
426  parSet[0] = ParameterSet( "Curv global scale", parStep[0], parMin[0], parMax[0] );
427  parSet[1] = ParameterSet( "deltaK bin1", parStep[1], parMin[1], parMax[1] );
428  parSet[2] = ParameterSet( "deltaK bin2", parStep[2], parMin[2], parMax[2] );
429  parSet[3] = ParameterSet( "deltaK bin3", parStep[3], parMin[3], parMax[3] );
430  parSet[4] = ParameterSet( "deltaK bin4", parStep[4], parMin[4], parMax[4] );
431  parSet[5] = ParameterSet( "deltaK bin5", parStep[5], parMin[5], parMax[5] );
432  parSet[6] = ParameterSet( "deltaK bin6", parStep[6], parMin[6], parMax[6] );
433  parSet[7] = ParameterSet( "deltaK bin7", parStep[7], parMin[7], parMax[7] );
434  parSet[8] = ParameterSet( "deltaK bin9", parStep[8], parMin[8], parMax[8] );
435  parSet[9] = ParameterSet( "deltaK bin10", parStep[9], parMin[9], parMax[9] );
436  parSet[10] = ParameterSet( "deltaK bin11", parStep[10], parMin[10], parMax[10] );
437  parSet[11] = ParameterSet( "deltaK bin12", parStep[11], parMin[11], parMax[11] );
438  parSet[12] = ParameterSet( "deltaK bin13", parStep[12], parMin[12], parMax[12] );
439  parSet[13] = ParameterSet( "deltaK bin14", parStep[13], parMin[13], parMax[13] );
440  parSet[14] = ParameterSet( "deltaK bin15", parStep[14], parMin[14], parMax[14] );
441  parSet[15] = ParameterSet( "deltaK bin17", parStep[15], parMin[15], parMax[15] );
442  parSet[16] = ParameterSet( "deltaK bin18", parStep[16], parMin[16], parMax[16] );
443  parSet[17] = ParameterSet( "deltaK bin19", parStep[17], parMin[17], parMax[17] );
444  parSet[18] = ParameterSet( "deltaK bin20", parStep[18], parMin[18], parMax[18] );
445  parSet[19] = ParameterSet( "deltaK bin21", parStep[19], parMin[19], parMax[19] );
446  parSet[20] = ParameterSet( "deltaK bin22", parStep[20], parMin[20], parMax[20] );
447  parSet[21] = ParameterSet( "deltaK bin23", parStep[21], parMin[21], parMax[21] );
448  parSet[22] = ParameterSet( "deltaK bin25", parStep[22], parMin[22], parMax[22] );
449  parSet[23] = ParameterSet( "deltaK bin26", parStep[23], parMin[23], parMax[23] );
450  parSet[24] = ParameterSet( "deltaK bin27", parStep[24], parMin[24], parMax[24] );
451  parSet[25] = ParameterSet( "deltaK bin28", parStep[25], parMin[25], parMax[25] );
452  parSet[26] = ParameterSet( "deltaK bin29", parStep[26], parMin[26], parMax[26] );
453  parSet[27] = ParameterSet( "deltaK bin30", parStep[27], parMin[27], parMax[27] );
454  parSet[28] = ParameterSet( "deltaK bin31", parStep[28], parMin[28], parMax[28] );
455  parSet[29] = ParameterSet( "deltaK bin33", parStep[29], parMin[29], parMax[29] );
456  parSet[30] = ParameterSet( "deltaK bin34", parStep[30], parMin[30], parMax[30] );
457  parSet[31] = ParameterSet( "deltaK bin35", parStep[31], parMin[31], parMax[31] );
458  parSet[32] = ParameterSet( "deltaK bin36", parStep[32], parMin[32], parMax[32] );
459  parSet[33] = ParameterSet( "deltaK bin37", parStep[33], parMin[33], parMax[33] );
460  parSet[34] = ParameterSet( "deltaK bin38", parStep[34], parMin[34], parMax[34] );
461  parSet[35] = ParameterSet( "deltaK bin39", parStep[35], parMin[35], parMax[35] );
462  parSet[36] = ParameterSet( "deltaK bin41", parStep[36], parMin[36], parMax[36] );
463  parSet[37] = ParameterSet( "deltaK bin42", parStep[37], parMin[37], parMax[37] );
464  parSet[38] = ParameterSet( "deltaK bin43", parStep[38], parMin[38], parMax[38] );
465  parSet[39] = ParameterSet( "deltaK bin44", parStep[39], parMin[39], parMax[39] );
466  parSet[40] = ParameterSet( "deltaK bin45", parStep[40], parMin[40], parMax[40] );
467  parSet[41] = ParameterSet( "deltaK bin46", parStep[41], parMin[41], parMax[41] );
468  parSet[42] = ParameterSet( "deltaK bin47", parStep[42], parMin[42], parMax[42] );
469  parSet[43] = ParameterSet( "deltaK bin49", parStep[43], parMin[43], parMax[43] );
470  parSet[44] = ParameterSet( "deltaK bin50", parStep[44], parMin[44], parMax[44] );
471  parSet[45] = ParameterSet( "deltaK bin51", parStep[45], parMin[45], parMax[45] );
472  parSet[46] = ParameterSet( "deltaK bin52", parStep[46], parMin[46], parMax[46] );
473  parSet[47] = ParameterSet( "deltaK bin53", parStep[47], parMin[47], parMax[47] );
474  parSet[48] = ParameterSet( "deltaK bin54", parStep[48], parMin[48], parMax[48] );
475  parSet[49] = ParameterSet( "deltaK bin55", parStep[49], parMin[49], parMax[49] );
476 
477  std::cout << "setting parameters" << std::endl;
478  for( int i=0; i<this->parNum_; ++i ) {
479  std::cout << "parStep["<<i<<"] = " << parStep[i]
480  << ", parMin["<<i<<"] = " << parMin[i]
481  << ", parMax["<<i<<"] = " << parMin[i] << std::endl;
482  }
483  this->setPar( Start, Step, Mini, Maxi, ind, parname, parScale, parScaleOrder, parSet );
484  }
485 
486 };
487 
488 
490 scaleFunctionBase<double * > * scaleFunctionService( const int identifier );
491 
493 scaleFunctionBase<std::vector<double> > * scaleFunctionVecService( const int identifier );
494 
495 // -------------- //
496 // Smear functors //
497 // -------------- //
498 
499 class smearFunctionBase {
500  public:
501  virtual void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) = 0;
502  smearFunctionBase() {
503  cotgth_ = 0.;
504  gRandom_ = new TRandom();
505  }
506  virtual ~smearFunctionBase() = 0;
507 protected:
508  void smearEta(double & eta) {
509  double theta;
510  if (cotgth_!=0) {
511  theta = atan(1/cotgth_);
512  } else {
513  theta = TMath::Pi()/2;
514  }
515  if (theta<0) theta += TMath::Pi();
516  eta = -log(tan(theta/2));
517  }
518  double cotgth_;
519  TRandom * gRandom_;
520 };
521 inline smearFunctionBase::~smearFunctionBase() { } // defined even though it's pure virtual; should be faster this way.
522 
523 // No smearing
524 // -----------
525 class smearFunctionType0 : public smearFunctionBase {
526  public:
527  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override { }
528 };
529 // The 3 parameters of smearType1 are: pt dependence of pt smear, phi smear and
530 // cotgtheta smear.
531 class smearFunctionType1 : public smearFunctionBase {
532  public:
533  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override {
534  pt = pt*(1.0+y[0]*parSmear[0]*pt);
535  phi = phi*(1.0+y[1]*parSmear[1]);
536  double tmp = 2*atan(exp(-eta));
537  cotgth_ = cos(tmp)/sin(tmp)*(1.0+y[2]*parSmear[2]);
538  smearEta(eta);
539  }
540 };
541 
542 class smearFunctionType2 : public smearFunctionBase {
543  public:
544  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override {
545  pt = pt*(1.0+y[0]*parSmear[0]*pt+y[1]*parSmear[1]*std::fabs(eta));
546  phi = phi*(1.0+y[2]*parSmear[2]);
547  double tmp = 2*atan(exp(-eta));
548  cotgth_ = cos(tmp)/sin(tmp)*(1.0+y[3]*parSmear[3]);
549  smearEta(eta);
550  }
551 };
552 
553 class smearFunctionType3 : public smearFunctionBase {
554  public:
555  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override {
556  pt = pt*(1.0+y[0]*parSmear[0]*pt+y[1]*parSmear[1]*std::fabs(eta));
557  phi = phi*(1.0+y[2]*parSmear[2]);
558  double tmp = 2*atan(exp(-eta));
559  cotgth_ = cos(tmp)/sin(tmp)*(1.0+y[3]*parSmear[3]+y[4]*parSmear[4]*std::fabs(eta));
560  smearEta(eta);
561  }
562 };
563 // The six parameters of SmearType=4 are respectively:
564 // Pt dep. of Pt res., |eta| dep. of Pt res., Phi res., |eta| res.,
565 // |eta| dep. of |eta| res., Pt^2 dep. of Pt res.
566 class smearFunctionType4 : public smearFunctionBase {
567  public:
568  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override {
569  pt = pt*(1.0+y[0]*parSmear[0]*pt+y[1]*parSmear[1]*std::fabs(eta)+y[5]*parSmear[5]*pow(pt,2));
570  phi = phi*(1.0+y[2]*parSmear[2]);
571  double tmp = 2*atan(exp(-eta));
572  cotgth_ = cos(tmp)/sin(tmp)*(1.0+y[3]*parSmear[3]+y[4]*parSmear[4]*std::fabs(eta));
573  smearEta(eta);
574  }
575 };
576 
577 class smearFunctionType5 : public smearFunctionBase {
578  public:
579  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override {
580  pt = pt*(1.0+y[0]*parSmear[0]*pt+y[1]*parSmear[1]*std::fabs(eta)+y[5]*parSmear[5]*pow(pt,2));
581  phi = phi*(1.0+y[2]*parSmear[2]+y[6]*parSmear[6]*pt);
582  double tmp = 2*atan(exp(-eta));
583  cotgth_ = cos(tmp)/sin(tmp)*(1.0+y[3]*parSmear[3]+y[4]*parSmear[4]*std::fabs(eta));
584  smearEta(eta);
585  }
586 };
587 
588 //Smearing for MC correction based on the resolution function Type 15 for misaligned MC
589 class smearFunctionType6 : public smearFunctionBase {
590  public:
591  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override {
592  double sigmaSmear = 0;
593  double sigmaPtAl = 0;
594  double sigmaPtMisal = 0;
595  double ptPart = parSmear[0] + parSmear[1]*1/pt + pt*parSmear[2];
596  double fabsEta = std::fabs(eta);
597 
598  sigmaPtAl = parSmear[14]*etaByPoints(eta, parSmear[15]);
599 
600  if (std::fabs(eta)<=1.4){
601  sigmaPtMisal = ptPart + parSmear[3] + parSmear[4]*std::fabs(eta) + parSmear[5]*eta*eta;
602  sigmaSmear = sqrt(std::fabs(pow(sigmaPtMisal,2)-pow(sigmaPtAl,2)));
603  pt = pt*gRandom_->Gaus(1,sigmaSmear);
604  }
605  else if (eta>1.4){//eta in right endcap
606  double par = parSmear[3] + parSmear[4]*1.4 + parSmear[5]*1.4*1.4 - (parSmear[6] + parSmear[7]*(1.4-parSmear[8]) + parSmear[9]*(1.4-parSmear[8])*(1.4-parSmear[8]));
607  sigmaPtMisal = par + ptPart + parSmear[6] + parSmear[7]*std::fabs((fabsEta-parSmear[8])) + parSmear[9]*(fabsEta-parSmear[8])*(fabsEta-parSmear[8]);
608  sigmaSmear = sqrt(std::fabs(pow(sigmaPtMisal,2)-pow(sigmaPtAl,2)));
609  pt = pt*gRandom_->Gaus(1,sigmaSmear);
610  }
611  else{//eta in left endcap
612  double par = parSmear[3] + parSmear[4]*1.4 + parSmear[5]*1.4*1.4 - (parSmear[10] + parSmear[11]*(1.4-parSmear[12]) + parSmear[13]*(1.4-parSmear[12])*(1.4-parSmear[12]));
613  sigmaPtMisal = par + ptPart + parSmear[10] + parSmear[11]*std::fabs((fabsEta-parSmear[12])) + parSmear[13]*(fabsEta-parSmear[12])*(fabsEta-parSmear[12]);
614  sigmaSmear = sqrt(std::fabs(pow(sigmaPtMisal,2)-pow(sigmaPtAl,2)));
615  pt = pt*gRandom_->Gaus(1,sigmaSmear);
616  }
617  }
618  protected:
623  double etaByPoints(const double & inEta, const double & border) {
624  Double_t eta = std::fabs(inEta);
625  if( 0. <= eta && eta <= 0.2 ) return 0.00942984;
626  else if( 0.2 < eta && eta <= 0.4 ) return 0.0104489;
627  else if( 0.4 < eta && eta <= 0.6 ) return 0.0110521;
628  else if( 0.6 < eta && eta <= 0.8 ) return 0.0117338;
629  else if( 0.8 < eta && eta <= 1.0 ) return 0.0138142;
630  else if( 1.0 < eta && eta <= 1.2 ) return 0.0165826;
631  else if( 1.2 < eta && eta <= 1.4 ) return 0.0183663;
632  else if( 1.4 < eta && eta <= 1.6 ) return 0.0169904;
633  else if( 1.6 < eta && eta <= 1.8 ) return 0.0173289;
634  else if( 1.8 < eta && eta <= 2.0 ) return 0.0205821;
635  else if( 2.0 < eta && eta <= 2.2 ) return 0.0250032;
636  else if( 2.2 < eta && eta <= 2.4 ) return 0.0339477;
637  // ATTENTION: This point has a big error and it is very displaced from the rest of the distribution.
638  else if( 2.4 < eta && eta <= 2.6 ) return border;
639  return ( 0. );
640  }
641 };
642 
643 class smearFunctionType7 : public smearFunctionBase
644 {
645  public:
646  void smear(double & pt, double & eta, double & phi, const double * y, const std::vector<double> & parSmear) override
647  {
648  double sigmaSquared = sigmaPtDiff.squaredDiff(eta);
649  TF1 G("G", "[0]*exp(-0.5*pow(x,2)/[1])", -5., 5.);
650  double norm = 1/(sqrt(2*TMath::Pi()*sigmaSquared));
651  G.SetParameter (0,norm);
652  G.SetParameter (1,sigmaSquared);
653  pt = pt*(1-G.GetRandom());
654  }
655  SigmaPtDiff sigmaPtDiff;
656 };
657 
659 smearFunctionBase * smearFunctionService( const int identifier );
660 
661 // // Defined globally...
662 // static smearFunctionBase * smearFunctionArray[] = {
663 // new smearFunctionType0,
664 // new smearFunctionType1,
665 // new smearFunctionType2,
666 // new smearFunctionType3,
667 // new smearFunctionType4,
668 // new smearFunctionType5
669 // };
670 
675 template <class T>
676 class resolutionFunctionBase {
677  public:
678  virtual double sigmaPt(const double & pt, const double & eta, const T & parval) = 0;
679  virtual double sigmaPtError(const double & pt, const double & eta, const T & parval, const T & parError)
680  {
681  return 0.;
682  }
683  virtual double sigmaPhi(const double & pt, const double & eta, const T & parval) = 0;
684  virtual double sigmaCotgTh(const double & pt, const double & eta, const T & parval) = 0;
685  virtual double covPt1Pt2(const double & pt1, const double & eta1, const double & pt2, const double & eta2, const T & parval)
686  {
687  return 0.;
688  }
689  resolutionFunctionBase() {}
690  virtual ~resolutionFunctionBase() = 0;
692  virtual void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
693  const T & parResol, const std::vector<int> & parResolOrder, const int muonType) {};
694  virtual void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
695  const T & parResol, const std::vector<int> & parResolOrder,
696  const std::vector<double> & parStep,
697  const std::vector<double> & parMin,
698  const std::vector<double> & parMax,
699  const int muonType)
700  {
701  std::cout << "The method setParameters must be implemented for this resolution function" << std::endl;
702  exit(1);
703  }
704  virtual int parNum() const { return parNum_; }
705  protected:
706  int parNum_;
708  virtual void setPar(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
709  TString* parname, const T & parResol, const std::vector<int> & parResolOrder,
710  double* thisStep, double* thisMini, double* thisMaxi, TString* thisParName ) {
711  for( int iPar=0; iPar<this->parNum_; ++iPar ) {
712  Start[iPar] = parResol[iPar];
713  Step[iPar] = thisStep[iPar];
714  Mini[iPar] = thisMini[iPar];
715  Maxi[iPar] = thisMaxi[iPar];
716  ind[iPar] = parResolOrder[iPar];
717  parname[iPar] = thisParName[iPar];
718  }
719  }
720  virtual void setPar(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
721  TString* parname, const T & parResol, const std::vector<int> & parResolOrder, const std::vector<ParameterSet> & parSet ) {
722  if( int(parSet.size()) != this->parNum_ ) {
723  std::cout << "Error: wrong number of parameter initializations = " << parSet.size() << ". Number of parameters is " << this->parNum_ << std::endl;
724  exit(1);
725  }
726  for( int iPar=0; iPar<this->parNum_; ++iPar ) {
727  Start[iPar] = parResol[iPar];
728  Step[iPar] = parSet[iPar].step;
729  Mini[iPar] = parSet[iPar].mini;
730  Maxi[iPar] = parSet[iPar].maxi;
731  ind[iPar] = parResolOrder[iPar];
732  parname[iPar] = parSet[iPar].name;
733  }
734  }
735 };
736 template <class T> inline resolutionFunctionBase<T>::~resolutionFunctionBase() { } // defined even though it's pure virtual; should be faster this way.
737 
738 // SLIMMED VERSION
739 // The full set of resolutionFunction developed in the past can be found at
740 // https://github.com/cms-sw/cmssw/blob/CMSSW_5_3_X/MuonAnalysis/MomentumScaleCalibration/interface/Functions.h#L3082-L5052
741 
742 // null
743 // --------
744 template <class T>
745 class resolutionFunctionType0 : public resolutionFunctionBase<T> {
746 public:
747  resolutionFunctionType0() {
748  // One of the two is required. This follows from when templates are used by the compiler and the names lookup rules in c++.
749  this->parNum_ = 0;
750  }
751  double sigmaPt(const double & pt, const double & eta, const T & parval) override { return 0; }
752  double sigmaCotgTh(const double & pt, const double & eta, const T & parval) override { return 0; }
753  double sigmaPhi(const double & pt, const double & eta, const T & parval) override { return 0.; }
754  // derivatives ---------------
755  double sigmaPtError(const double & pt, const double & eta, const T & parval, const T & parError) override { return 0; }
756  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const T & parResol, const std::vector<int> & parResolOrder, const int muonType) override {}
757  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
758  const T & parResol, const std::vector<int> & parResolOrder,
759  const std::vector<double> & parStep,
760  const std::vector<double> & parMin,
761  const std::vector<double> & parMax,
762  const int muonType) override {}
763 
764 };
765 
766 // Binned in eta to fit the Z (parametrization as linear sum)
767 template <class T>
768 class resolutionFunctionType45 : public resolutionFunctionBase<T> {
769  public:
770  int etaBin(const double & eta)
771  {
772  // 12 bins from -2.4 to 2.4
773  // std::cout << "for eta = " << eta << ", bin = " << bin << std::endl;
774  if( eta < -2.0 ) return 1;
775  if( eta < -1.8 ) return 2;
776  if( eta < -1.6 ) return 3;
777  if( eta < -1.2 ) return 4;
778  if( eta < -0.8 ) return 5;
779  if( eta < 0. ) return 6;
780  if( eta < 0.8 ) return 7;
781  if( eta < 1.2 ) return 8;
782  if( eta < 1.6 ) return 9;
783  if( eta < 1.8 ) return 10;
784  if( eta < 2.0 ) return 11;
785  return 12;
786  }
787 
788  resolutionFunctionType45() { this->parNum_ = 13; }
789 
790  double sigmaPt(const double & pt, const double & eta, const T & parval) override
791  {
792  return (parval[0]*pt + parval[etaBin(eta)]);
793  }
794  double sigmaCotgTh(const double & pt, const double & eta, const T & parval) override {
795  return 0;
796  }
797  double sigmaPhi(const double & pt, const double & eta, const T & parval) override {
798  return 0.;
799  }
800 
801  // derivatives ---------------
802  double sigmaPtError(const double & pt, const double & eta, const T & parval, const T & parError) override
803  {
804  // Use the etaBin function to select the right bin for the parameter
805  return sqrt( pow(pt*parError[0], 2) + pow(parError[etaBin(eta)], 2));
806  }
807 
808  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
809  TString* parname, const T & parResol, const std::vector<int> & parResolOrder,
810  const int muonType) override
811  {
812  std::vector<ParameterSet> parSet(this->parNum_);
813  // name, step, mini, maxi
814  parSet[0] = ParameterSet( "Pt res. sc.", 0.002, -0.1, 0.1 );
815  parSet[1] = ParameterSet( "eta bin 1", 0.00002, -0.01, 0.01 );
816  parSet[2] = ParameterSet( "eta bin 2", 0.00002, -0.01, 0.01 );
817  parSet[3] = ParameterSet( "eta bin 3", 0.00002, -0.01, 0.01 );
818  parSet[4] = ParameterSet( "eta bin 4", 0.00002, -0.01, 0.01 );
819  parSet[5] = ParameterSet( "eta bin 5", 0.00002, -0.01, 0.01 );
820  parSet[6] = ParameterSet( "eta bin 6", 0.00002, -0.01, 0.01 );
821  parSet[7] = ParameterSet( "eta bin 7", 0.00002, -0.01, 0.01 );
822  parSet[8] = ParameterSet( "eta bin 8", 0.00002, -0.01, 0.01 );
823  parSet[9] = ParameterSet( "eta bin 9", 0.00002, -0.01, 0.01 );
824  parSet[10] = ParameterSet( "eta bin 10", 0.00002, -0.01, 0.01 );
825  parSet[11] = ParameterSet( "eta bin 11", 0.00002, -0.01, 0.01 );
826  parSet[12] = ParameterSet( "eta bin 12", 0.00002, -0.01, 0.01 );
827 
828  std::cout << "setting parameters" << std::endl;
829  this->setPar( Start, Step, Mini, Maxi, ind, parname, parResol, parResolOrder, parSet );
830  }
831 
832  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
833  const T & parResol, const std::vector<int> & parResolOrder,
834  const std::vector<double> & parStep,
835  const std::vector<double> & parMin,
836  const std::vector<double> & parMax,
837  const int muonType) override
838  {
839  if( (int(parStep.size()) != this->parNum_) || (int(parMin.size()) != this->parNum_) || (int(parMax.size()) != this->parNum_) ) {
840  std::cout << "Error: par step or min or max do not match with number of parameters" << std::endl;
841  std::cout << "parNum = " << this->parNum_ << std::endl;
842  std::cout << "parStep.size() = " << parStep.size() << std::endl;
843  std::cout << "parMin.size() = " << parMin.size() << std::endl;
844  std::cout << "parMax.size() = " << parMax.size() << std::endl;
845  exit(1);
846  }
847  std::vector<ParameterSet> parSet(this->parNum_);
848  // name, step, mini, maxi
849  parSet[0] = ParameterSet( "Pt res. sc.", parStep[0], parMin[0], parMax[0] );
850  parSet[1] = ParameterSet( "eta bin 1", parStep[1], parMin[1], parMax[1] );
851  parSet[2] = ParameterSet( "eta bin 2", parStep[2], parMin[2], parMax[2] );
852  parSet[3] = ParameterSet( "eta bin 3", parStep[3], parMin[3], parMax[3] );
853  parSet[4] = ParameterSet( "eta bin 4", parStep[4], parMin[4], parMax[4] );
854  parSet[5] = ParameterSet( "eta bin 5", parStep[5], parMin[5], parMax[5] );
855  parSet[6] = ParameterSet( "eta bin 6", parStep[6], parMin[6], parMax[6] );
856  parSet[7] = ParameterSet( "eta bin 7", parStep[7], parMin[7], parMax[7] );
857  parSet[8] = ParameterSet( "eta bin 8", parStep[8], parMin[8], parMax[8] );
858  parSet[9] = ParameterSet( "eta bin 9", parStep[9], parMin[9], parMax[9] );
859  parSet[10] = ParameterSet( "eta bin 10", parStep[10], parMin[10], parMax[10] );
860  parSet[11] = ParameterSet( "eta bin 11", parStep[11], parMin[11], parMax[11] );
861  parSet[12] = ParameterSet( "eta bin 12", parStep[12], parMin[12], parMax[12] );
862 
863  std::cout << "setting parameters" << std::endl;
864  for( int i=0; i<this->parNum_; ++i ) {
865  std::cout << "parStep["<<i<<"] = " << parStep[i]
866  << ", parMin["<<i<<"] = " << parMin[i]
867  << ", parMax["<<i<<"] = " << parMin[i] << std::endl;
868  }
869  this->setPar( Start, Step, Mini, Maxi, ind, parname, parResol, parResolOrder, parSet );
870  }
871 };
872 
873 
874 // Binned in eta to fit the Z (parametrization as sum in quadrature)
875 template <class T>
876 class resolutionFunctionType46 : public resolutionFunctionBase<T> {
877  public:
878  int etaBin(const double & eta)
879  {
880  // 12 bins from -2.4 to 2.4
881  // std::cout << "for eta = " << eta << ", bin = " << bin << std::endl;
882  if( eta < -2.0 ) return 1;
883  if( eta < -1.8 ) return 2;
884  if( eta < -1.6 ) return 3;
885  if( eta < -1.2 ) return 4;
886  if( eta < -0.8 ) return 5;
887  if( eta < 0. ) return 6;
888  if( eta < 0.8 ) return 7;
889  if( eta < 1.2 ) return 8;
890  if( eta < 1.6 ) return 9;
891  if( eta < 1.8 ) return 10;
892  if( eta < 2.0 ) return 11;
893  return 12;
894  }
895 
896  resolutionFunctionType46() { this->parNum_ = 13; }
897 
898  double sigmaPt(const double & pt, const double & eta, const T & parval) override
899  {
900  return sqrt(pow(parval[0]*pt,2) + pow(parval[etaBin(eta)],2));
901  }
902  double sigmaCotgTh(const double & pt, const double & eta, const T & parval) override {
903  return 0;
904  }
905  double sigmaPhi(const double & pt, const double & eta, const T & parval) override {
906  return 0.;
907  }
908 
909  // derivatives ---------------
910  double sigmaPtError(const double & pt, const double & eta, const T & parval, const T & parError) override
911  {
912  // Use the etaByBin function to select the right bin for the parameter
913  double r = sqrt(pow(parval[0]*pt,2) + pow(parval[etaBin(eta)],2));
914  return sqrt( pow(pt*pt*parval[0]*parError[0],2) + pow(parval[etaBin(eta)]*parError[etaBin(eta)],2) )/r;
915  }
916 
917  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
918  TString* parname, const T & parResol, const std::vector<int> & parResolOrder,
919  const int muonType) override
920  {
921  std::vector<ParameterSet> parSet(this->parNum_);
922  // name, step, mini, maxi
923  parSet[0] = ParameterSet( "Pt res. sc.", 0.0002, 0., 0.1 );
924  parSet[1] = ParameterSet( "eta bin 1", 0.00002, 0., 0.01 );
925  parSet[2] = ParameterSet( "eta bin 2", 0.00002, 0., 0.01 );
926  parSet[3] = ParameterSet( "eta bin 3", 0.00002, 0., 0.01 );
927  parSet[4] = ParameterSet( "eta bin 4", 0.00002, 0., 0.01 );
928  parSet[5] = ParameterSet( "eta bin 5", 0.00002, 0., 0.01 );
929  parSet[6] = ParameterSet( "eta bin 6", 0.00002, 0., 0.01 );
930  parSet[7] = ParameterSet( "eta bin 7", 0.00002, 0., 0.01 );
931  parSet[8] = ParameterSet( "eta bin 8", 0.00002, 0., 0.01 );
932  parSet[9] = ParameterSet( "eta bin 9", 0.00002, 0., 0.01 );
933  parSet[10] = ParameterSet( "eta bin 10", 0.00002, 0., 0.01 );
934  parSet[11] = ParameterSet( "eta bin 11", 0.00002, 0., 0.01 );
935  parSet[12] = ParameterSet( "eta bin 12", 0.00002, 0., 0.01 );
936 
937  std::cout << "setting parameters" << std::endl;
938  this->setPar( Start, Step, Mini, Maxi, ind, parname, parResol, parResolOrder, parSet );
939  }
940 
941  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
942  const T & parResol, const std::vector<int> & parResolOrder,
943  const std::vector<double> & parStep,
944  const std::vector<double> & parMin,
945  const std::vector<double> & parMax,
946  const int muonType) override
947  {
948  if( (int(parStep.size()) != this->parNum_) || (int(parMin.size()) != this->parNum_) || (int(parMax.size()) != this->parNum_) ) {
949  std::cout << "Error: par step or min or max do not match with number of parameters" << std::endl;
950  std::cout << "parNum = " << this->parNum_ << std::endl;
951  std::cout << "parStep.size() = " << parStep.size() << std::endl;
952  std::cout << "parMin.size() = " << parMin.size() << std::endl;
953  std::cout << "parMax.size() = " << parMax.size() << std::endl;
954  exit(1);
955  }
956  std::vector<ParameterSet> parSet(this->parNum_);
957  // name, step, mini, maxi
958  parSet[0] = ParameterSet( "Pt res. sc.", parStep[0], parMin[0], parMax[0] );
959  parSet[1] = ParameterSet( "eta bin 1", parStep[1], parMin[1], parMax[1] );
960  parSet[2] = ParameterSet( "eta bin 2", parStep[2], parMin[2], parMax[2] );
961  parSet[3] = ParameterSet( "eta bin 3", parStep[3], parMin[3], parMax[3] );
962  parSet[4] = ParameterSet( "eta bin 4", parStep[4], parMin[4], parMax[4] );
963  parSet[5] = ParameterSet( "eta bin 5", parStep[5], parMin[5], parMax[5] );
964  parSet[6] = ParameterSet( "eta bin 6", parStep[6], parMin[6], parMax[6] );
965  parSet[7] = ParameterSet( "eta bin 7", parStep[7], parMin[7], parMax[7] );
966  parSet[8] = ParameterSet( "eta bin 8", parStep[8], parMin[8], parMax[8] );
967  parSet[9] = ParameterSet( "eta bin 9", parStep[9], parMin[9], parMax[9] );
968  parSet[10] = ParameterSet( "eta bin 10", parStep[10], parMin[10], parMax[10] );
969  parSet[11] = ParameterSet( "eta bin 11", parStep[11], parMin[11], parMax[11] );
970  parSet[12] = ParameterSet( "eta bin 12", parStep[12], parMin[12], parMax[12] );
971 
972  std::cout << "setting parameters" << std::endl;
973  for( int i=0; i<this->parNum_; ++i ) {
974  std::cout << "parStep["<<i<<"] = " << parStep[i]
975  << ", parMin["<<i<<"] = " << parMin[i]
976  << ", parMax["<<i<<"] = " << parMin[i] << std::endl;
977  }
978  this->setPar( Start, Step, Mini, Maxi, ind, parname, parResol, parResolOrder, parSet );
979  }
980 };
981 
982 // Binned in eta to fit the Z (parametrization as sum in quadrature) and including an overall covariance
983 template <class T>
984 class resolutionFunctionType47 : public resolutionFunctionBase<T> {
985  public:
986  int etaBin(const double & eta)
987  {
988  // 12 bins from -2.4 to 2.4
989  // std::cout << "for eta = " << eta << ", bin = " << bin << std::endl;
990  if( eta < -2.0 ) return 1;
991  if( eta < -1.8 ) return 2;
992  if( eta < -1.6 ) return 3;
993  if( eta < -1.2 ) return 4;
994  if( eta < -0.8 ) return 5;
995  if( eta < 0. ) return 6;
996  if( eta < 0.8 ) return 7;
997  if( eta < 1.2 ) return 8;
998  if( eta < 1.6 ) return 9;
999  if( eta < 1.8 ) return 10;
1000  if( eta < 2.0 ) return 11;
1001  return 12;
1002  }
1003 
1004  resolutionFunctionType47() { this->parNum_ = 14; }
1005 
1006  double sigmaPt(const double & pt, const double & eta, const T & parval) override
1007  {
1008  return sqrt(pow(parval[0]*pt,2) + pow(parval[etaBin(eta)],2) + pow(parval[13],2));
1009  }
1010  double sigmaCotgTh(const double & pt, const double & eta, const T & parval) override {
1011  return 0;
1012  }
1013  double sigmaPhi(const double & pt, const double & eta, const T & parval) override {
1014  return 0.;
1015  }
1016 
1017  double covPt1Pt2(const double & pt1, const double & eta1, const double & pt2, const double & eta2, const T & parval) override
1018  {
1019  return parval[14];
1020  }
1021 
1022  // derivatives ---------------
1023  double sigmaPtError(const double & pt, const double & eta, const T & parval, const T & parError) override
1024  {
1025  // Use the etaByBin function to select the right bin for the parameter
1026  double r = sqrt(pow(parval[0]*pt,2) + pow(parval[etaBin(eta)],2) + pow(parval[13],2));
1027  return sqrt( pow(pt*pt*parval[0]*parError[0],2) + pow(parval[etaBin(eta)]*parError[etaBin(eta)],2) + pow(parval[13]*parError[13],2) )/r;
1028  }
1029 
1030  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind,
1031  TString* parname, const T & parResol, const std::vector<int> & parResolOrder,
1032  const int muonType) override
1033  {
1034  std::vector<ParameterSet> parSet(this->parNum_);
1035  // name, step, mini, maxi
1036  parSet[0] = ParameterSet( "Pt res. sc.", 0.0002, 0., 0.1 );
1037  parSet[1] = ParameterSet( "eta bin 1", 0.00002, 0., 0.01 );
1038  parSet[2] = ParameterSet( "eta bin 2", 0.00002, 0., 0.01 );
1039  parSet[3] = ParameterSet( "eta bin 3", 0.00002, 0., 0.01 );
1040  parSet[4] = ParameterSet( "eta bin 4", 0.00002, 0., 0.01 );
1041  parSet[5] = ParameterSet( "eta bin 5", 0.00002, 0., 0.01 );
1042  parSet[6] = ParameterSet( "eta bin 6", 0.00002, 0., 0.01 );
1043  parSet[7] = ParameterSet( "eta bin 7", 0.00002, 0., 0.01 );
1044  parSet[8] = ParameterSet( "eta bin 8", 0.00002, 0., 0.01 );
1045  parSet[9] = ParameterSet( "eta bin 9", 0.00002, 0., 0.01 );
1046  parSet[10] = ParameterSet( "eta bin 10", 0.00002, 0., 0.01 );
1047  parSet[11] = ParameterSet( "eta bin 11", 0.00002, 0., 0.01 );
1048  parSet[12] = ParameterSet( "eta bin 12", 0.00002, 0., 0.01 );
1049  parSet[13] = ParameterSet( "cov(pt1,pt2)", 0.00002, 0., 0.01 );
1050 
1051  std::cout << "setting parameters" << std::endl;
1052  this->setPar( Start, Step, Mini, Maxi, ind, parname, parResol, parResolOrder, parSet );
1053  }
1054 
1055  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
1056  const T & parResol, const std::vector<int> & parResolOrder,
1057  const std::vector<double> & parStep,
1058  const std::vector<double> & parMin,
1059  const std::vector<double> & parMax,
1060  const int muonType) override
1061  {
1062  if( (int(parStep.size()) != this->parNum_) || (int(parMin.size()) != this->parNum_) || (int(parMax.size()) != this->parNum_) ) {
1063  std::cout << "Error: par step or min or max do not match with number of parameters" << std::endl;
1064  std::cout << "parNum = " << this->parNum_ << std::endl;
1065  std::cout << "parStep.size() = " << parStep.size() << std::endl;
1066  std::cout << "parMin.size() = " << parMin.size() << std::endl;
1067  std::cout << "parMax.size() = " << parMax.size() << std::endl;
1068  exit(1);
1069  }
1070  std::vector<ParameterSet> parSet(this->parNum_);
1071  // name, step, mini, maxi
1072  parSet[0] = ParameterSet( "Pt res. sc.", parStep[0], parMin[0], parMax[0] );
1073  parSet[1] = ParameterSet( "eta bin 1", parStep[1], parMin[1], parMax[1] );
1074  parSet[2] = ParameterSet( "eta bin 2", parStep[2], parMin[2], parMax[2] );
1075  parSet[3] = ParameterSet( "eta bin 3", parStep[3], parMin[3], parMax[3] );
1076  parSet[4] = ParameterSet( "eta bin 4", parStep[4], parMin[4], parMax[4] );
1077  parSet[5] = ParameterSet( "eta bin 5", parStep[5], parMin[5], parMax[5] );
1078  parSet[6] = ParameterSet( "eta bin 6", parStep[6], parMin[6], parMax[6] );
1079  parSet[7] = ParameterSet( "eta bin 7", parStep[7], parMin[7], parMax[7] );
1080  parSet[8] = ParameterSet( "eta bin 8", parStep[8], parMin[8], parMax[8] );
1081  parSet[9] = ParameterSet( "eta bin 9", parStep[9], parMin[9], parMax[9] );
1082  parSet[10] = ParameterSet( "eta bin 10", parStep[10], parMin[10], parMax[10] );
1083  parSet[11] = ParameterSet( "eta bin 11", parStep[11], parMin[11], parMax[11] );
1084  parSet[12] = ParameterSet( "eta bin 12", parStep[12], parMin[12], parMax[12] );
1085  parSet[13] = ParameterSet( "cov(pt1,pt2)",parStep[13], parMin[13], parMax[13] );
1086 
1087  std::cout << "setting parameters" << std::endl;
1088  for( int i=0; i<this->parNum_; ++i ) {
1089  std::cout << "parStep["<<i<<"] = " << parStep[i]
1090  << ", parMin["<<i<<"] = " << parMin[i]
1091  << ", parMax["<<i<<"] = " << parMin[i] << std::endl;
1092  }
1093  this->setPar( Start, Step, Mini, Maxi, ind, parname, parResol, parResolOrder, parSet );
1094  }
1095 };
1096 
1097 
1098 // ------------ ATTENTION ----------- //
1099 // Other functions are not in for now //
1100 // ---------------------------------- //
1101 
1103 resolutionFunctionBase<double *> * resolutionFunctionService( const int identifier );
1104 
1106 resolutionFunctionBase<std::vector<double> > * resolutionFunctionVecService( const int identifier );
1107 
1132 class backgroundFunctionBase {
1133  public:
1134  backgroundFunctionBase(const double & lowerLimit, const double & upperLimit) :
1135  lowerLimit_(lowerLimit), upperLimit_(upperLimit) {}
1136  virtual ~backgroundFunctionBase()
1137  {
1138  delete functionForIntegral_;
1139  };
1140  virtual double operator()( const double * parval, const double & mass, const double & eta ) const = 0;
1141  virtual double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const
1142  {
1143  return operator()(parval, mass, eta1);
1144  }
1145  virtual int parNum() const { return parNum_; }
1147  virtual void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) = 0;
1148  virtual TF1* functionForIntegral(const std::vector<double>::const_iterator & parBgrIt) const
1149  {
1150  functionForIntegral_ = new FunctionForIntegral(this, parBgrIt);
1151  TF1 * backgroundFunctionForIntegral = new TF1("backgroundFunctionForIntegral", functionForIntegral_,
1152  lowerLimit_, upperLimit_, this->parNum_);
1153  return( backgroundFunctionForIntegral );
1154  }
1155  virtual double fracVsEta(const double * parval, const double & eta1, const double & eta2) const { return 1.; }
1156 
1157 protected:
1158  int parNum_;
1159  double lowerLimit_;
1160  double upperLimit_;
1162  virtual void setPar(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname,
1163  const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt,
1164  double* thisStep, double* thisMini, double* thisMaxi, TString* thisParName ) {
1165  for( int iPar=0; iPar<this->parNum_; ++iPar ) {
1166  Start[iPar] = *(parBgrIt+iPar);
1167  Step[iPar] = thisStep[iPar];
1168  Mini[iPar] = thisMini[iPar];
1169  Maxi[iPar] = thisMaxi[iPar];
1170  ind[iPar] = *(parBgrOrderIt+iPar);
1171  // EM 2012.05.22 this line is crashing cmsRun (need to be fixed) parname[iPar] = thisParName[iPar];
1172  }
1173  }
1174  class FunctionForIntegral
1175  {
1176  public:
1177  FunctionForIntegral( const backgroundFunctionBase * function,
1178  const std::vector<double>::const_iterator & parBgrIt ) :
1179  function_(function)
1180  {
1181  parval_ = new double[function_->parNum()];
1182  for( int i=0; i < function_->parNum(); ++i ) {
1183  parval_[i] = *(parBgrIt+i);
1184  }
1185  }
1186  ~FunctionForIntegral()
1187  {
1188  delete parval_;
1189  }
1190  double operator()(const double * mass, const double *) const
1191  {
1192  // FIXME: this is a gross approximation. The function should be integrated in eta over the sample.
1193  return( (*function_)(parval_, *mass, 0.) );
1194  }
1195  protected:
1196  const backgroundFunctionBase * function_;
1197  double * parval_;
1198  };
1199  mutable FunctionForIntegral * functionForIntegral_;
1200 };
1201 
1203 // -------
1204 class backgroundFunctionType1 : public backgroundFunctionBase
1205 {
1206  public:
1213  backgroundFunctionType1(const double & lowerLimit, const double & upperLimit) :
1214  backgroundFunctionBase(lowerLimit, upperLimit)
1215  { this->parNum_ = 2; }
1216  double operator()( const double * parval, const double & mass, const double & eta ) const override
1217  {
1218  double a = 1.;
1219  double b = parval[1];
1220 
1221  double norm = -(a*lowerLimit_ + b*lowerLimit_*lowerLimit_/2.);
1222 
1223  if( -a/b > upperLimit_ ) norm += a*upperLimit_ + b*upperLimit_*upperLimit_/2.;
1224  else norm += -a*a/(2*b);
1225 
1226  if( mass < -a/b && norm != 0 ) return (a + b*mass)/norm;
1227  else return 0;
1228  }
1229  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1230  double thisStep[] = {0.01, 0.01};
1231  TString thisParName[] = {"Constant", "Linear"};
1232  if( muonType == 1 ) {
1233  double thisMini[] = {0.0, -300.};
1234  double thisMaxi[] = {1.0, 0.};
1235  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1236  } else {
1237  double thisMini[] = {0.0, -300.};
1238  double thisMaxi[] = {1.0, 0.};
1239  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1240  }
1241  }
1242 };
1244 // ------------
1245 class backgroundFunctionType2 : public backgroundFunctionBase {
1246  public:
1252  backgroundFunctionType2(const double & lowerLimit, const double & upperLimit) :
1253  backgroundFunctionBase(lowerLimit, upperLimit)
1254  { this->parNum_ = 2; }
1255  double operator()( const double * parval, const double & mass, const double & eta ) const override
1256  {
1257  double Bgrp2 = parval[1];
1258  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
1259  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
1260  else return 0.;
1261  }
1262  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1263  double thisStep[] = {0.01, 0.01};
1264  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
1265  if( muonType == 1 ) {
1266  double thisMini[] = {0.0, 0.};
1267  double thisMaxi[] = {1.0, 10.};
1268  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1269  } else {
1270  double thisMini[] = {0.0, 0.};
1271  double thisMaxi[] = {1.0, 10.};
1272  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1273  }
1274  }
1275 
1276 
1277 
1278 
1279  // virtual double fracVsEta(const double * parval, const double & resEta) const
1280  // {
1281  // // return( 0.6120 - 0.0225*eta*eta );
1282  // return( 1. - 0.0225*resEta*resEta ); // so that a = 1 for eta = 0.
1283  // }
1284 
1285 
1286 
1287 
1288 };
1289 
1290 
1291 
1292 
1293 
1294 
1295 
1300 //
1301 //class backgroundFunctionType3 : public backgroundFunctionBase {
1302 // public:
1303 // // pass parval[shift]
1304 // backgroundFunctionType3(const double & lowerLimit, const double & upperLimit) :
1305 // backgroundFunctionBase(lowerLimit, upperLimit)
1306 // { this->parNum_ = 3; }
1307 // virtual double operator()( const double * parval, const int resTotNum, const int ires, const bool * resConsidered,
1308 // const double * ResMass, const double ResHalfWidth[], const int MuonType, const double & mass, const int nbins ) {
1309 // double PB = 0.;
1310 // double Bgrp2 = parval[1];
1311 // double Bgrp3 = parval[2];
1312 // for (int ires=0; ires<resTotNum; ires++) {
1313 // // In this case, by integrating between A and B, we get for f=exp(a-bx)+k:
1314 // // INT = exp(a)/b*(exp(-bA)-exp(-bB))+k*(B-A) so our function, which in 1000 bins between A and B
1315 // // gets a total of 1, is f = (exp(a-bx)+k)*(B-A)/nbins / (INT)
1316 // // ----------------------------------------------------------------------------------------------
1317 // if (resConsidered[ires]) {
1318 // if (exp(-Bgrp2*(ResMass[ires]-ResHalfWidth[ires]))-exp(-Bgrp2*(ResMass[ires]+ResHalfWidth[ires]))>0) {
1319 // PB += (exp(-Bgrp2*mass)+Bgrp3) *
1320 // 2*ResHalfWidth[ires]/(double)nbins /
1321 // ( (exp(-Bgrp2*(ResMass[ires]-ResHalfWidth[ires]))-exp(-Bgrp2*(ResMass[ires]+ResHalfWidth[ires])))/
1322 // Bgrp2 + Bgrp3*2*ResHalfWidth[ires] );
1323 // } else {
1324 // std::cout << "Impossible to compute Background probability! - some fix needed - Bgrp2=" << Bgrp2 << std::endl;
1325 // }
1326 // }
1327 // }
1328 // return PB;
1329 // }
1330 // virtual void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) {
1331 // double thisStep[] = {0.1, 0.001, 0.1};
1332 // TString thisParName[] = {"Bgr fraction", "Bgr slope", "Bgr constant"};
1333 // if( muonType == 1 ) {
1334 // double thisMini[] = {0.0, 0.000000001, 0.0};
1335 // double thisMaxi[] = {1.0, 0.2, 1000};
1336 // this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1337 // } else {
1338 // double thisMini[] = {0.0, 0.000000001, 0.0};
1339 // double thisMaxi[] = {1.0, 0.2, 1000};
1340 // this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1341 // }
1342 // }
1343 // virtual TF1* functionForIntegral(const std::vector<double>::const_iterator & parBgrIt) const {return 0;};
1344 //};
1345 
1346 
1347 
1349 // --------------------------------
1350 class backgroundFunctionType4 : public backgroundFunctionBase
1351 {
1352  public:
1358  backgroundFunctionType4(const double & lowerLimit, const double & upperLimit) :
1359  backgroundFunctionBase(lowerLimit, upperLimit)
1360  { this->parNum_ = 4; }
1361  double operator()( const double * parval, const double & mass, const double & eta ) const override
1362  {
1363  double Bgrp2 = parval[1] + parval[2]*eta*eta;
1364  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
1365  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
1366  else return 0.;
1367  }
1368  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1369  double thisStep[] = {0.01, 0.01, 0.01, 0.01};
1370  TString thisParName[] = {"Bgr fraction", "Bgr slope", "Bgr slope eta^2 dependence", "background fraction eta dependence"};
1371  if( muonType == 1 ) {
1372  double thisMini[] = {0.0, 0., 0., -1.};
1373  double thisMaxi[] = {1.0, 10., 10., 1.};
1374  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1375  } else {
1376  double thisMini[] = {0.0, 0., -1., -1.};
1377  double thisMaxi[] = {1.0, 10., 1., 1.};
1378  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1379  }
1380  }
1381  /* virtual double fracVsEta(const double * parval, const double & resEta) const */
1382  /* { */
1383  /* return( 1. - parval[3]*resEta*resEta ); // so that a = 1 for eta = 0. */
1384  /* } */
1385 };
1386 
1388 // ---------------------------
1389 class backgroundFunctionType5 : public backgroundFunctionBase
1390 {
1391  public:
1396  backgroundFunctionType5(const double & lowerLimit, const double & upperLimit) :
1397  backgroundFunctionBase(lowerLimit, upperLimit)
1398  { this->parNum_ = 3; }
1399  double operator()( const double * parval, const double & mass, const double & eta ) const override
1400  {
1401  double b = parval[1];
1402  // double c = parval[2];
1403  double a = 1 + parval[2]*eta*eta;
1404 
1405  double norm = -(a*lowerLimit_ + b*lowerLimit_*lowerLimit_/2.);
1406 
1407  if( -a/b > upperLimit_ ) norm += a*upperLimit_ + b*upperLimit_*upperLimit_/2.;
1408  else norm += -a*a/(2*b);
1409 
1410  if( mass < -a/b && norm != 0 ) return (a + b*mass)/norm;
1411  else return 0;
1412  }
1413  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1414  double thisStep[] = {0.01, 0.01, 0.01};
1415  TString thisParName[] = {"Bgr fraction", "Constant", "Linear"};
1416  if( muonType == 1 ) {
1417  double thisMini[] = {0.0, 0., -300.};
1418  double thisMaxi[] = {1.0, 300., 0.};
1419  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1420  } else {
1421  double thisMini[] = {0.0, 0., -300.};
1422  double thisMaxi[] = {1.0, 300., 0.};
1423  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1424  }
1425  }
1426 };
1427 
1428 
1430 // --------------------------
1431 class backgroundFunctionType6 : public backgroundFunctionBase {
1432  public:
1433  backgroundFunctionType6(const double & lowerLimit, const double & upperLimit) :
1434  backgroundFunctionBase(lowerLimit, upperLimit)
1435  {
1436  this->parNum_ = 2;
1437  }
1438  double operator()( const double * parval, const double & mass, const double & eta ) const override {return 0.;}
1439  double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const override
1440  {
1441  double Bgrp2 = 0.;
1442  if( fabs(eta1) <= 1.3 && fabs(eta2) <= 1.3 ) {
1443  Bgrp2 = -1.20528;
1444  }
1445  else if( (fabs(eta1) <= 1.6 && fabs(eta1) > 1.3) && (fabs(eta2) <= 1.6 && fabs(eta2) > 1.3) ) {
1446  Bgrp2 = 0.234713;
1447  }
1448  else if( fabs(eta1) > 1.6 && fabs(eta2) > 1.6 ) {
1449  Bgrp2 = -0.667103;
1450  }
1451  else if( (fabs(eta1) <= 1.3 && (fabs(eta2) > 1.3 && fabs(eta2) <= 1.6)) ||
1452  (fabs(eta2) <= 1.3 && (fabs(eta1) > 1.3 && fabs(eta1) <= 1.6)) ) {
1453  Bgrp2 = -0.656904;
1454  }
1455  else if( (fabs(eta1) <= 1.3 && fabs(eta2) > 1.6) ||
1456  (fabs(eta2) <= 1.3 && fabs(eta1) > 1.6) ) {
1457  Bgrp2 = 0.155328;
1458  }
1459  else if( ((fabs(eta1) > 1.3 && fabs(eta1) <= 1.6) && fabs(eta2) > 1.6) ||
1460  ((fabs(eta2) > 1.3 && fabs(eta2) <= 1.6) && fabs(eta1) > 1.6) ) {
1461  Bgrp2 = -0.177154;
1462  }
1463  else {
1464  std::cout << "WARNING: this should not happen for eta1 = " << eta1 << " and eta2 = " << eta2 << std::endl;
1465  Bgrp2 = -0.667103;
1466  }
1467  Bgrp2*=-1.;
1468  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
1469  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
1470  else return 0.;
1471 
1472  }
1473  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1474  double thisStep[] = {0.01, 0.01};
1475  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
1476  if( muonType == 1 ) {
1477  double thisMini[] = {0.0, 0.};
1478  double thisMaxi[] = {1.0, 10.};
1479  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1480  } else {
1481  double thisMini[] = {0.0, 0.};
1482  double thisMaxi[] = {1.0, 10.};
1483  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1484  }
1485  }
1486 
1487  double fracVsEta(const double * parval, const double & eta1, const double & eta2) const override
1488  {
1489  if( fabs(eta1) <= 1.3 && fabs(eta2) <= 1.3 ) {
1490  return (1.-0.910903);
1491  }
1492  else if( (fabs(eta1) <= 1.6 && fabs(eta1) > 1.3) && (fabs(eta2) <= 1.6 && fabs(eta2) > 1.3) ) {
1493  return (1.-0.801469);
1494  }
1495  else if( fabs(eta1) > 1.6 && fabs(eta2) > 1.6 ) {
1496  return (1.-0.658196);
1497  }
1498  else if( (fabs(eta1) <= 1.3 && (fabs(eta2) > 1.3 && fabs(eta2) <= 1.6)) ||
1499  (fabs(eta2) <= 1.3 && (fabs(eta1) > 1.3 && fabs(eta1) <= 1.6)) ) {
1500  return (1.-0.873411);
1501  }
1502  else if( (fabs(eta1) <= 1.3 && fabs(eta2) > 1.6) ||
1503  (fabs(eta2) <= 1.3 && fabs(eta1) > 1.6) ) {
1504  return (1.-0.784674);
1505  }
1506  else if( ((fabs(eta1) > 1.3 && fabs(eta1) <= 1.6) && fabs(eta2) > 1.6) ||
1507  ((fabs(eta2) > 1.3 && fabs(eta2) <= 1.6) && fabs(eta1) > 1.6) ) {
1508  return (1.-0.714398);
1509  }
1510  else {
1511  std::cout << "WARNING: this should not happen for eta1 = " << eta1 << " and eta2 = " << eta2 << std::endl;
1512  return (1.-0.658196);
1513  }
1514  }
1515 };
1516 
1518 // ---------------------------------------------------------
1519 class backgroundFunctionType7 : public backgroundFunctionBase {
1520  public:
1521  backgroundFunctionType7(const double & lowerLimit, const double & upperLimit) :
1522  backgroundFunctionBase(lowerLimit, upperLimit)
1523  {
1524  this->parNum_ = 2;
1525  }
1526  double operator()( const double * parval, const double & mass, const double & eta ) const override {return 0.;}
1527  double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const override
1528  {
1529  double Bgrp2 = 0.;
1530  if( (fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 0. && fabs(eta2) < 0.9) ) {
1531  Bgrp2 = (-1.42465);
1532  }
1533  else if( (fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) ) {
1534  Bgrp2 = (-1.38576);
1535  }
1536  else if( (fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) ) {
1537  Bgrp2 = (-0.333728);
1538  }
1539  else if( (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) ) {
1540  Bgrp2 = (0.94066);
1541  }
1542  else if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) ) {
1543  Bgrp2 = (0.371026);
1544  }
1545  else if( (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) ) {
1546  Bgrp2 = (-0.959101);
1547  }
1548  else if( (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9) ) {
1549  Bgrp2 = (-1.13829);
1550  }
1551  else if( (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0) ) {
1552  Bgrp2 = (-0.921581);
1553  }
1554  else if( (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.) ) {
1555  Bgrp2 = (-0.664338);
1556  }
1557  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 0.9 && fabs(eta2) < 1.3)) ||
1558  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 0.9 && fabs(eta1) < 1.3)) ) {
1559  Bgrp2 = (-1.07581);
1560  }
1561  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.3 && fabs(eta2) < 1.5)) ||
1562  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.3 && fabs(eta1) < 1.5)) ) {
1563  Bgrp2 = (-0.250272);
1564  }
1565  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6)) ||
1566  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6)) ) {
1567  Bgrp2 = (0.101785);
1568  }
1569  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1570  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1571  Bgrp2 = (0.360397);
1572  }
1573  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1574  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1575  Bgrp2 = (0.689136);
1576  }
1577  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1578  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1579  Bgrp2 = (0.860723);
1580  }
1581  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1582  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1583  Bgrp2 = (1.21908);
1584  }
1585  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1586  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1587  Bgrp2 = (2.4453);
1588  }
1589  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.3 && fabs(eta2) < 1.5)) ||
1590  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.3 && fabs(eta1) < 1.5)) ) {
1591  Bgrp2 = (-1.14152);
1592  }
1593  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6)) ||
1594  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6)) ) {
1595  Bgrp2 = (-0.77241);
1596  }
1597  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1598  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1599  Bgrp2 = (-0.516479);
1600  }
1601  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1602  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1603  Bgrp2 = (-0.361401);
1604  }
1605  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1606  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1607  Bgrp2 = (-0.33143);
1608  }
1609  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1610  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1611  Bgrp2 = (-0.20813);
1612  }
1613  else if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1614  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1615  Bgrp2 = (0.158002);
1616  }
1617  else if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6)) ||
1618  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6)) ) {
1619  Bgrp2 = (0.273222);
1620  }
1621  else if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1622  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1623  Bgrp2 = (0.247639);
1624  }
1625  else if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1626  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1627  Bgrp2 = (-0.148616);
1628  }
1629  else if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1630  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1631  Bgrp2 = (-0.413175);
1632  }
1633  else if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1634  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1635  Bgrp2 = (-0.230031);
1636  }
1637  else if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1638  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1639  Bgrp2 = (-0.122756);
1640  }
1641  else if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1642  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1643  Bgrp2 = (0.650851);
1644  }
1645  else if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1646  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1647  Bgrp2 = (-0.0985001);
1648  }
1649  else if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1650  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1651  Bgrp2 = (-0.402548);
1652  }
1653  else if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1654  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1655  Bgrp2 = (-0.27401);
1656  }
1657  else if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1658  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1659  Bgrp2 = (-0.22863);
1660  }
1661  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1662  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1663  Bgrp2 = (-0.436959);
1664  }
1665  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1666  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1667  Bgrp2 = (-0.506041);
1668  }
1669  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1670  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1671  Bgrp2 = (-0.31618);
1672  }
1673  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1674  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1675  Bgrp2 = (-0.365653);
1676  }
1677  else if( ((fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1678  ((fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1679  Bgrp2 = (-1.16783);
1680  }
1681  else if( ((fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1682  ((fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1683  Bgrp2 = (-0.730701);
1684  }
1685  else if( ((fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1686  ((fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1687  Bgrp2 = (-0.5271);
1688  }
1689  else if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 1.9) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1690  ((fabs(eta2) >= 1.8 && fabs(eta2) < 1.9) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1691  Bgrp2 = (-0.99893);
1692  }
1693  else if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 1.9) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1694  ((fabs(eta2) >= 1.8 && fabs(eta2) < 1.9) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1695  Bgrp2 = (-0.687263);
1696  }
1697  else if( ((fabs(eta1) >= 1.9 && fabs(eta1) < 2.0) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1698  ((fabs(eta2) >= 1.9 && fabs(eta2) < 2.0) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1699  Bgrp2 = (-0.722394);
1700  }
1701  else {
1702  std::cout << "WARNING: this should not happen for eta1 = " << eta1 << " and eta2 = " << eta2 << std::endl;
1703  Bgrp2 = -0.664338;
1704  }
1705  Bgrp2*=-1.;
1706  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
1707  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
1708  else return 0.;
1709 
1710  }
1711  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1712  double thisStep[] = {0.01, 0.01};
1713  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
1714  if( muonType == 1 ) {
1715  double thisMini[] = {0.0, 0.};
1716  double thisMaxi[] = {1.0, 10.};
1717  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1718  } else {
1719  double thisMini[] = {0.0, 0.};
1720  double thisMaxi[] = {1.0, 10.};
1721  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1722  }
1723  }
1724 
1725  double fracVsEta(const double * parval, const double & eta1, const double & eta2) const override
1726  {
1727  if( (fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 0. && fabs(eta2) < 0.9) ) {
1728  return (1.-0.915365);
1729  }
1730  if( (fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) ) {
1731  return (1.-0.914149);
1732  }
1733  if( (fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) ) {
1734  return (1.-0.855918);
1735  }
1736  if( (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) ) {
1737  return (1.-0.70221);
1738  }
1739  if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) ) {
1740  return (1.-0.701489);
1741  }
1742  if( (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) ) {
1743  return (1.-0.651162);
1744  }
1745  if( (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9) ) {
1746  return (1.-0.639839);
1747  }
1748  if( (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0) ) {
1749  return (1.-0.64915);
1750  }
1751  if( (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.) ) {
1752  return (1.-0.687878);
1753  }
1754  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 0.9 && fabs(eta2) < 1.3)) ||
1755  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 0.9 && fabs(eta1) < 1.3)) ) {
1756  return (1.-0.903486);
1757  }
1758  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.3 && fabs(eta2) < 1.5)) ||
1759  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.3 && fabs(eta1) < 1.5)) ) {
1760  return (1.-0.882516);
1761  }
1762  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6)) ||
1763  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6)) ) {
1764  return (1.-0.85477);
1765  }
1766  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1767  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1768  return (1.-0.804919);
1769  }
1770  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1771  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1772  return (1.-0.75411);
1773  }
1774  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1775  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1776  return (1.-0.714128);
1777  }
1778  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1779  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1780  return (1.-0.645403);
1781  }
1782  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.9) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1783  ((fabs(eta2) >= 0. && fabs(eta2) < 0.9) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1784  return (1.-0.588049);
1785  }
1786  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.3 && fabs(eta2) < 1.5)) ||
1787  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.3 && fabs(eta1) < 1.5)) ) {
1788  return (1.-0.901123);
1789  }
1790  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6)) ||
1791  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6)) ) {
1792  return (1.-0.87852);
1793  }
1794  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1795  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1796  return (1.-0.862266);
1797  }
1798  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1799  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1800  return (1.-0.846385);
1801  }
1802  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1803  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1804  return (1.-0.825401);
1805  }
1806  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1807  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1808  return (1.-0.812449);
1809  }
1810  if( ((fabs(eta1) >= 0.9 && fabs(eta1) < 1.3) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1811  ((fabs(eta2) >= 0.9 && fabs(eta2) < 1.3) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1812  return (1.-0.753754);
1813  }
1814  if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.5 && fabs(eta2) < 1.6)) ||
1815  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.5 && fabs(eta1) < 1.6)) ) {
1816  return (1.-0.794143);
1817  }
1818  if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1819  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1820  return (1.-0.761375);
1821  }
1822  if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1823  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1824  return (1.-0.765572);
1825  }
1826  if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1827  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1828  return (1.-0.749438);
1829  }
1830  if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1831  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1832  return (1.-0.750941);
1833  }
1834  if( ((fabs(eta1) >= 1.3 && fabs(eta1) < 1.5) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1835  ((fabs(eta2) >= 1.3 && fabs(eta2) < 1.5) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1836  return (1.-0.722832);
1837  }
1838  if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.7)) ||
1839  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.7)) ) {
1840  return (1.-0.699723);
1841  }
1842  if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1843  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1844  return (1.-0.734044);
1845  }
1846  if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1847  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1848  return (1.-0.719434);
1849  }
1850  if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1851  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1852  return (1.-0.718889);
1853  }
1854  if( ((fabs(eta1) >= 1.5 && fabs(eta1) < 1.6) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1855  ((fabs(eta2) >= 1.5 && fabs(eta2) < 1.6) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1856  return (1.-0.689382);
1857  }
1858  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.7 && fabs(eta2) < 1.8)) ||
1859  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 1.7 && fabs(eta1) < 1.8)) ) {
1860  return (1.-0.681106);
1861  }
1862  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1863  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1864  return (1.-0.685783);
1865  }
1866  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1867  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1868  return (1.-0.695924);
1869  }
1870  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.7) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1871  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.7) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1872  return (1.-0.670977);
1873  }
1874  if( ((fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.8 && fabs(eta2) < 1.9)) ||
1875  ((fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) && (fabs(eta1) >= 1.8 && fabs(eta1) < 1.9)) ) {
1876  return (1.-0.654816);
1877  }
1878  if( ((fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1879  ((fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1880  return (1.-0.670969);
1881  }
1882  if( ((fabs(eta1) >= 1.7 && fabs(eta1) < 1.8) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1883  ((fabs(eta2) >= 1.7 && fabs(eta2) < 1.8) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1884  return (1.-0.659082);
1885  }
1886  if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 1.9) && (fabs(eta2) >= 1.9 && fabs(eta2) < 2.0)) ||
1887  ((fabs(eta2) >= 1.8 && fabs(eta2) < 1.9) && (fabs(eta1) >= 1.9 && fabs(eta1) < 2.0)) ) {
1888  return (1.-0.648371);
1889  }
1890  if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 1.9) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1891  ((fabs(eta2) >= 1.8 && fabs(eta2) < 1.9) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1892  return (1.-0.659114);
1893  }
1894  if( ((fabs(eta1) >= 1.9 && fabs(eta1) < 2.0) && (fabs(eta2) >= 2.0 && fabs(eta2) < 1000.)) ||
1895  ((fabs(eta2) >= 1.9 && fabs(eta2) < 2.0) && (fabs(eta1) >= 2.0 && fabs(eta1) < 1000.)) ) {
1896  return (1.-0.660482);
1897  }
1898  else {
1899  std::cout << "WARNING: this should not happen for eta1 = " << eta1 << " and eta2 = " << eta2 << std::endl;
1900  return (1.-0.687878);
1901  }
1902  }
1903 };
1904 //Function 8
1905 class backgroundFunctionType8 : public backgroundFunctionBase {
1906  public:
1907  backgroundFunctionType8(const double & lowerLimit, const double & upperLimit) :
1908  backgroundFunctionBase(lowerLimit, upperLimit)
1909  {
1910  this->parNum_ = 2;
1911  }
1912  double operator()( const double * parval, const double & mass, const double & eta ) const override {return 0.;}
1913  double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const override
1914  {
1915  double Bgrp2 = 0.;
1916 if( (fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0. && fabs(eta2) < 0.85) ) {
1917  Bgrp2 = (-2.17047);
1918 }
1919 else if( (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) ) {
1920  Bgrp2 = (-2.12913);
1921 }
1922 else if( (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) ) {
1923  Bgrp2 = (-2.19963);
1924 }
1925 else if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1000) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000) ) {
1926  Bgrp2 = (-0.386394);
1927 }
1928 else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25)) ||
1929  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25)) ) {
1930  Bgrp2 = (-1.71339);
1931 }
1932 else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
1933  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
1934  Bgrp2 = (-0.206566);
1935 }
1936 else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
1937  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
1938  Bgrp2 = (4.4815);
1939 }
1940 else if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
1941  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
1942  Bgrp2 = (-1.87985);
1943 }
1944 else if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
1945  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
1946  Bgrp2 = (-0.163569);
1947 }
1948 else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
1949  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
1950  Bgrp2 = (-1.67935);
1951 }
1952  Bgrp2*=-1.;
1953  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
1954  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
1955  else return 0.;
1956 
1957  }
1958  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
1959  double thisStep[] = {0.01, 0.01};
1960  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
1961  if( muonType == 1 ) {
1962  double thisMini[] = {0.0, 0.};
1963  double thisMaxi[] = {1.0, 10.};
1964  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1965  } else {
1966  double thisMini[] = {0.0, 0.};
1967  double thisMaxi[] = {1.0, 10.};
1968  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
1969  }
1970  }
1971 
1972  double fracVsEta(const double * parval, const double & eta1, const double & eta2) const override
1973  {
1974 if( (fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0. && fabs(eta2) < 0.85) ) {
1975  return (1.-0.907727);
1976 }
1977 if( (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) ) {
1978  return (1.-0.907715);
1979 }
1980 if( (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) ) {
1981  return (1.-0.912233);
1982 }
1983 if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1000) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000) ) {
1984  return (1.-0.876776);
1985 }
1986 if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25)) ||
1987  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25)) ) {
1988  return (1.-0.913046);
1989 }
1990 if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
1991  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
1992  return (1.-0.916765);
1993 }
1994 if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
1995  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
1996  return (1.-0.6);
1997 }
1998 if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
1999  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2000  return (1.-0.907471);
2001 }
2002 if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2003  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2004  return (1.-0.899253);
2005 }
2006 if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2007  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2008  return (1.-0.879108);
2009 }
2010  return 0;}
2011 };
2012 
2014 backgroundFunctionBase * backgroundFunctionService( const int identifier, const double & lowerLimit, const double & upperLimit );
2015 
2016 //Function Type 9
2017 
2018 class backgroundFunctionType9 : public backgroundFunctionBase {
2019  public:
2020  backgroundFunctionType9(const double & lowerLimit, const double & upperLimit) :
2021  backgroundFunctionBase(lowerLimit, upperLimit)
2022  {
2023  this->parNum_ = 2;
2024  }
2025  double operator()( const double * parval, const double & mass, const double & eta ) const override {return 0.;}
2026  double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const override
2027  {
2028  double Bgrp2 = 0.;
2029 
2030 if( (fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0. && fabs(eta2) < 0.85) ) {
2031  Bgrp2 = (-1.80833);
2032 }
2033 else if( (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) ) {
2034  Bgrp2 = (-1.98281);
2035 }
2036 else if( (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) ) {
2037  Bgrp2 = (-1.79632);
2038 }
2039 else if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1000) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000) ) {
2040  Bgrp2 = (-1.14645);
2041 }
2042 else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25)) ||
2043  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25)) ) {
2044  Bgrp2 = (-1.55747);
2045 }
2046 else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2047  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2048  Bgrp2 = (-0.337598);
2049 }
2050 else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2051  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2052  Bgrp2 = (5.36513);
2053 }
2054 else if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2055  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2056  Bgrp2 = (-1.44363);
2057 }
2058 else if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2059  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2060  Bgrp2 = (-0.54614);
2061 }
2062 else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2063  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2064  Bgrp2 = (-1.41059);
2065 }
2066 else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2067  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2068  Bgrp2 = (-1.41059);
2069 }
2070  Bgrp2*=-1.;
2071  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
2072  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
2073  else return 0.;
2074 
2075  }
2076  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
2077  double thisStep[] = {0.01, 0.01};
2078  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
2079  if( muonType == 1 ) {
2080  double thisMini[] = {0.0, 0.};
2081  double thisMaxi[] = {1.0, 10.};
2082  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
2083  } else {
2084  double thisMini[] = {0.0, 0.};
2085  double thisMaxi[] = {1.0, 10.};
2086  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
2087  }
2088  }
2089 
2090  double fracVsEta(const double * parval, const double & eta1, const double & eta2) const override
2091 
2092  {
2093  // std::cout << "Eta1 " << eta1 << " eta2 = " << eta2 << std::endl;
2094 if( (fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0. && fabs(eta2) < 0.85) ) {
2095  return (1.-0.893683);
2096 }
2097 if( (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) ) {
2098  return (1.-0.888968);
2099 }
2100 if( (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) ) {
2101  return (1.-0.885926);
2102 }
2103 if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1000) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000) ) {
2104  return (1.-0.866615);
2105 }
2106 if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25)) ||
2107  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25)) ) {
2108  return (1.-0.892856);
2109 }
2110 if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2111  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2112  return (1.-0.884864);
2113 }
2114 if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2115  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2116  return (1.-0.6);
2117 }
2118 if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2119  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2120  return (1.-0.894739);
2121 }
2122 if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2123  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2124  return (1.-0.880597);
2125 }
2126 if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2127  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2128  return (1.-0.869165);
2129 }
2130 if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2131  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2132  return (1.-0.869165);
2133 }
2134  return 0;}
2135 };
2136 
2137 
2138 class backgroundFunctionType10 : public backgroundFunctionBase {
2139  public:
2140  backgroundFunctionType10(const double & lowerLimit, const double & upperLimit) :
2141  backgroundFunctionBase(lowerLimit, upperLimit)
2142  {
2143  this->parNum_ = 2;
2144  }
2145  double operator()( const double * parval, const double & mass, const double & eta ) const override {return 0.;}
2146  double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const override
2147  {
2148  double Bgrp2 = 0.;
2149  if( (fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0. && fabs(eta2) < 0.85) ) {
2150  Bgrp2 = (-1.80833);
2151  }
2152  else if( (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) ) {
2153  Bgrp2 = (-1.98281);
2154  }
2155  else if( (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) ) {
2156  Bgrp2 = (-1.79632);
2157  }
2158  else if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) ) {
2159  Bgrp2 = (-1.50855);
2160  }
2161  else if( (fabs(eta1) >= 1.8 && fabs(eta1) < 2.0) && (fabs(eta2) >= 1.8 && fabs(eta2) < 2.0) ) {
2162  Bgrp2 = (-0.498511);
2163  }
2164  else if( (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2) ) {
2165  Bgrp2 = (-0.897031);
2166  }
2167  else if( (fabs(eta1) >= 2.2 && fabs(eta1) < 1000) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000) ) {
2168  Bgrp2 = (-0.75954);
2169  }
2170  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25)) ||
2171  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25)) ) {
2172  Bgrp2 = (-1.55747);
2173  }
2174  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2175  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2176  Bgrp2 = (-0.337598);
2177  }
2178  else if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2179  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2180  Bgrp2 = (3.5163);
2181  }
2182  else if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2183  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2184  Bgrp2 = (-1.44363);
2185  }
2186  else if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2187  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2188  Bgrp2 = (-0.54614);
2189  }
2190  else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.8)) ||
2191  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.8)) ) {
2192  Bgrp2 = (-1.36442);
2193  }
2194  else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.8 && fabs(eta2) < 2.0)) ||
2195  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.8 && fabs(eta1) < 2.0)) ) {
2196  Bgrp2 = (-1.66202);
2197  }
2198  else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2)) ||
2199  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2)) ) {
2200  Bgrp2 = (-0.62038);
2201  }
2202  else if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2203  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2204  Bgrp2 = (0.662449);
2205  }
2206  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.8 && fabs(eta2) < 2.0)) ||
2207  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) && (fabs(eta1) >= 1.8 && fabs(eta1) < 2.0)) ) {
2208  Bgrp2 = (-0.723325);
2209  }
2210  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2)) ||
2211  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) && (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2)) ) {
2212  Bgrp2 = (-1.54405);
2213  }
2214  else if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2215  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2216  Bgrp2 = (-1.1104);
2217  }
2218  else if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 2.0) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2)) ||
2219  ((fabs(eta2) >= 1.8 && fabs(eta2) < 2.0) && (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2)) ) {
2220  Bgrp2 = (-1.56277);
2221  }
2222  else if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 2.0) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2223  ((fabs(eta2) >= 1.8 && fabs(eta2) < 2.0) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2224  Bgrp2 = (-1.0827);
2225  }
2226  else if( ((fabs(eta1) >= 2.0 && fabs(eta1) < 2.2) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2227  ((fabs(eta2) >= 2.0 && fabs(eta2) < 2.2) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2228  Bgrp2 = (-1.05662);
2229  }
2230  Bgrp2*=-1.;
2231  double norm = -(exp(-Bgrp2*upperLimit_) - exp(-Bgrp2*lowerLimit_))/Bgrp2;
2232  if( norm != 0 ) return exp(-Bgrp2*mass)/norm;
2233  else return 0.;
2234  }
2235  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
2236  double thisStep[] = {0.01, 0.01};
2237  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
2238  if( muonType == 1 ) {
2239  double thisMini[] = {0.0, 0.};
2240  double thisMaxi[] = {1.0, 10.};
2241  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
2242  } else {
2243  double thisMini[] = {0.0, 0.};
2244  double thisMaxi[] = {1.0, 10.};
2245  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
2246  }
2247  }
2248 
2249  double fracVsEta(const double * parval, const double & eta1, const double & eta2) const override
2250  {
2251  if( (fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0. && fabs(eta2) < 0.85) ) {
2252  return (1.-0.893683);
2253  }
2254  if( (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) ) {
2255  return (1.-0.888968);
2256  }
2257  if( (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) ) {
2258  return (1.-0.885926);
2259  }
2260  if( (fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) ) {
2261  return (1.-0.892823);
2262  }
2263  if( (fabs(eta1) >= 1.8 && fabs(eta1) < 2.0) && (fabs(eta2) >= 1.8 && fabs(eta2) < 2.0) ) {
2264  return (1.-0.888735);
2265  }
2266  if( (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2) ) {
2267  return (1.-0.87497);
2268  }
2269  if( (fabs(eta1) >= 2.2 && fabs(eta1) < 1000) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000) ) {
2270  return (1.-0.895275);
2271  }
2272  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 0.85 && fabs(eta2) < 1.25)) ||
2273  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 0.85 && fabs(eta1) < 1.25)) ) {
2274  return (1.-0.892856);
2275  }
2276  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2277  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2278  return (1.-0.884864);
2279  }
2280  if( ((fabs(eta1) >= 0. && fabs(eta1) < 0.85) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2281  ((fabs(eta2) >= 0. && fabs(eta2) < 0.85) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2282  return (1.-0.834572);
2283  }
2284  if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.25 && fabs(eta2) < 1.6)) ||
2285  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.25 && fabs(eta1) < 1.6)) ) {
2286  return (1.-0.894739);
2287  }
2288  if( ((fabs(eta1) >= 0.85 && fabs(eta1) < 1.25) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1000)) ||
2289  ((fabs(eta2) >= 0.85 && fabs(eta2) < 1.25) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1000)) ) {
2290  return (1.-0.880597);
2291  }
2292  if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.6 && fabs(eta2) < 1.8)) ||
2293  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.6 && fabs(eta1) < 1.8)) ) {
2294  return (1.-0.892911);
2295  }
2296  if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 1.8 && fabs(eta2) < 2.0)) ||
2297  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 1.8 && fabs(eta1) < 2.0)) ) {
2298  return (1.-0.880506);
2299  }
2300  if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2)) ||
2301  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2)) ) {
2302  return (1.-0.885718);
2303  }
2304  if( ((fabs(eta1) >= 1.25 && fabs(eta1) < 1.6) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2305  ((fabs(eta2) >= 1.25 && fabs(eta2) < 1.6) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2306  return (1.-0.853141);
2307  }
2308  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 1.8 && fabs(eta2) < 2.0)) ||
2309  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) && (fabs(eta1) >= 1.8 && fabs(eta1) < 2.0)) ) {
2310  return (1.-0.88822);
2311  }
2312  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2)) ||
2313  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) && (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2)) ) {
2314  return (1.-0.87028);
2315  }
2316  if( ((fabs(eta1) >= 1.6 && fabs(eta1) < 1.8) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2317  ((fabs(eta2) >= 1.6 && fabs(eta2) < 1.8) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2318  return (1.-0.869603);
2319  }
2320  if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 2.0) && (fabs(eta2) >= 2.0 && fabs(eta2) < 2.2)) ||
2321  ((fabs(eta2) >= 1.8 && fabs(eta2) < 2.0) && (fabs(eta1) >= 2.0 && fabs(eta1) < 2.2)) ) {
2322  return (1.-0.877922);
2323  }
2324  if( ((fabs(eta1) >= 1.8 && fabs(eta1) < 2.0) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2325  ((fabs(eta2) >= 1.8 && fabs(eta2) < 2.0) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2326  return (1.-0.865997);
2327  }
2328  if( ((fabs(eta1) >= 2.0 && fabs(eta1) < 2.2) && (fabs(eta2) >= 2.2 && fabs(eta2) < 1000)) ||
2329  ((fabs(eta2) >= 2.0 && fabs(eta2) < 2.2) && (fabs(eta1) >= 2.2 && fabs(eta1) < 1000)) ) {
2330  return (1.-0.886109);
2331  }
2332  return 0;
2333  }
2334 };
2335 
2336 
2338 // --------------------------
2339 class backgroundFunctionType11 : public backgroundFunctionBase {
2340  public:
2341  backgroundFunctionType11(const double & lowerLimit, const double & upperLimit) :
2342  backgroundFunctionBase(lowerLimit, upperLimit)
2343  {
2344  this->parNum_ = 2;
2345  }
2346  double operator()( const double * parval, const double & mass, const double & eta ) const override {return 0.;}
2347  double operator()( const double * parval, const double & mass, const double & eta1, const double & eta2 ) const override
2348  {
2349  double Bgrp2 = 0.;
2350  if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= -100. && eta2 < -0.8) ) {
2351  Bgrp2 = (-0.0512353);
2352  }
2353  else if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= -0.8 && eta2 < 0.) ) {
2354  Bgrp2 = (-0.0448482);
2355  }
2356  else if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= 0. && eta2 < 0.8) ) {
2357  Bgrp2 = (-0.0193726);
2358  }
2359  else if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= 0.8 && eta2 < 100.) ) {
2360  Bgrp2 = (0.0225765);
2361  }
2362  else if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= -100. && eta2 < -0.8) ) {
2363  Bgrp2 = (-0.0822936);
2364  }
2365  else if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= -0.8 && eta2 < 0.) ) {
2366  Bgrp2 = (-0.0676357);
2367  }
2368  else if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= 0. && eta2 < 0.8) ) {
2369  Bgrp2 = (-0.0591544);
2370  }
2371  else if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= 0.8 && eta2 < 100.) ) {
2372  Bgrp2 = (-0.0235858);
2373  }
2374  else if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= -100. && eta2 < -0.8) ) {
2375  Bgrp2 = (-0.0317051);
2376  }
2377  else if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= -0.8 && eta2 < 0.) ) {
2378  Bgrp2 = (-0.06139);
2379  }
2380  else if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= 0. && eta2 < 0.8) ) {
2381  Bgrp2 = (-0.0747737);
2382  }
2383  else if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= 0.8 && eta2 < 100.) ) {
2384  Bgrp2 = (-0.0810139);
2385  }
2386  else if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= -100. && eta2 < -0.8) ) {
2387  Bgrp2 = (0.0229602);
2388  }
2389  else if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= -0.8 && eta2 < 0.) ) {
2390  Bgrp2 = (-0.0224212);
2391  }
2392  else if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= 0. && eta2 < 0.8) ) {
2393  Bgrp2 = (-0.0446273);
2394  }
2395  else if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= 0.8 && eta2 < 100.) ) {
2396  Bgrp2 = (-0.0554561);
2397  }
2398  else {
2399  std::cout << "WARNING, backgroundFunctionType11: this should not happen for eta1 = " << eta1 << " and eta2 = " << eta2 << std::endl;
2400  return (-0.05);
2401  }
2402  double norm = (exp(Bgrp2*upperLimit_) - exp(Bgrp2*lowerLimit_))/Bgrp2;
2403  if( norm != 0 ) return exp(Bgrp2*mass)/norm;
2404  else return 0.;
2405 
2406  }
2407  void setParameters(double* Start, double* Step, double* Mini, double* Maxi, int* ind, TString* parname, const std::vector<double>::const_iterator & parBgrIt, const std::vector<int>::const_iterator & parBgrOrderIt, const int muonType) override {
2408  double thisStep[] = {0.01, 0.01};
2409  TString thisParName[] = {"Bgr fraction", "Bgr slope"};
2410  if( muonType == 1 ) {
2411  double thisMini[] = {-1.0, 10.};
2412  double thisMaxi[] = {1.0 , 10.};
2413  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
2414  } else {
2415  double thisMini[] = {-1.0, 10.};
2416  double thisMaxi[] = { 1.0, 10.};
2417  this->setPar( Start, Step, Mini, Maxi, ind, parname, parBgrIt, parBgrOrderIt, thisStep, thisMini, thisMaxi, thisParName );
2418  }
2419  }
2420 
2421  double fracVsEta(const double * parval, const double & eta1, const double & eta2) const override
2422  {
2423  if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= -100. && eta2 < -0.8) ) {
2424  return (1.-0.966316);
2425  }
2426  if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= -0.8 && eta2 < 0.) ) {
2427  return (1.-0.966875);
2428  }
2429  if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= 0. && eta2 < 0.8) ) {
2430  return (1.-0.955311);
2431  }
2432  if( (eta1 >= -100. && eta1 < -0.8) && (eta2 >= 0.8 && eta2 < 100.) ) {
2433  return (1.-0.928771);
2434  }
2435  if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= -100. && eta2 < -0.8) ) {
2436  return (1.-0.983255);
2437  }
2438  if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= -0.8 && eta2 < 0.) ) {
2439  return (1.-0.982203);
2440  }
2441  if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= 0. && eta2 < 0.8) ) {
2442  return (1.-0.972127);
2443  }
2444  if( (eta1 >= -0.8 && eta1 < 0.) && (eta2 >= 0.8 && eta2 < 100.) ) {
2445  return (1.-0.962929);
2446  }
2447  if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= -100. && eta2 < -0.8) ) {
2448  return (1.-0.965597);
2449  }
2450  if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= -0.8 && eta2 < 0.) ) {
2451  return (1.-0.969461);
2452  }
2453  if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= 0. && eta2 < 0.8) ) {
2454  return (1.-0.979922);
2455  }
2456  if( (eta1 >= 0. && eta1 < 0.8) && (eta2 >= 0.8 && eta2 < 100.) ) {
2457  return (1.-0.984247);
2458  }
2459  if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= -100. && eta2 < -0.8) ) {
2460  return (1.-0.934252);
2461  }
2462  if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= -0.8 && eta2 < 0.) ) {
2463  return (1.-0.952914);
2464  }
2465  if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= 0. && eta2 < 0.8) ) {
2466  return (1.-0.960191);
2467  }
2468  if( (eta1 >= 0.8 && eta1 < 100.) && (eta2 >= 0.8 && eta2 < 100.) ) {
2469  return (1.-0.966175);
2470  }
2471  else {
2472  std::cout << "WARNING, backgroundFunctionType11: this should not happen for eta1 = " << eta1 << " and eta2 = " << eta2 << std::endl;
2473  return (1.-0.97);
2474  }
2475  }
2476 };
2477 
2478 
2479 #endif // FUNCTIONS_H
Pi
const double Pi
Definition: CosmicMuonParameters.h:18
backgroundFunctionType4::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1368
scaleFunctionType50::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parScale, const std::vector< int > &parScaleOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:167
backgroundFunctionType9::backgroundFunctionType9
backgroundFunctionType9(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:2020
scaleFunctionService
scaleFunctionBase< double * > * scaleFunctionService(const int identifier)
Service to build the scale functor corresponding to the passed identifier.
Definition: Functions.cc:3
resolutionFunctionType0::sigmaCotgTh
double sigmaCotgTh(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:752
backgroundFunctionType1::backgroundFunctionType1
backgroundFunctionType1(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1213
ParameterSet::name
TString name
Definition: Functions.h:27
Scenarios_cff.scale
scale
Definition: Scenarios_cff.py:2186
ParameterSet::mini
double mini
Definition: Functions.h:28
mps_fire.i
i
Definition: mps_fire.py:338
SigmaPtDiff
Returns ( sigmaPt/Pt(data)^2 - sigmaPt/Pt(MC)^2 )
Definition: SigmaPtDiff.h:71
backgroundFunctionType7
Exponential binned in eta, much finer binning then type6.
Definition: Functions.h:1519
scaleFunctionType64::scale
double scale(const double &pt, const double &eta, const double &phi, const int chg, const T &parScale) const override
Definition: Functions.h:268
resolutionFunctionType45::sigmaPtError
double sigmaPtError(const double &pt, const double &eta, const T &parval, const T &parError) override
Definition: Functions.h:802
backgroundFunctionType11::operator()
double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const override
Definition: Functions.h:2347
backgroundFunctionBase::functionForIntegral
virtual TF1 * functionForIntegral(const std::vector< double >::const_iterator &parBgrIt) const
Definition: Functions.h:1148
backgroundFunctionType8::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1912
backgroundFunctionBase::FunctionForIntegral::operator()
double operator()(const double *mass, const double *) const
Definition: Functions.h:1190
smearFunctionType6::etaByPoints
double etaByPoints(const double &inEta, const double &border)
Definition: Functions.h:623
smearFunctionType4::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:568
chg
const float chg[109]
Definition: CoreSimTrack.cc:5
backgroundFunctionType10::fracVsEta
double fracVsEta(const double *parval, const double &eta1, const double &eta2) const override
Definition: Functions.h:2249
resolutionFunctionType45::sigmaCotgTh
double sigmaCotgTh(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:794
backgroundFunctionType10::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:2235
smearFunctionBase::gRandom_
TRandom * gRandom_
Definition: Functions.h:519
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
resolutionFunctionType47::etaBin
int etaBin(const double &eta)
Definition: Functions.h:986
backgroundFunctionType5::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1413
resolutionFunctionType45::sigmaPhi
double sigmaPhi(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:797
backgroundFunctionType11::backgroundFunctionType11
backgroundFunctionType11(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:2341
theta
Geom::Theta< T > theta() const
Definition: Basic3DVectorLD.h:179
resolutionFunctionType0::sigmaPt
double sigmaPt(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:751
backgroundFunctionType7::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1526
backgroundFunctionType9::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:2076
scaleFunctionBase::setPar
virtual void setPar(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< ParameterSet > &parSet)
Definition: Functions.h:82
backgroundFunctionType7::fracVsEta
double fracVsEta(const double *parval, const double &eta1, const double &eta2) const override
Definition: Functions.h:1725
smearFunctionType7::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:646
resolutionFunctionBase::setPar
virtual void setPar(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, double *thisStep, double *thisMini, double *thisMaxi, TString *thisParName)
This method sets the parameters.
Definition: Functions.h:708
scaleFunctionBase::~scaleFunctionBase
virtual ~scaleFunctionBase()=0
Definition: Functions.h:103
backgroundFunctionType4::backgroundFunctionType4
backgroundFunctionType4(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1358
resolutionFunctionBase::setParameters
virtual void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType)
Definition: Functions.h:694
smearFunctionType0::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:527
smearFunctionBase::smearEta
void smearEta(double &eta)
Definition: Functions.h:508
backgroundFunctionType6::fracVsEta
double fracVsEta(const double *parval, const double &eta1, const double &eta2) const override
Definition: Functions.h:1487
resolutionFunctionType45::etaBin
int etaBin(const double &eta)
Definition: Functions.h:770
backgroundFunctionType5::backgroundFunctionType5
backgroundFunctionType5(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1396
backgroundFunctionType11::fracVsEta
double fracVsEta(const double *parval, const double &eta1, const double &eta2) const override
Definition: Functions.h:2421
smearFunctionType6::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:591
backgroundFunctionType2::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1255
backgroundFunctionType2::backgroundFunctionType2
backgroundFunctionType2(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1252
resolutionFunctionType45::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:808
scaleFunctionType0::resetParameters
void resetParameters(std::vector< double > *scaleVec) const override
This method is used to reset the scale parameters to neutral values (useful for iterations > 0) ...
Definition: Functions.h:116
backgroundFunctionType10::backgroundFunctionType10
backgroundFunctionType10(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:2140
backgroundFunctionBase::operator()
virtual double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const
Definition: Functions.h:1141
backgroundFunctionType1::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1229
backgroundFunctionBase::~backgroundFunctionBase
virtual ~backgroundFunctionBase()
Definition: Functions.h:1136
backgroundFunctionBase::FunctionForIntegral::function_
const backgroundFunctionBase * function_
Definition: Functions.h:1196
backgroundFunctionBase::fracVsEta
virtual double fracVsEta(const double *parval, const double &eta1, const double &eta2) const
Definition: Functions.h:1155
scaleFunctionType50::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parScale, const std::vector< int > &parScaleOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType) override
Definition: Functions.h:200
resolutionFunctionVecService
resolutionFunctionBase< std::vector< double > > * resolutionFunctionVecService(const int identifier)
Service to build the resolution functor corresponding to the passed identifier when receiving a std::...
Definition: Functions.cc:50
backgroundFunctionType4
Exponential with eta dependence.
Definition: Functions.h:1350
resolutionFunctionBase::parNum
virtual int parNum() const
Definition: Functions.h:704
backgroundFunctionType8::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1958
resolutionFunctionType45::sigmaPt
double sigmaPt(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:790
hiDetachedQuadStep_cff.pt2
pt2
Definition: hiDetachedQuadStep_cff.py:70
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:18
conversionPostprocessing_cfi.etaBin
etaBin
Definition: conversionPostprocessing_cfi.py:18
scaleFunctionBase::setPar
virtual void setPar(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parScale, const std::vector< int > &parScaleOrder, double *thisStep, double *thisMini, double *thisMaxi, TString *thisParName)
This method sets the parameters.
Definition: Functions.h:70
backgroundFunctionType9::fracVsEta
double fracVsEta(const double *parval, const double &eta1, const double &eta2) const override
Definition: Functions.h:2090
smearFunctionType2::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:544
scaleFunctionType0::scale
double scale(const double &pt, const double &eta, const double &phi, const int chg, const T &parScale) const override
Definition: Functions.h:115
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
resolutionFunctionBase::setPar
virtual void setPar(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< ParameterSet > &parSet)
Definition: Functions.h:720
scaleFunctionType64::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parScale, const std::vector< int > &parScaleOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType) override
Definition: Functions.h:409
funct::tan
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
backgroundFunctionType10::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:2145
resolutionFunctionBase::covPt1Pt2
virtual double covPt1Pt2(const double &pt1, const double &eta1, const double &pt2, const double &eta2, const T &parval)
Definition: Functions.h:685
backgroundFunctionBase::backgroundFunctionBase
backgroundFunctionBase(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1134
resolutionFunctionType46::sigmaPtError
double sigmaPtError(const double &pt, const double &eta, const T &parval, const T &parError) override
Definition: Functions.h:910
scaleFunctionBase::setParameters
virtual void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType)
Definition: Functions.h:56
scaleFunctionType50::resetParameters
void resetParameters(std::vector< double > *scaleVec) const override
This method is used to reset the scale parameters to neutral values (useful for iterations > 0) ...
Definition: Functions.h:161
resolutionFunctionType47::sigmaPtError
double sigmaPtError(const double &pt, const double &eta, const T &parval, const T &parError) override
Definition: Functions.h:1023
smearFunctionType1::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:533
ParameterSet::ParameterSet
ParameterSet(const TString &inputName, const double &inputStep, const double &inputMini, const double &inputMaxi)
Definition: Functions.h:19
smearFunctionService
smearFunctionBase * smearFunctionService(const int identifier)
Service to build the smearing functor corresponding to the passed identifier.
Definition: Functions.cc:23
backgroundFunctionType7::operator()
double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const override
Definition: Functions.h:1527
resolutionFunctionType46::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType) override
Definition: Functions.h:941
hiDetachedQuadStep_cff.pt1
pt1
Definition: hiDetachedQuadStep_cff.py:70
backgroundFunctionType6::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1438
backgroundFunctionType11::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:2407
resolutionFunctionType46::sigmaPhi
double sigmaPhi(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:905
ParameterSet::step
double step
Definition: Functions.h:28
resolutionFunctionType0::sigmaPtError
double sigmaPtError(const double &pt, const double &eta, const T &parval, const T &parError) override
Definition: Functions.h:755
scaleFunctionVecService
scaleFunctionBase< std::vector< double > > * scaleFunctionVecService(const int identifier)
Service to build the scale functor corresponding to the passed identifier when receiving a std::vecto...
Definition: Functions.cc:13
backgroundFunctionType5::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1399
resolutionFunctionType0::sigmaPhi
double sigmaPhi(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:753
scaleFunctionType0::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parScale, const std::vector< int > &parScaleOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:117
resolutionFunctionType0::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:756
resolutionFunctionType47::covPt1Pt2
double covPt1Pt2(const double &pt1, const double &eta1, const double &pt2, const double &eta2, const T &parval) override
Definition: Functions.h:1017
resolutionFunctionType46::sigmaCotgTh
double sigmaCotgTh(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:902
resolutionFunctionType47::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1030
backgroundFunctionType4::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1361
backgroundFunctionType9::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:2025
backgroundFunctionType8::fracVsEta
double fracVsEta(const double *parval, const double &eta1, const double &eta2) const override
Definition: Functions.h:1972
b
double b
Definition: hdecay.h:120
scaleFunctionType64::resetParameters
void resetParameters(std::vector< double > *scaleVec) const override
This method is used to reset the scale parameters to neutral values (useful for iterations > 0) ...
Definition: Functions.h:344
backgroundFunctionBase::FunctionForIntegral::FunctionForIntegral
FunctionForIntegral(const backgroundFunctionBase *function, const std::vector< double >::const_iterator &parBgrIt)
Definition: Functions.h:1177
resolutionFunctionType46::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:917
ParameterSet::maxi
double maxi
Definition: Functions.h:28
smearFunctionType5::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:579
smearFunctionType7::sigmaPtDiff
SigmaPtDiff sigmaPtDiff
Definition: Functions.h:655
backgroundFunctionType2::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1262
backgroundFunctionBase::functionForIntegral_
FunctionForIntegral * functionForIntegral_
Definition: Functions.h:1199
resolutionFunctionType46::etaBin
int etaBin(const double &eta)
Definition: Functions.h:878
tmp
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
backgroundFunctionType5
Linear with eta dependence.
Definition: Functions.h:1389
backgroundFunctionBase::parNum
virtual int parNum() const
Definition: Functions.h:1145
a
double a
Definition: hdecay.h:121
backgroundFunctionType10::operator()
double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const override
Definition: Functions.h:2146
scaleFunctionType64::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parScale, const std::vector< int > &parScaleOrder, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:350
backgroundFunctionType1::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:1216
resolutionFunctionService
resolutionFunctionBase< double * > * resolutionFunctionService(const int identifier)
Service to build the resolution functor corresponding to the passed identifier.
Definition: Functions.cc:38
backgroundFunctionType11
Exponential binned in eta (Z, Run2012C PromptReco-v1 + PromptReco-v2)
Definition: Functions.h:2339
backgroundFunctionType6::operator()
double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const override
Definition: Functions.h:1439
resolutionFunctionBase::setParameters
virtual void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const int muonType)
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:692
resolutionFunctionType47::sigmaCotgTh
double sigmaCotgTh(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:1010
resolutionFunctionType45::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType) override
Definition: Functions.h:832
backgroundFunctionType7::backgroundFunctionType7
backgroundFunctionType7(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1521
backgroundFunctionType6::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1473
resolutionFunctionType47::sigmaPt
double sigmaPt(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:1006
resolutionFunctionType0::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType) override
Definition: Functions.h:757
gather_cfg.cout
cout
Definition: gather_cfg.py:144
step
step
Definition: StallMonitor.cc:94
backgroundFunctionType11::operator()
double operator()(const double *parval, const double &mass, const double &eta) const override
Definition: Functions.h:2346
JetChargeProducer_cfi.exp
exp
Definition: JetChargeProducer_cfi.py:6
backgroundFunctionType9::operator()
double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const override
Definition: Functions.h:2026
resolutionFunctionType47::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const T &parResol, const std::vector< int > &parResolOrder, const std::vector< double > &parStep, const std::vector< double > &parMin, const std::vector< double > &parMax, const int muonType) override
Definition: Functions.h:1055
scaleFunctionBase::resetParameters
virtual void resetParameters(std::vector< double > *scaleVec) const
This method is used to reset the scale parameters to neutral values (useful for iterations > 0) ...
Definition: Functions.h:49
smearFunctionBase::~smearFunctionBase
virtual ~smearFunctionBase()=0
Definition: Functions.h:521
backgroundFunctionService
backgroundFunctionBase * backgroundFunctionService(const int identifier, const double &lowerLimit, const double &upperLimit)
Service to build the background functor corresponding to the passed identifier.
Definition: Functions.cc:62
T
long double T
Definition: Basic3DVectorLD.h:66
resolutionFunctionType46::sigmaPt
double sigmaPt(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:898
backgroundFunctionType8::operator()
double operator()(const double *parval, const double &mass, const double &eta1, const double &eta2) const override
Definition: Functions.h:1913
scaleFunctionBase::parNum
virtual int parNum() const
Definition: Functions.h:66
resolutionFunctionBase::sigmaPtError
virtual double sigmaPtError(const double &pt, const double &eta, const T &parval, const T &parError)
Definition: Functions.h:679
backgroundFunctionBase::setPar
virtual void setPar(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, double *thisStep, double *thisMini, double *thisMaxi, TString *thisParName)
This method sets the parameters.
Definition: Functions.h:1162
backgroundFunctionType7::setParameters
void setParameters(double *Start, double *Step, double *Mini, double *Maxi, int *ind, TString *parname, const std::vector< double >::const_iterator &parBgrIt, const std::vector< int >::const_iterator &parBgrOrderIt, const int muonType) override
This method is used to differentiate parameters among the different functions.
Definition: Functions.h:1711
resolutionFunctionBase::~resolutionFunctionBase
virtual ~resolutionFunctionBase()=0
Definition: Functions.h:736
smearFunctionType3::smear
void smear(double &pt, double &eta, double &phi, const double *y, const std::vector< double > &parSmear) override
Definition: Functions.h:555
backgroundFunctionType6::backgroundFunctionType6
backgroundFunctionType6(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1433
backgroundFunctionType8::backgroundFunctionType8
backgroundFunctionType8(const double &lowerLimit, const double &upperLimit)
Definition: Functions.h:1907
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
scaleFunctionType50::scale
double scale(const double &pt, const double &eta, const double &phi, const int chg, const T &parScale) const override
Definition: Functions.h:127
resolutionFunctionType47::sigmaPhi
double sigmaPhi(const double &pt, const double &eta, const T &parval) override
Definition: Functions.h:1013
backgroundFunctionType6
Exponential binned in eta.
Definition: Functions.h:1431
Generated for CMSSW Reference Manual by doxygen 1.8.11