ATLAS Offline Software
Reconstruction
egamma
egammaLayerRecalibTool
Root
corr_HV_EMBPS.cxx
Go to the documentation of this file.
1
/*
2
Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
3
*/
4
5
#include "
egammaLayerRecalibTool/corr_HV_EMBPS.h
"
6
#include <cmath>
7
#include <cstdio>
8
#include <cstdlib>
9
#include <iostream>
10
#include <string>
11
12
#include <TProfile2D.h>
13
#include <TFile.h>
14
#include <TSystem.h>
15
16
#include "
PathResolver/PathResolver.h
"
17
18
19
corr_HV_EMBPS::corr_HV_EMBPS
(){
20
std::string
filename
=
PathResolverFindCalibFile
(
"egammaLayerRecalibTool/v1/HVmaps_embps_2012.root"
);
21
m_file
= TFile::Open(
filename
.c_str());
22
if
(not
m_file
or
m_file
->IsZombie()) {
23
std::cerr <<
"FATAL: cannot open "
<<
filename
<< std::endl;
24
}
25
else
{
26
for
(
int
iperiod=0;iperiod<6;iperiod++) {
27
for
(
int
igap=0;igap<2;igap++) {
28
char
name
[12];
29
snprintf(
name
,
sizeof
(
name
),
"h%d%d"
,iperiod+1,igap);
30
m_hHV
[iperiod][igap]=(
TProfile2D
*) (
m_file
->Get(
name
));
31
if
(not
m_hHV
[iperiod][igap]) {
32
std::cerr <<
"FATAL: cannot find "
<<
name
<< std::endl;
33
}
34
}
35
}
36
}
37
}
38
//===============================================================================
39
corr_HV_EMBPS::~corr_HV_EMBPS
()
40
{
41
m_file
->Close();
42
}
43
44
//===============================================================================
45
float
corr_HV_EMBPS::getCorr
(
int
run
,
float
eta
,
float
phi
)
const
46
{
47
if
(std::fabs(
eta
) > 1.5)
return
1.;
48
if
(
run
<200804)
return
1.;
49
int
iperiod;
50
if
(
run
<204932) iperiod=0;
51
else
if
(
run
<205112) iperiod=1;
52
else
if
(
run
<211670) iperiod=2;
53
else
if
(
run
<212619) iperiod=3;
54
else
if
(
run
<212809) iperiod=4;
55
else
iperiod=5;
56
//std::cout << "period " << iperiod << std::endl;
57
58
int
ieta=(
int
)((
eta
+1.6)/0.4);
59
if
(ieta<0 || ieta>7)
return
1.;
60
if
(
phi
<0.)
phi
=
phi
+2*
M_PI
;
61
int
iphi=(
int
)(
phi
*(32./(2*
M_PI
)));
62
if
(iphi<0) iphi=0;
63
if
(iphi>31) iphi=31;
64
//std::cout << " ieta,iphi " << ieta << " " << iphi << std::endl;
65
66
float
hv1 =
m_hHV
[iperiod][0]->
GetBinContent
(ieta+1,iphi+1);
67
float
hv2 =
m_hHV
[iperiod][1]->
GetBinContent
(ieta+1,iphi+1);
68
69
float
c1
=
getRecoCorrection
(hv1,
eta
);
70
float
c2
=
getRecoCorrection
(hv2,
eta
);
71
float
delta1 =
getDataCorrection
(hv1,
eta
);
72
float
delta2 =
getDataCorrection
(hv2,
eta
);
73
74
float
extra1=1.;
75
if
(
run
<211670 && hv1<1300.) extra1 =
getExtraScaleFactor
();
76
float
extra2=1.;
77
if
(
run
<211670 && hv2<1300.) extra2 =
getExtraScaleFactor
();
78
79
if
((
c1
+
c2
)<0.01)
return
1.;
80
float
oldScale = (
c1
+
c2
);
81
float
newScale = (
c1
*delta1*extra1 +
c2
*delta2*extra2);
82
83
float
newCorr = oldScale/newScale;
84
85
//std::cout << " run,eta,phi " << run << " " << eta << " " << phi << " " << " period "<< iperiod << std::endl;
86
//std::cout << " HV values " << hv1 << " " << hv2 << std::endl;
87
//std::cout << " reco scales " << c1 << " " << c2 << std::endl;
88
//std::cout << " delta scale " << delta1 << " " << delta2 << std::endl;
89
//std::cout << " extra " << extra1 << " " << extra2 << std::endl;
90
//std::cout << " old/new scales " << oldScale << " " << newScale << std::endl;
91
//std::cout << " newCorr " << newCorr << std::endl;
92
93
return
newCorr;
94
95
}
96
97
//===============================================================================
98
// return scale factor of response vs HV used in reconstruction
99
float
corr_HV_EMBPS::getRecoCorrection
(
float
hv,
float
eta
)
100
{
101
float
nominal
= 2000.;
102
float
T
= 88.37;
103
104
int
ieta=(std::fabs(
eta
)/0.025);
105
float
d
;
106
if
(ieta>=0 && ieta<16)
d
= 0.196;
//cm
107
else
if
(ieta>=16 && ieta<32)
d
= 0.193;
//cm
108
else
if
(ieta>=32 && ieta<48)
d
= 0.2;
//cm
109
else
d
= 0.190;
//cm
110
111
if
(hv>(
nominal
-2.))
return
1.;
112
double
efield=std::fabs(hv)/(
d
*1
e
+3);
113
double
enominal=
nominal
/(
d
*1
e
+3);
114
double
scale
=
Respo
(efield, enominal,
T
);
115
116
return
scale
;
117
}
118
119
//===============================================================================
120
float
corr_HV_EMBPS::Respo
(
float
e
,
float
e_nominal,
float
tempe)
121
{
122
if
(
e
< -999.)
return
1.;
123
if
(
e
< 0.01)
return
0;
124
if
(
e
> e_nominal )
return
1;
125
float
resp = (
InvCharge
(e_nominal)*
vdrift
(
e
,tempe))/(
InvCharge
(
e
)*
vdrift
(e_nominal,tempe));
126
return
resp;
127
}
128
129
//===============================================================================
130
float
corr_HV_EMBPS::InvCharge
(
float
e
)
131
{
132
float
q
= 1.;
133
if
(
e
> 2.)
q
=(1.+0.36/
e
);
134
else
if
(
e
>1
e
-6)
q
=(1.04+0.25/
e
);
135
return
q
;
136
}
137
138
//===============================================================================
139
float
corr_HV_EMBPS::vdrift
(
float
e
,
float
tempe)
140
{
141
const
float
T
= tempe;
142
static
const
float
P[6] = {-0.01481,-0.0075,0.141,12.4,1.627,0.317};
143
if
(
e
< -999.)
return
0.;
144
float
vd = (P[0]*
T
+1)*( P[2]*
e
*
std::log
(1+ (P[3]/
e
)) + P[4]*
std::pow
(
e
,P[5])) + P[1]*
T
;
// vdrift formula walcowialk mm/micro_s
145
return
vd;
146
}
147
148
149
//===============================================================================
150
float
corr_HV_EMBPS::getDataCorrection
(
float
hv,
float
eta
)
151
//
152
// values are measured 1 + <Eraw0_end of run>-<Eraw0 beginning of run> / <Eraw0> (i.e <Eraw end>/<Eraw begin>)
153
// for the sectors with both sides going from 1600->1200V from beginning to end of run
154
// as a function of eta (possibly)
155
{
156
if
(hv>1300.)
return
1.;
157
if
(hv>1000.) {
158
if
(
eta
<-1.2)
return
0.9925;
159
if
(
eta
<-0.8)
return
0.9918;
160
if
(
eta
<-0.4)
return
0.9889;
161
if
(
eta
<0.)
return
0.9935;
162
if
(
eta
<0.4)
return
0.9908;
163
if
(
eta
<0.8)
return
0.99197;
164
if
(
eta
<1.2)
return
0.9974;
165
return
0.9971;
166
}
167
return
1.015;
168
}
169
170
//===============================================================================
171
float
corr_HV_EMBPS::getExtraScaleFactor
()
172
{
173
return
1.013;
174
}
corr_HV_EMBPS::getCorr
float getCorr(int run, float eta, float phi) const
get correction factor to apply to raw EMBPS energy : corrected raw EMBPS energy = correction factor *...
Definition:
corr_HV_EMBPS.cxx:45
phi
Scalar phi() const
phi method
Definition:
AmgMatrixBasePlugin.h:64
CaloCellPos2Ntuple.int
int
Definition:
CaloCellPos2Ntuple.py:24
corr_HV_EMBPS::~corr_HV_EMBPS
~corr_HV_EMBPS()
Definition:
corr_HV_EMBPS.cxx:39
eta
Scalar eta() const
pseudorapidity method
Definition:
AmgMatrixBasePlugin.h:79
hist_file_dump.d
d
Definition:
hist_file_dump.py:137
extractSporadic.c1
c1
Definition:
extractSporadic.py:134
conifer::pow
constexpr int pow(int x)
Definition:
conifer.h:20
corr_HV_EMBPS::m_file
TFile * m_file
Definition:
corr_HV_EMBPS.h:64
TProfile2D
Definition:
rootspy.cxx:531
M_PI
#define M_PI
Definition:
ActiveFraction.h:11
corr_HV_EMBPS::Respo
static float Respo(float e, float e_nominal, float tempe)
Definition:
corr_HV_EMBPS.cxx:120
corr_HV_EMBPS::getDataCorrection
static float getDataCorrection(float hv, float eta)
Definition:
corr_HV_EMBPS.cxx:150
yodamerge_tmp.scale
scale
Definition:
yodamerge_tmp.py:138
corr_HV_EMBPS::getRecoCorrection
static float getRecoCorrection(float hv, float eta)
Definition:
corr_HV_EMBPS.cxx:99
corr_HV_EMBPS::vdrift
static float vdrift(float e, float tempe)
Definition:
corr_HV_EMBPS.cxx:139
top::nominal
@ nominal
Definition:
ScaleFactorRetriever.h:29
run
Definition:
run.py:1
corr_HV_EMBPS.h
corr_HV_EMBPS::corr_HV_EMBPS
corr_HV_EMBPS()
constructor (initialization done there reading a root file for the HV maps per period
Definition:
corr_HV_EMBPS.cxx:19
PathResolver.h
name
std::string name
Definition:
Control/AthContainers/Root/debug.cxx:192
compileRPVLLRates.c2
c2
Definition:
compileRPVLLRates.py:361
TProfile2D::GetBinContent
double GetBinContent(int) const
Definition:
rootspy.cxx:546
PathResolverFindCalibFile
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
Definition:
PathResolver.cxx:431
corr_HV_EMBPS::m_hHV
TProfile2D * m_hHV[6][2]
Definition:
corr_HV_EMBPS.h:63
corr_HV_EMBPS::InvCharge
static float InvCharge(float e)
Definition:
corr_HV_EMBPS.cxx:130
DiTauMassTools::MaxHistStrategyV2::e
e
Definition:
PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:26
python.CaloCondTools.log
log
Definition:
CaloCondTools.py:20
CaloCellTimeCorrFiller.filename
filename
Definition:
CaloCellTimeCorrFiller.py:24
extractSporadic.q
list q
Definition:
extractSporadic.py:98
corr_HV_EMBPS::getExtraScaleFactor
static float getExtraScaleFactor()
Definition:
corr_HV_EMBPS.cxx:171
TSU::T
unsigned long long T
Definition:
L1TopoDataTypes.h:35
Generated on Sun Jun 30 2024 21:13:01 for ATLAS Offline Software by
1.8.18