ATLAS Offline Software
Functions
muCombUtil Namespace Reference

Functions

void setMuFastRes (std::vector< double > &, double, double, double, double, double, double)
 utility vector set More...
 
void setIDSCANRes (std::vector< double > &, double, double)
 utility vector set More...
 
double getMuFastRes (const std::vector< double > &vec, const MuonFeature *feature)
 
int whichECRegion (const float, const float)
 utility function for getMuFastRes More...
 
double getIDSCANRes (std::vector< double >, std::vector< double >, std::vector< double >, std::vector< double >, std::vector< double >, double, double)
 Get parametrized IDSCAN 1/pt resolution. More...
 
double getMuFastEtaRes (const MuonFeature *)
 Get parametrized muFast Eta resolution (extrapolated) More...
 
double getMuFastPhiRes (const MuonFeature *)
 Get parametrized muFast Phi resolution (extrapolated) More...
 
double getG4ExtEtaRes (double, double)
 Get parametrized Geant4 Eta resolution (extrapolated) More...
 
double getG4ExtPhiRes (double, double)
 Get parametrized Geant4 Phi resolution (extrapolated) More...
 
double getDeltaPhi (double, double)
 Get DeltaPhi. More...
 
double getDeltaEta (double, double)
 Get DeltaEta. More...
 
double getDeltaR (double, double, double, double)
 Get DeltaR. More...
 
double getCombinedAverage (double, double, double, double)
 Get weighted mean. More...
 
double getCombinedAverageSigma (double, double)
 Get sigma of weighted mean. More...
 
double getChi2 (int &, double, double, double, double, double, double, double, double, double, double, double, double, double, bool)
 Get OLD style (i.e. muFast time) Chi2. More...
 
double getStdChi2 (int &, double, double, double, double, double, double, double, double, double, double, double, double, bool)
 Get Std Chi2. More...
 
double getMuFastRes (std::vector< double >, const MuonFeature *)
 Get parametrized muFast 1/pt resolution (extrapolated) More...
 

Function Documentation

◆ getChi2()

double muCombUtil::getChi2 ( int &  ndof,
double  ipt,
double  eta1,
double  seta1,
double  phi1,
double  sphi1,
double  ipt1,
double  sipt1,
double  eta2,
double  seta2,
double  phi2,
double  sphi2,
double  ipt2,
double  sipt2,
bool  useAbsPt 
)

Get OLD style (i.e. muFast time) Chi2.

Definition at line 374 of file muCombUtil.cxx.

376  {
377 
378  double deta = getDeltaEta(eta1,eta2);
379  double sdeta = seta1*seta1+seta2*seta2;
380  double dphi = getDeltaPhi(phi1,phi2);
381  double sdphi = sphi1*sphi1+sphi2*sphi2;
382  double dipt_1 = ipt - ipt1;
383  if (useAbsPt) dipt_1 = fabs(ipt) - fabs(ipt1);
384  double sdipt_1 = sipt1*sipt1;
385  double dipt_2 = ipt - ipt2;
386  if (useAbsPt) dipt_2 = fabs(ipt) - fabs(ipt2);
387  double sdipt_2 = sipt2*sipt2;
388 
389  double chi2 = 0.0;
390  ndof = 0;
391  if (sdeta != 0) { chi2 += deta*deta/sdeta; ndof++; }
392  if (sdphi != 0) { chi2 += dphi*dphi/sdphi; ndof++; }
393  if (sdipt_1 != 0) { chi2 += dipt_1*dipt_1/sdipt_1; ndof++; }
394  if (sdipt_2 != 0) { chi2 += dipt_2*dipt_2/sdipt_2; ndof++; }
395 
396  if (ndof == 0) return 1.0e30;
397  else return chi2;
398  }

◆ getCombinedAverage()

double muCombUtil::getCombinedAverage ( double  p1,
double  sp1,
double  p2,
double  sp2 
)

Get weighted mean.

Definition at line 362 of file muCombUtil.cxx.

362  {
363  if (sp1 != 0 && sp2 != 0) return (sp2*sp2*fabs(p1) + sp1*sp1*fabs(p2))/(sp1*sp1 + sp2*sp2);
364  else if (sp1 != 0 && sp2 == 0) return fabs(p2);
365  else if (sp1 == 0 && sp2 != 0) return fabs(p1);
366  else return (fabs(p1)+fabs(p2))*0.5;
367  }

