CMS 3D CMS Logo
Functions
ZdcSimpleRecAlgoImpl Namespace Reference

Functions

template<class Digi , class RecHit >
RecHit reco1 (const Digi &digi, const HcalCoder &coder, const HcalCalibrations &calibs, const std::vector< unsigned int > &myNoiseTS, const std::vector< unsigned int > &mySignalTS, int lowGainOffset, double lowGainFrac, bool slewCorrect, const HcalPulseContainmentCorrection *corr, HcalTimeSlew::BiasSetting slewFlavor)
 
template<class Digi , class RecHit >
RecHit reco2 (const Digi &digi, const HcalCoder &coder, const HcalCalibrations &calibs, const std::vector< unsigned int > &myNoiseTS, const std::vector< unsigned int > &mySignalTS, int lowGainOffset, double lowGainFrac, bool slewCorrect, const HcalPulseContainmentCorrection *corr, HcalTimeSlew::BiasSetting slewFlavor)
 

Function Documentation

template<class Digi , class RecHit >
RecHit ZdcSimpleRecAlgoImpl::reco1 ( const Digi &  digi,
const HcalCoder &  coder,
const HcalCalibrations &  calibs,
const std::vector< unsigned int > &  myNoiseTS,
const std::vector< unsigned int > &  mySignalTS,
int  lowGainOffset,
double  lowGainFrac,
bool  slewCorrect,
const HcalPulseContainmentCorrection *  corr,
HcalTimeSlew::BiasSetting  slewFlavor 
)
inline

Definition at line 31 of file ZdcSimpleRecAlgo.cc.

References HcalCoder::adc2fC(), HcalPulseContainmentCorrection::getCorrection(), mps_fire::i, LogDebug, HcalHitReconstructor_zdc_cfi::lowGainFrac, HcalHitReconstructor_zdc_cfi::lowGainOffset, findQualityFiles::maxI, dqmiodumpmetadata::n, HcalCalibrations::pedestal(), HcalCalibrations::respcorrgain(), CaloSamples::size(), and ntuplemaker::time.

