CMS 3D CMS Logo

HcalDeterministicFit.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <cmath>
3 #include <climits>
5 
6 constexpr float HcalDeterministicFit::invGpar[3];
7 constexpr float HcalDeterministicFit::negThresh[2];
8 constexpr int HcalDeterministicFit::HcalRegion[2];
9 constexpr float HcalDeterministicFit::rCorr[2];
10 constexpr float HcalDeterministicFit::rCorrSiPM[2];
11 
12 using namespace std;
13 
15 
17 
20  bool iApplyTimeSlew,
21  double respCorr) {
22  fTimeSlew_ = tsParam;
23  fTimeSlewBias_ = bias;
24  applyTimeSlew_ = iApplyTimeSlew;
25  frespCorr_ = respCorr;
26 }
27 
28 constexpr float HcalDeterministicFit::landauFrac[];
29 // Landau function integrated in 1 ns intervals
30 //Landau pulse shape from https://indico.cern.ch/event/345283/contribution/3/material/slides/0.pdf
31 //Landau turn on by default at left edge of time slice
32 // normalized to 1 on [0,10000]
33 void HcalDeterministicFit::getLandauFrac(float tStart, float tEnd, float &sum) const {
34  if (std::abs(tStart - tEnd - tsWidth) < 0.1f) {
35  sum = 0.f;
36  return;
37  }
38  sum = landauFrac[int(ceil(tStart + tsWidth))];
39  return;
40 }
41 
42 constexpr float HcalDeterministicFit::siPM205Frac[];
43 void HcalDeterministicFit::get205Frac(float tStart, float tEnd, float &sum) const {
44  if (std::abs(tStart - tEnd - tsWidth) < 0.1f) {
45  sum = 0.f;
46  return;
47  }
48  sum = siPM205Frac[int(ceil(tStart + tsWidth))];
49  return;
50 }
51 
52 constexpr float HcalDeterministicFit::siPM206Frac[];
53 void HcalDeterministicFit::get206Frac(float tStart, float tEnd, float &sum) const {
54  if (std::abs(tStart - tEnd - tsWidth) < 0.1f) {
55  sum = 0.f;
56  return;
57  }
58  sum = siPM206Frac[int(ceil(tStart + tsWidth))];
59  return;
60 }
61 
62 constexpr float HcalDeterministicFit::siPM207Frac[];
63 void HcalDeterministicFit::get207Frac(float tStart, float tEnd, float &sum) const {
64  if (std::abs(tStart - tEnd - tsWidth) < 0.1f) {
65  sum = 0.f;
66  return;
67  }
68  sum = siPM207Frac[int(ceil(tStart + tsWidth))];
69  return;
70 }
71 
72 void HcalDeterministicFit::getFrac(float tStart, float tEnd, float &sum, FType fType) const {
73  switch (fType) {
74  case shape205:
75  get205Frac(tStart, tEnd, sum);
76  break;
77  case shape206:
78  get206Frac(tStart, tEnd, sum);
79  break;
80  case shape207:
81  get207Frac(tStart, tEnd, sum);
82  break;
83  case shapeLandau:
84  getLandauFrac(tStart, tEnd, sum);
85  break;
86  }
87 }
88 
90  float &reconstructedEnergy,
91  float &reconstructedTime,
92  const HcalTimeSlew *hcalTimeSlew_delay) const {
93  unsigned int soi = channelData.soi();
94 
95  std::vector<double> corrCharge;
96  corrCharge.reserve(channelData.nSamples());
97  std::vector<double> inputCharge;
98  inputCharge.reserve(channelData.nSamples());
99  std::vector<double> inputPedestal;
100  inputPedestal.reserve(channelData.nSamples());
101  std::vector<double> inputNoise;
102  inputNoise.reserve(channelData.nSamples());
103 
104  double gainCorr = 0;
105  double respCorr = 0;
106 
107  for (unsigned int ip = 0; ip < channelData.nSamples(); ip++) {
108  double charge = channelData.tsRawCharge(ip);
109  double ped = channelData.tsPedestal(ip);
110  double noise = channelData.tsPedestalWidth(ip);
111  double gain = channelData.tsGain(ip);
112 
113  gainCorr = gain;
114  inputCharge.push_back(charge);
115  inputPedestal.push_back(ped);
116  inputNoise.push_back(noise);
117  }
118 
119  fPedestalSubFxn_.calculate(inputCharge, inputPedestal, inputNoise, corrCharge, soi, channelData.nSamples());
120 
121  if (fTimeSlew_ == 0)
122  respCorr = 1.0;
123  else if (fTimeSlew_ == 1)
124  channelData.hasTimeInfo() ? respCorr = rCorrSiPM[0] : respCorr = rCorr[0];
125  else if (fTimeSlew_ == 2)
126  channelData.hasTimeInfo() ? respCorr = rCorrSiPM[1] : respCorr = rCorr[1];
127  else if (fTimeSlew_ == 3)
128  respCorr = frespCorr_;
129 
130  float tsShift3, tsShift4, tsShift5;
131  tsShift3 = 0.f, tsShift4 = 0.f, tsShift5 = 0.f;
132 
133  if (applyTimeSlew_) {
134  tsShift3 = hcalTimeSlew_delay->delay(inputCharge[soi - 1], fTimeSlew_, fTimeSlewBias_, !channelData.hasTimeInfo());
135  tsShift4 = hcalTimeSlew_delay->delay(inputCharge[soi], fTimeSlew_, fTimeSlewBias_, !channelData.hasTimeInfo());
136  tsShift5 = hcalTimeSlew_delay->delay(inputCharge[soi + 1], fTimeSlew_, fTimeSlewBias_, !channelData.hasTimeInfo());
137  }
138 
139  float ch3, ch4, ch5, i3, n3, nn3, i4, n4, i5, n5;
140  ch4 = 0.f, i3 = 0.f, n3 = 0.f, nn3 = 0.f, i4 = 0.f, n4 = 0.f, i5 = 0.f, n5 = 0.f;
141 
142  FType fType;
143  if (channelData.hasTimeInfo() && channelData.recoShape() == 205)
144  fType = shape205;
145  else if (channelData.hasTimeInfo() && channelData.recoShape() == 206)
146  fType = shape206;
147  else if (channelData.hasTimeInfo() && channelData.recoShape() == 207)
148  fType = shape207;
149  else
150  fType = shapeLandau;
151 
152  getFrac(-tsShift3, -tsShift3 + tsWidth, i3, fType);
153  getFrac(-tsShift3 + tsWidth, -tsShift3 + tsWidth * 2, n3, fType);
154  getFrac(-tsShift3 + tsWidth * 2, -tsShift3 + tsWidth * 3, nn3, fType);
155 
156  getFrac(-tsShift4, -tsShift4 + tsWidth, i4, fType);
157  getFrac(-tsShift4 + tsWidth, -tsShift4 + tsWidth * 2, n4, fType);
158 
159  getFrac(-tsShift5, -tsShift5 + tsWidth, i5, fType);
160  getFrac(-tsShift5 + tsWidth, -tsShift5 + tsWidth * 2, n5, fType);
161 
162  if (i3 != 0 && i4 != 0 && i5 != 0) {
163  ch3 = corrCharge[soi - 1] / i3;
164  ch4 = (i3 * corrCharge[soi] - n3 * corrCharge[soi - 1]) / (i3 * i4);
165  ch5 = (n3 * n4 * corrCharge[soi - 1] - i4 * nn3 * corrCharge[soi - 1] - i3 * n4 * corrCharge[soi] +
166  i3 * i4 * corrCharge[soi + 1]) /
167  (i3 * i4 * i5);
168 
169  if (ch3 < negThresh[0]) {
170  ch3 = negThresh[0];
171  ch4 = corrCharge[soi] / i4;
172  ch5 = (i4 * corrCharge[soi + 1] - n4 * corrCharge[soi]) / (i4 * i5);
173  }
174  if (ch5 < negThresh[0] && ch4 > negThresh[1]) {
175  double ratio = (corrCharge[soi] - ch3 * i3) / (corrCharge[soi + 1] - negThresh[0] * i5);
176  if (ratio < 5 && ratio > 0.5) {
177  double invG = invGpar[0] + invGpar[1] * std::sqrt(2 * std::log(invGpar[2] / ratio));
178  float iG = 0.f;
179  getFrac(-invG, -invG + tsWidth, iG, fType);
180  if (iG != 0) {
181  ch4 = (corrCharge[soi] - ch3 * n3) / (iG);
182  tsShift4 = invG;
183  }
184  }
185  }
186  }
187 
188  if (ch4 < 1) {
189  ch4 = 0.f;
190  }
191 
192  reconstructedEnergy = ch4 * gainCorr * respCorr;
193  reconstructedTime = tsShift4;
194 }
HcalDeterministicFit.h
HcalDeterministicFit::siPM206Frac
static constexpr float siPM206Frac[125]
Definition: HcalDeterministicFit.h:82
testProducerWithPsetDescEmpty_cfi.i3
i3
Definition: testProducerWithPsetDescEmpty_cfi.py:47
HcalDeterministicFit::~HcalDeterministicFit
~HcalDeterministicFit()
Definition: HcalDeterministicFit.cc:16
HBHEChannelInfo::hasTimeInfo
constexpr bool hasTimeInfo() const
Definition: HBHEChannelInfo.h:139
HcalTimeSlew::ParaSource
ParaSource
Definition: HcalTimeSlew.h:70
HBHEChannelInfo::nSamples
constexpr unsigned nSamples() const
Definition: HBHEChannelInfo.h:136
HcalDeterministicFit::get206Frac
void get206Frac(float tStart, float tEnd, float &sum) const
Definition: HcalDeterministicFit.cc:53
HcalDeterministicFit::getFrac
void getFrac(float, float, float &, FType) const
Definition: HcalDeterministicFit.cc:72
f
double f[11][100]
Definition: MuScleFitUtils.cc:78
HcalDeterministicFit::get207Frac
void get207Frac(float tStart, float tEnd, float &sum) const
Definition: HcalDeterministicFit.cc:63
HBHEChannelInfo::tsPedestalWidth
constexpr double tsPedestalWidth(const unsigned ts) const
Definition: HBHEChannelInfo.h:167
HcalDeterministicFit::phase1Apply
void phase1Apply(const HBHEChannelInfo &channelData, float &reconstructedEnergy, float &reconstructedTime, const HcalTimeSlew *hcalTimeSlew_delay) const
Definition: HcalDeterministicFit.cc:89
HBHEChannelInfo::tsGain
constexpr double tsGain(const unsigned ts) const
Definition: HBHEChannelInfo.h:168
HcalDeterministicFit::landauFrac
static constexpr float landauFrac[125]
Definition: HcalDeterministicFit.h:51
HcalDeterministicFit::getLandauFrac
void getLandauFrac(float tStart, float tEnd, float &sum) const
Definition: HcalDeterministicFit.cc:33
HcalDeterministicFit::FType
FType
Definition: HcalDeterministicFit.h:21
HcalDeterministicFit::siPM207Frac
static constexpr float siPM207Frac[125]
Definition: HcalDeterministicFit.h:97
HcalTimeSlew::BiasSetting
BiasSetting
Definition: HcalTimeSlew.h:71
HcalTimeSlew
Definition: HcalTimeSlew.h:19
HcalTimeSlew::delay
float delay(float fC, BiasSetting bias=Medium) const
Returns the amount (ns) by which a pulse of the given number of fC will be delayed by the timeslew ef...
Definition: HcalTimeSlew.cc:20
reco::ceil
constexpr int32_t ceil(float num)
Definition: constexpr_cmath.h:7
HcalDeterministicFit::rCorr
static constexpr float rCorr[2]
Definition: HcalDeterministicFit.h:48
HcalDeterministicFit::get205Frac
void get205Frac(float tStart, float tEnd, float &sum) const
Definition: HcalDeterministicFit.cc:43
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
HBHEChannelInfo::recoShape
constexpr int recoShape() const
Definition: HBHEChannelInfo.h:134
HBHEChannelInfo::tsPedestal
constexpr double tsPedestal(const unsigned ts) const
Definition: HBHEChannelInfo.h:166
particleFlowDisplacedVertex_cfi.ratio
ratio
Definition: particleFlowDisplacedVertex_cfi.py:93
ALCARECOTkAlJpsiMuMu_cff.charge
charge
Definition: ALCARECOTkAlJpsiMuMu_cff.py:47
HcalDeterministicFit::siPM205Frac
static constexpr float siPM205Frac[125]
Definition: HcalDeterministicFit.h:67
HcalDeterministicFit::init
void init(HcalTimeSlew::ParaSource tsParam, HcalTimeSlew::BiasSetting bias, bool iApplyTimeSlew, double respCorr)
Definition: HcalDeterministicFit.cc:18
HcalDeterministicFit::rCorrSiPM
static constexpr float rCorrSiPM[2]
Definition: HcalDeterministicFit.h:49
createfilelist.int
int
Definition: createfilelist.py:10
gpuVertexFinder::noise
__shared__ int noise
Definition: gpuFitVertices.h:50
PedestalClient_cfi.gain
gain
Definition: PedestalClient_cfi.py:37
HcalDeterministicFit::HcalDeterministicFit
HcalDeterministicFit()
Definition: HcalDeterministicFit.cc:14
HcalDeterministicFit::invGpar
static constexpr float invGpar[3]
Definition: HcalDeterministicFit.h:47
std
Definition: JetResolutionObject.h:76
HBHEChannelInfo::tsRawCharge
constexpr double tsRawCharge(const unsigned ts) const
Definition: HBHEChannelInfo.h:165
HcalDeterministicFit::HcalRegion
static constexpr int HcalRegion[2]
Definition: HcalDeterministicFit.h:44
dqm-mbProfile.log
log
Definition: dqm-mbProfile.py:17
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
HBHEChannelInfo::soi
constexpr unsigned soi() const
Definition: HBHEChannelInfo.h:137
HcalDeterministicFit::negThresh
static constexpr float negThresh[2]
Definition: HcalDeterministicFit.h:46
HBHEChannelInfo
Definition: HBHEChannelInfo.h:15