◆ getCombinedAverageSigma()

double muCombUtil::getCombinedAverageSigma ( double  sp1,
double  sp2 
)

Get sigma of weighted mean.

Definition at line 369 of file muCombUtil.cxx.

369  {
370  if (sp1 != 0 && sp2 != 0) return sqrt((sp1*sp1*sp2*sp2)/(sp1*sp1 + sp2*sp2));
371  else return 0.0;
372  }

◆ getDeltaEta()

double muCombUtil::getDeltaEta ( double  eta1,
double  eta2 
)

Get DeltaEta.

Definition at line 352 of file muCombUtil.cxx.

352  {
353  return fabs(eta1-eta2);
354  }

◆ getDeltaPhi()

double muCombUtil::getDeltaPhi ( double  phi1,
double  phi2 
)

Get DeltaPhi.

Definition at line 345 of file muCombUtil.cxx.

345  {
346  double dphi = phi1 - phi2;
347  if (dphi > M_PI) dphi -= 2*M_PI;
348  if (dphi < -M_PI) dphi += 2*M_PI;
349  return fabs(dphi);
350  }

◆ getDeltaR()

double muCombUtil::getDeltaR ( double  eta1,
double  phi1,
double  eta2,
double  phi2 
)

Get DeltaR.

Definition at line 356 of file muCombUtil.cxx.

356  {
357  double deta = getDeltaEta(eta1,eta2);
358  double dphi = getDeltaPhi(phi1,phi2);
359  return sqrt(deta*deta + dphi*dphi);
360  }

◆ getG4ExtEtaRes()

double muCombUtil::getG4ExtEtaRes ( double  pt,
double  eta 
)

Get parametrized Geant4 Eta resolution (extrapolated)

Definition at line 288 of file muCombUtil.cxx.

288  {
289 
290  //double pt = feature->pt();
291  //double eta = feature->eta();
292  if (pt < 4. ) pt = 4.;
293  if (pt > 40.) pt = 40.;
294  // bool ts = false;
295  // if (feature->radius() <= 10.) ts = true;
296 
297  //G4 extrapolator eta resolution parametrized as A/pt^2 + B/pt + C + D*pt
298  // Resolution evaluated in regions of eta: 0.0,0.5,1.0,1.5,2.0,above...
299  if (fabs(eta) < 0.5) {
300  return 6.0e-3/(pt*pt) + 4.4e-2/pt + 1.2e-2 + (-3.2e-5*pt);
301  }
302  else if (fabs(eta) >= 0.5 && fabs(eta) < 1.0) {
303  return 5.6e-3/(pt*pt) + 3.9e-2/pt + 1.3e-2 + (-3.7e-5*pt);
304  }
305  else if (fabs(eta) >= 1.0 && fabs(eta) < 1.5) {
306  return 2.1e-1/(pt*pt) + 8.8e-2/pt + 7.1e-4 + (1.7e-4*pt);
307  }
308  else if (fabs(eta) >= 1.5 && fabs(eta) < 2.0) {
309  return 4.5e-3/(pt*pt) + 1.4e-2/pt + 1.1e-2 + (-5.7e-5*pt);
310  }
311  else { // if (fabs(eta) >= 2.0) {
312  return 3.1e-2/(pt*pt) + 4.9e-1/pt + (-2.6e-3) + (3.2e-5*pt);
313  }
314  }

◆ getG4ExtPhiRes()

double muCombUtil::getG4ExtPhiRes ( double  pt,
double  eta 
)

Get parametrized Geant4 Phi resolution (extrapolated)

Definition at line 317 of file muCombUtil.cxx.