40  {
41  CaloSamples tool;
42  coder.adc2fC(digi, tool);
43  int ifirst = mySignalTS[0];
44  int n = mySignalTS.size();
45  double ampl = 0;
46  int maxI = -1;
47  double maxA = -1e10;
48  double ta = 0;
49  double fc_ampl = 0;
50  double lowGEnergy = 0;
51  double lowGfc_ampl = 0;
52  double TempLGAmp = 0;
53  // TS increment for regular energy to lowGainEnergy
54  // Signal in higher TS (effective "low Gain") has a fraction of the whole signal
55  // This constant for fC --> GeV is dervied from 2010 PbPb analysis of single neutrons
56  // assumed similar fraction for EM and HAD sections
57  // this variable converts from current assumed TestBeam values for fC--> GeV
58  // to the lowGain TS region fraction value (based on 1N Had, assume EM same response)
59  // regular energy
60  for (int i = ifirst; i < tool.size() && i < n + ifirst; i++) {
61  int capid = digi[i].capid();
62  ta = (tool[i] - calibs.pedestal(capid)); // pedestal subtraction
63  fc_ampl += ta;
64  ta *= calibs.respcorrgain(capid); // fC --> GeV
65  ampl += ta;
66  if (ta > maxA) {
67  maxA = ta;
68  maxI = i;
69  }
70  }
71  // calculate low Gain Energy (in 2010 PbPb, signal TS 4,5,6, lowGain TS: 6,7,8)
72  int topLowGain = 10;
73  if ((n + ifirst + lowGainOffset) <= 10) {
74  topLowGain = n + ifirst + lowGainOffset;
75  } else {
76  topLowGain = 10;
77  }
78  for (int iLG = (ifirst + lowGainOffset); iLG < tool.size() && iLG < topLowGain; iLG++) {
79  int capid = digi[iLG].capid();
80  TempLGAmp = (tool[iLG] - calibs.pedestal(capid)); // pedestal subtraction
81  lowGfc_ampl += TempLGAmp;
82  TempLGAmp *= calibs.respcorrgain(capid); // fC --> GeV
83  TempLGAmp *= lowGainFrac; // TS (signalRegion) --> TS (lowGainRegion)
84  lowGEnergy += TempLGAmp;
85  }
86  double time = -9999;
87  // Time based on regular energy (lowGainEnergy signal assumed to happen at same Time)
89  if (maxI == 0 || maxI == (tool.size() - 1)) {
90  LogDebug("HCAL Pulse") << "ZdcSimpleRecAlgo::reco1 :"
91  << " Invalid max amplitude position, "
92  << " max Amplitude: " << maxI << " first: " << ifirst << " last: " << (tool.size() - 1)
93  << std::endl;
94  } else {
95  int capid = digi[maxI - 1].capid();
96  double Energy0 = ((tool[maxI - 1]) * calibs.respcorrgain(capid));
97  // if any of the energies used in the weight are negative, make them 0 instead
98  // these are actually QIE values, not energy
99  if (Energy0 < 0) {
100  Energy0 = 0.;
101  }
102  capid = digi[maxI].capid();
103  double Energy1 = ((tool[maxI]) * calibs.respcorrgain(capid));
104  if (Energy1 < 0) {
105  Energy1 = 0.;
106  }
107  capid = digi[maxI + 1].capid();
108  double Energy2 = ((tool[maxI + 1]) * calibs.respcorrgain(capid));
109  if (Energy2 < 0) {
110  Energy2 = 0.;
111  }
112  //
113  double TSWeightEnergy = ((maxI - 1) * Energy0 + maxI * Energy1 + (maxI + 1) * Energy2);
114  double EnergySum = Energy0 + Energy1 + Energy2;
115  double AvgTSPos = 0.;
116  if (EnergySum != 0)
117  AvgTSPos = TSWeightEnergy / EnergySum;
118  // If time is zero, set it to the "nonsensical" -99
119  // Time should be between 75ns and 175ns (Timeslices 3-7)
120  if (AvgTSPos == 0) {
121  time = -99;
122  } else {
123  time = (AvgTSPos * 25.0);
124  }
125  if (corr != nullptr) {
126  // Apply phase-based amplitude correction:
127  ampl *= corr->getCorrection(fc_ampl);
128  }
129  }
130  return RecHit(digi.id(), ampl, time, lowGEnergy);
131  }
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:670
mps_fire.i
i
Definition: mps_fire.py:341
HcalPulseContainmentCorrection::getCorrection
double getCorrection(double fc_ampl) const
Definition: HcalPulseContainmentCorrection.cc:51
HcalHitReconstructor_zdc_cfi.lowGainOffset
lowGainOffset
Definition: HcalHitReconstructor_zdc_cfi.py:14
HcalHitReconstructor_zdc_cfi.lowGainFrac
lowGainFrac
Definition: HcalHitReconstructor_zdc_cfi.py:15
HcalCoder::adc2fC
virtual void adc2fC(const HBHEDataFrame &df, CaloSamples &lf) const =0
HcalCalibrations::pedestal
constexpr double pedestal(int fCapId) const
get pedestal for capid=0..3
Definition: HcalCalibrations.h:46
CaloSamples::size
int size() const
get the size
Definition: CaloSamples.h:24
HcalCalibrations::respcorrgain
constexpr double respcorrgain(int fCapId) const
get response corrected gain for capid=0..3
Definition: HcalCalibrations.h:42
template<class Digi , class RecHit >
RecHit ZdcSimpleRecAlgoImpl::reco2 ( const Digi &  digi,
const HcalCoder &  coder,
const HcalCalibrations &  calibs,
const std::vector< unsigned int > &  myNoiseTS,
const std::vector< unsigned int > &  mySignalTS,
int  lowGainOffset,
double  lowGainFrac,
bool  slewCorrect,
const HcalPulseContainmentCorrection *  corr,
HcalTimeSlew::BiasSetting  slewFlavor 
)
inline

Definition at line 136 of file ZdcSimpleRecAlgo.cc.

References HcalCoder::adc2fC(), HcalPulseContainmentCorrection::getCorrection(), LogDebug, HcalHitReconstructor_zdc_cfi::lowGainFrac, HcalHitReconstructor_zdc_cfi::lowGainOffset, findQualityFiles::maxI, hgcalDigitizer_cfi::noise, HcalCalibrations::respcorrgain(), CaloSamples::size(), and ntuplemaker::time.