317  {
318 
319  //double pt = feature->pt();
320  //double eta = feature->eta();
321  if (pt < 4. ) pt = 4.;
322  if (pt > 40.) pt = 40.;
323  // bool ts = false;
324  // if (feature->radius() <= 10.) ts = true;
325 
326  //G4 extrapolator phi resolution parametrized as A/pt^2 + B/pt + C + D*pt
327  // Resolution evaluated in regions of eta: 0.0,0.5,1.0,1.5,2.0,above...
328  if (fabs(eta) < 0.5) {
329  return 4.1e-2/(pt*pt) + 1.5e-1/pt + 1.4e-3 + (5.8e-5*pt);
330  }
331  else if (fabs(eta) >= 0.5 && fabs(eta) < 1.0) {
332  return 9.4e-3/(pt*pt) + 8.8e-2/pt + 1.1e-2 + (-1.e-4*pt);
333  }
334  else if (fabs(eta) >= 1.0 && fabs(eta) < 1.5) {
335  return 3.5e-1/(pt*pt) + 9.4e-2/pt + 8.2e-3 + (-3.3e-7*pt);
336  }
337  else if (fabs(eta) >= 1.5 && fabs(eta) < 2.0) {
338  return 2.6e-2/(pt*pt) + 1.2e-1/pt + 3.8e-3 + (2.2e-5*pt);
339  }
340  else { // if (fabs(eta) >= 2.0) {
341  return 2.7e-1/(pt*pt) + 3.2e-2/pt + (9.0e-3) + (-2.6e-5*pt);
342  }
343  }

◆ getIDSCANRes()

double muCombUtil::getIDSCANRes ( std::vector< double >  barrelvec,
std::vector< double >  ec1vec,
std::vector< double >  ec2vec,
std::vector< double >  ec3vec,
std::vector< double >  ec4vec,
double  pt_id,
double  eta_id 
)

Get parametrized IDSCAN 1/pt resolution.

Definition at line 190 of file muCombUtil.cxx.

192  {
193 
194  double pt = pt_id;
195  if (pt == 0) return 1.0e33;
196 
197  double AbsPtInv = fabs(1./pt);
198  double AbsEta = fabs(eta_id);
199 
200  std::vector<double> vec;
201  if (AbsEta < 1.05) vec = std::move(barrelvec);
202  else if (AbsEta>=1.05 && AbsEta<1.35) vec = std::move(ec1vec);
203  else if (AbsEta>=1.35 && AbsEta<1.65) vec = std::move(ec2vec);
204  else if (AbsEta>=1.65 && AbsEta<2.0) vec = std::move(ec3vec);
205  else vec = std::move(ec4vec);
206 
207  return vec[0]*AbsPtInv+vec[1]/1000.;
208  }

◆ getMuFastEtaRes()

double muCombUtil::getMuFastEtaRes ( const MuonFeature feature)

Get parametrized muFast Eta resolution (extrapolated)

Definition at line 211 of file muCombUtil.cxx.

211  {
212 
213  double pt = feature->pt();
214  double eta = feature->eta();
215  if (pt < 4. ) pt = 4.;
216  if (pt > 40.) pt = 40.;
217  bool ts = false;
218  if (feature->radius() <= 10.) ts = true;
219 
220  if (fabs(eta) < 1.) {//barrel
221  return 7.75e-2/pt + 8.1e-3;
222  }
223  else {
224  if (ts) {//trigger station
225  if (fabs(eta) >= 1. && fabs(eta) < 1.2) {//ec0
226  return 2.0e-1/pt + 1.5e-3;
227  }
228  else if (fabs(eta) >= 1.2 && fabs(eta) < 1.5) {//ec1
229  return 2.5e-1/pt + 0.1e-3;
230  }
231  else if (fabs(eta) >= 1.5 && fabs(eta) < 2.0) {//ec2
232  return 1.3e-1/pt + 6.0e-3;
233  }
234  else {//ec3
235  return 1.3e-1/pt + 8.0e-3;
236  }
237  }
238  else {
239  if (fabs(eta) >= 1. && fabs(eta) < 1.6) {//ec0/1
240  return 0.071/pt + 0.0055;
241  }
242  else {//ec2/3
243  return 0.055/pt + 0.0063;
244  }
245  }
246  }
247  }

◆ getMuFastPhiRes()

double muCombUtil::getMuFastPhiRes ( const MuonFeature feature)

Get parametrized muFast Phi resolution (extrapolated)

Definition at line 250 of file muCombUtil.cxx.

250  {
251 
252  double pt = feature->pt();
253  double eta = feature->eta();
254  if (pt < 4. ) pt = 4.;
255  if (pt > 40.) pt = 40.;
256  bool ts = false;
257  if (feature->radius() <= 10.) ts = true;
258 
259  if (fabs(eta) < 1.) {//barrel
260  return 1.9e-1/pt - 3.4e-4;
261  }
262  else {
263  if (ts) {//trigger station
264  if (fabs(eta) >= 1. && fabs(eta) < 1.2) {//ec0
265  return 2.0e-1/pt + 2.0e-3;
266  }
267  else if (fabs(eta) >= 1.2 && fabs(eta) < 1.5) {//ec1
268  return 1.8e-1/pt + 3.0e-3;
269  }
270  else if (fabs(eta) >= 1.5 && fabs(eta) < 2.0) {//ec2
271  return 1.5e-1/pt + 8.0e-3;
272  }
273  else {//ec3
274  return 1.0e-1/pt + 1.8e-2;
275  }
276  }
277  else {
278  //if (fabs(eta) >= 1. && fabs(eta) < 1.6) {//ec0/1
279  return 0.08/pt + 0.0025;
280  //}
281  //else {//ec2/3
282  // return 0.08/pt + 0.0025;
283  //}
284  }
285  }
286  }

◆ getMuFastRes() [1/2]

double muCombUtil::getMuFastRes ( const std::vector< double > &  vec,
const MuonFeature feature 
)

Definition at line 42 of file muCombUtil.cxx.

42  {
43 
44  double ptGev = feature->pt();
45  double pt = ptGev*1000.; //muFast Pt(MeV)
46 
47  if (pt == 0) return 1.0e33;
48 
49  double AbsPtInv = fabs(1./pt); //muFast 1/Pt (1/MeV)
50  int add = feature->saddress();
51  double AbsEta = fabs(feature->eta());
52  double phi = feature->phi();
53 
54  if ( add != -1) {
55  if (AbsPtInv < 0.000186) {
56  return vec[0]*AbsPtInv + vec[1]/1000.;
57  }
58  else {
59  double AbsPtInv3 = AbsPtInv*AbsPtInv*AbsPtInv;
60  double AbsPtInv2 = AbsPtInv*AbsPtInv;
61  return vec[2]*AbsPtInv3/(1000.*1000.) +
62  vec[3]*AbsPtInv2/(1000.) +
63  vec[4]*AbsPtInv +
64  vec[5]/1000.;
65  }
66  }
67  else {//Takuya/Kunihiro updated numbers
68 
69  const int N_PARAMS = 5;
70  const double vparEC1[N_PARAMS] = {0.291483, -6.11348, 65.1099, -285.664, 440.041};
71  const double vparEC2[N_PARAMS] = {0.286307, -4.6759, 43.2815, -163.185, 210.786};
72  const double vparEC3[N_PARAMS] = {0.330699, -6.70755, 70.4725, -291.85, 408.739};
73  const double vparEC4[N_PARAMS] = {0.261738, -4.69971, 47.4762, -183.98, 236.813};
74  const double vparEC5[N_PARAMS] = {0.196301, -3.57276, 38.3744, -159.808, 228.256};
75  const double vparEC6[N_PARAMS] = {0.172939, -3.10788, 33.3823, -142.996, 212.957};
76  const double vparEC7[N_PARAMS] = {0.233017, -4.377, 42.5691, -171.752, 245.702};
77  const double vparEC8[N_PARAMS] = {0.22389, -4.16259, 40.1369, -162.824, 236.39};
78  const double vparEC9[N_PARAMS] = {0.197992, -3.52117, 33.5997, -136.014, 197.474};
79  const double vparECA[N_PARAMS] = {0.417289, -0.852254,-31.9257, 308.873, -719.591};
80  const double vparECB[N_PARAMS] = {0.526612, -8.04087, 82.1906, -336.87, 462.973};
81 
82  double AbsPtInvGeV = AbsPtInv * 1000;
83 
84  const double AbsPtInvGeVMin = 5e-3; // 200 GeV
85  const double AbsPtInvGeVMax = 0.25; // 4 GeV
86  if( AbsPtInvGeV < AbsPtInvGeVMin ) AbsPtInvGeV = AbsPtInvGeVMin;
87  if( AbsPtInvGeV > AbsPtInvGeVMax ) AbsPtInvGeV = AbsPtInvGeVMax;
88 
89  //cout << "...AbsPtInvGeV=" << AbsPtInv << endl;
90 
91  const double* vpar;
92  int spReg = whichECRegion(AbsEta,phi);
93  //cout << "...spReg=" << spReg << endl;
94  if ( spReg==1 ) { vpar = vparECA; }
95  else if( spReg==2 ) { vpar = vparECB; }
96  else {
97  if ( AbsEta < 1.20) { vpar = vparEC1; }
98  else if( AbsEta < 1.35) { vpar = vparEC2; }
99  else if( AbsEta < 1.50) { vpar = vparEC3; }
100  else if( AbsEta < 1.65) { vpar = vparEC4; }
101  else if( AbsEta < 1.80) { vpar = vparEC5; }
102  else if( AbsEta < 1.95) { vpar = vparEC6; }
103  else if( AbsEta < 2.10) { vpar = vparEC7; }
104  else if( AbsEta < 2.35) { vpar = vparEC8; }
105  else { vpar = vparEC9; }
106  }
107 
108  double fracRes = vpar[0] + vpar[1]*AbsPtInvGeV
109  + vpar[2]*AbsPtInvGeV*AbsPtInvGeV
110  + vpar[3]*AbsPtInvGeV*AbsPtInvGeV*AbsPtInvGeV
111  + vpar[4]*AbsPtInvGeV*AbsPtInvGeV*AbsPtInvGeV*AbsPtInvGeV;
112 
113  return fabs(fracRes * AbsPtInv);
114  }
115 }