145  {
146  CaloSamples tool;
147  coder.adc2fC(digi, tool);
148  // Reads noiseTS and signalTS from database
149  int ifirst = mySignalTS[0];
150  // int n = mySignalTS.size();
151  double ampl = 0;
152  int maxI = -1;
153  double maxA = -1e10;
154  double ta = 0;
155  double fc_ampl = 0;
156  double lowGEnergy = 0;
157  double lowGfc_ampl = 0;
158  double TempLGAmp = 0;
159  // TS increment for regular energy to lowGainEnergy
160  // Signal in higher TS (effective "low Gain") has a fraction of the whole signal
161  // This constant for fC --> GeV is dervied from 2010 PbPb analysis of single neutrons
162  // assumed similar fraction for EM and HAD sections
163  // this variable converts from current assumed TestBeam values for fC--> GeV
164  // to the lowGain TS region fraction value (based on 1N Had, assume EM same response)
165  double Allnoise = 0;
166  int noiseslices = 0;
167  int CurrentTS = 0;
168  double noise = 0;
169  // regular energy (both use same noise)
170  for (unsigned int iv = 0; iv < myNoiseTS.size(); ++iv) {
171  CurrentTS = myNoiseTS[iv];
172  if (CurrentTS >= digi.size())
173  continue;
174  Allnoise += tool[CurrentTS];
175  noiseslices++;
176  }
177  if (noiseslices != 0) {
178  noise = (Allnoise) / double(noiseslices);
179  } else {
180  noise = 0;
181  }
182  for (unsigned int ivs = 0; ivs < mySignalTS.size(); ++ivs) {
183  CurrentTS = mySignalTS[ivs];
184  if (CurrentTS >= digi.size())
185  continue;
186  int capid = digi[CurrentTS].capid();
187  // if(noise<0){
188  // // flag hit as having negative noise, and don't subtract anything, because
189  // // it will falsely increase the energy
190  // noisefactor=0.;
191  // }
192  ta = tool[CurrentTS] - noise;
193  fc_ampl += ta;
194  ta *= calibs.respcorrgain(capid); // fC --> GeV
195  ampl += ta;
196  if (ta > maxA) {
197  maxA = ta;
198  maxI = CurrentTS;
199  }
200  }
201  // calculate low Gain Energy (in 2010 PbPb, signal TS 4,5,6, lowGain TS: 6,7,8)
202  for (unsigned int iLGvs = 0; iLGvs < mySignalTS.size(); ++iLGvs) {
203  CurrentTS = mySignalTS[iLGvs] + lowGainOffset;
204  if (CurrentTS >= digi.size())
205  continue;
206  int capid = digi[CurrentTS].capid();
207  TempLGAmp = tool[CurrentTS] - noise;
208  lowGfc_ampl += TempLGAmp;
209  TempLGAmp *= calibs.respcorrgain(capid); // fC --> GeV
210  TempLGAmp *= lowGainFrac; // TS (signalRegion) --> TS (lowGainRegion)
211  lowGEnergy += TempLGAmp;
212  }
213  // if(ta<0){
214  // // flag hits that have negative energy
215  // }
216 
217  double time = -9999;
218  // Time based on regular energy (lowGainEnergy signal assumed to happen at same Time)
220  if (maxI == 0 || maxI == (tool.size() - 1)) {
221  LogDebug("HCAL Pulse") << "ZdcSimpleRecAlgo::reco2 :"
222  << " Invalid max amplitude position, "
223  << " max Amplitude: " << maxI << " first: " << ifirst << " last: " << (tool.size() - 1)
224  << std::endl;
225  } else {
226  int capid = digi[maxI - 1].capid();
227  double Energy0 = ((tool[maxI - 1]) * calibs.respcorrgain(capid));
228  // if any of the energies used in the weight are negative, make them 0 instead
229  // these are actually QIE values, not energy
230  if (Energy0 < 0) {
231  Energy0 = 0.;
232  }
233  capid = digi[maxI].capid();
234  double Energy1 = ((tool[maxI]) * calibs.respcorrgain(capid));
235  if (Energy1 < 0) {
236  Energy1 = 0.;
237  }
238  capid = digi[maxI + 1].capid();
239  double Energy2 = ((tool[maxI + 1]) * calibs.respcorrgain(capid));
240  if (Energy2 < 0) {
241  Energy2 = 0.;
242  }
243  //
244  double TSWeightEnergy = ((maxI - 1) * Energy0 + maxI * Energy1 + (maxI + 1) * Energy2);
245  double EnergySum = Energy0 + Energy1 + Energy2;
246  double AvgTSPos = 0.;
247  if (EnergySum != 0)
248  AvgTSPos = TSWeightEnergy / EnergySum;
249  // If time is zero, set it to the "nonsensical" -99
250  // Time should be between 75ns and 175ns (Timeslices 3-7)
251  if (AvgTSPos == 0) {
252  time = -99;
253  } else {
254  time = (AvgTSPos * 25.0);
255  }
256  if (corr != nullptr) {
257  // Apply phase-based amplitude correction:
258  ampl *= corr->getCorrection(fc_ampl);
259  }
260  }
261  return RecHit(digi.id(), ampl, time, lowGEnergy);
262  }
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:670
HcalPulseContainmentCorrection::getCorrection
double getCorrection(double fc_ampl) const
Definition: HcalPulseContainmentCorrection.cc:51
HcalHitReconstructor_zdc_cfi.lowGainOffset
lowGainOffset
Definition: HcalHitReconstructor_zdc_cfi.py:14
hgcalDigitizer_cfi.noise
noise
Definition: hgcalDigitizer_cfi.py:150
HcalHitReconstructor_zdc_cfi.lowGainFrac
lowGainFrac
Definition: HcalHitReconstructor_zdc_cfi.py:15
HcalCoder::adc2fC
virtual void adc2fC(const HBHEDataFrame &df, CaloSamples &lf) const =0
CaloSamples::size
int size() const
get the size
Definition: CaloSamples.h:24
HcalCalibrations::respcorrgain
constexpr double respcorrgain(int fCapId) const
get response corrected gain for capid=0..3
Definition: HcalCalibrations.h:42
Generated for CMSSW Reference Manual by doxygen 1.8.11