◆ getMuFastRes() [2/2]

double muCombUtil::getMuFastRes ( std::vector< double >  ,
const MuonFeature  
)

Get parametrized muFast 1/pt resolution (extrapolated)

◆ getStdChi2()

double muCombUtil::getStdChi2 ( int &  ndof,
double  eta1,
double  seta1,
double  phi1,
double  sphi1,
double  qOvpt1,
double  sqOvpt1,
double  eta2,
double  seta2,
double  phi2,
double  sphi2,
double  qOvpt2,
double  sqOvpt2,
bool  useAbsPt 
)

Get Std Chi2.

Definition at line 401 of file muCombUtil.cxx.

403  {
404 
405  double deta = getDeltaEta(eta1,eta2);
406  double sdeta = seta1*seta1+seta2*seta2;
407  double dphi = getDeltaPhi(phi1,phi2);
408  double sdphi = sphi1*sphi1+sphi2*sphi2;
409  double dipt = qOvpt1 - qOvpt2;
410  if (useAbsPt) dipt = fabs(qOvpt1) - fabs(qOvpt2);
411  double sdipt = sqOvpt1*sqOvpt1+sqOvpt2*sqOvpt2;
412 
413  double chi2 = 0.0;
414  ndof = 0;
415  if (sdeta != 0) { chi2 += deta*deta/sdeta; ndof++; }
416  if (sdphi != 0) { chi2 += dphi*dphi/sdphi; ndof++; }
417  if (sdipt != 0) { chi2 += dipt*dipt/sdipt; ndof++; }
418 
419  if (ndof == 0) return 1.0e30;
420  else return chi2;
421  }

◆ setIDSCANRes()

void muCombUtil::setIDSCANRes ( std::vector< double > &  vec,
double  p1,
double  p2 
)

utility vector set

Definition at line 35 of file muCombUtil.cxx.

35  {
36  vec.clear();
37  vec.push_back(p1);
38  vec.push_back(p2);
39  }

◆ setMuFastRes()

void muCombUtil::setMuFastRes ( std::vector< double > &  vec,
double  p1,
double  p2,
double  p3,
double  p4,
double  p5,
double  p6 
)

utility vector set

Definition at line 24 of file muCombUtil.cxx.

25  {
26  vec.clear();
27  vec.push_back(p1);
28  vec.push_back(p2);
29  vec.push_back(p3);
30  vec.push_back(p4);
31  vec.push_back(p5);
32  vec.push_back(p6);
33  }

◆ whichECRegion()

int muCombUtil::whichECRegion ( const float  eta,
const float  phi 
)

utility function for getMuFastRes

Definition at line 117 of file muCombUtil.cxx.

117  {
118  // 0: bulk
119  // 1: WeakBfield A
120  // 2: WeakBfield B
121 
122  float absEta = fabs(eta);
123 
124  if( ( 1.3 <= absEta && absEta < 1.45) &&
125  ( (0 <= fabs(phi) && fabs(phi) < M_PI/48. ) ||
126  (M_PI*11./48. <= fabs(phi) && fabs(phi) < M_PI*13./48. ) ||
127  (M_PI*23./48. <= fabs(phi) && fabs(phi) < M_PI*25./48. ) ||
128  (M_PI*35./48. <= fabs(phi) && fabs(phi) < M_PI*37./48. ) ||
129  (M_PI*47./48. <= fabs(phi) && fabs(phi) < M_PI )
130  )
131  ) return 1;
132 
133  else if( ( 1.5 <= absEta && absEta < 1.65 ) &&
134  ( (M_PI*3./32. <= fabs(phi) && fabs(phi) < M_PI*5./32. ) ||
135  (M_PI*11./32. <= fabs(phi) && fabs(phi) < M_PI*13./32.) ||
136  (M_PI*19./32. <= fabs(phi) && fabs(phi) < M_PI*21./32.) ||
137  (M_PI*27./32. <= fabs(phi) && fabs(phi) < M_PI*29./32.)
138  )
139  ) return 2;
140 
141  else return 0;
142 }
IDTPM::ndof
float ndof(const U &p)
Definition: TrackParametersHelper.h:142
ParticleGun_SamplingFraction.eta2
eta2
Definition: ParticleGun_SamplingFraction.py:96
phi
Scalar phi() const
phi method
Definition: AmgMatrixBasePlugin.h:64
muCombUtil::getDeltaEta
double getDeltaEta(double eta1, double eta2)
Get DeltaEta.
Definition: muCombUtil.cxx:352
MuonFeature::radius
float radius(void) const
Definition: MuonFeature.h:49
eta
Scalar eta() const
pseudorapidity method
Definition: AmgMatrixBasePlugin.h:79
test_pyathena.pt
pt
Definition: test_pyathena.py:11
MuonFeature::eta
float eta(void) const
Definition: MuonFeature.h:50
M_PI
#define M_PI
Definition: ActiveFraction.h:11
xAOD::eta1
setEt setPhi setE277 setWeta2 eta1
Definition: TrigEMCluster_v1.cxx:41
vec
std::vector< size_t > vec
Definition: CombinationsGeneratorTest.cxx:12
MuonFeature::phi
float phi(void) const
Definition: MuonFeature.h:51
muCombUtil::getDeltaPhi
double getDeltaPhi(double phi1, double phi2)
Get DeltaPhi.
Definition: muCombUtil.cxx:345
ParticleJetTools::p3
Amg::Vector3D p3(const xAOD::TruthVertex *p)
Definition: ParticleJetLabelCommon.cxx:55
CxxUtils::vec
typename vecDetail::vec_typedef< T, N >::type vec
Define a nice alias for the vectorized type.
Definition: vec.h:207
MuonFeature::saddress
int saddress(void) const
Definition: MuonFeature.h:47
muCombUtil::whichECRegion
int whichECRegion(const float eta, const float phi)
utility function for getMuFastRes
Definition: muCombUtil.cxx:117
chi2
double chi2(TH1 *h0, TH1 *h1)
Definition: comparitor.cxx:522
add
bool add(const std::string &hname, TKey *tobj)
Definition: fastadd.cxx:55
DiTauMassTools::MaxHistStrategyV2::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:26
TauGNNUtils::Variables::absEta
bool absEta(const xAOD::TauJet &tau, double &out)
Definition: TauGNNUtils.cxx:232
MuonFeature::pt
float pt(void) const
Definition: MuonFeature.h:48
python.CaloScaleNoiseConfig.ts
ts
Definition: CaloScaleNoiseConfig.py:86
Prompt::Def::AbsEta
@ AbsEta
Definition: VarHolder.h:77