CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
EcalUncalibRecHitWorkerMultiFit.cc
Go to the documentation of this file.
2 
8 
12 
15  // get the BX for the pulses to be activated
16  std::vector<int32_t> activeBXs = ps.getParameter<std::vector<int32_t>>("activeBXs");
17  activeBX.resize(activeBXs.size());
18  for (unsigned int ibx = 0; ibx < activeBXs.size(); ++ibx) {
19  activeBX.coeffRef(ibx) = activeBXs[ibx];
20  }
21 
22  // uncertainty calculation (CPU intensive)
23  ampErrorCalculation_ = ps.getParameter<bool>("ampErrorCalculation");
24  useLumiInfoRunHeader_ = ps.getParameter<bool>("useLumiInfoRunHeader");
25 
26  if (useLumiInfoRunHeader_) {
27  bunchSpacing_ = c.consumes<unsigned int>(edm::InputTag("bunchSpacingProducer"));
29  } else {
30  bunchSpacingManual_ = ps.getParameter<int>("bunchSpacing");
31  }
32 
33  doPrefitEB_ = ps.getParameter<bool>("doPrefitEB");
34  doPrefitEE_ = ps.getParameter<bool>("doPrefitEE");
35 
36  prefitMaxChiSqEB_ = ps.getParameter<double>("prefitMaxChiSqEB");
37  prefitMaxChiSqEE_ = ps.getParameter<double>("prefitMaxChiSqEE");
38 
39  dynamicPedestalsEB_ = ps.getParameter<bool>("dynamicPedestalsEB");
40  dynamicPedestalsEE_ = ps.getParameter<bool>("dynamicPedestalsEE");
41  mitigateBadSamplesEB_ = ps.getParameter<bool>("mitigateBadSamplesEB");
42  mitigateBadSamplesEE_ = ps.getParameter<bool>("mitigateBadSamplesEE");
43  gainSwitchUseMaxSampleEB_ = ps.getParameter<bool>("gainSwitchUseMaxSampleEB");
44  gainSwitchUseMaxSampleEE_ = ps.getParameter<bool>("gainSwitchUseMaxSampleEE");
45  selectiveBadSampleCriteriaEB_ = ps.getParameter<bool>("selectiveBadSampleCriteriaEB");
46  selectiveBadSampleCriteriaEE_ = ps.getParameter<bool>("selectiveBadSampleCriteriaEE");
47  addPedestalUncertaintyEB_ = ps.getParameter<double>("addPedestalUncertaintyEB");
48  addPedestalUncertaintyEE_ = ps.getParameter<double>("addPedestalUncertaintyEE");
49  simplifiedNoiseModelForGainSwitch_ = ps.getParameter<bool>("simplifiedNoiseModelForGainSwitch");
61 
62  // algorithm to be used for timing
63  auto const& timeAlgoName = ps.getParameter<std::string>("timealgo");
64  if (timeAlgoName == "RatioMethod")
66  else if (timeAlgoName == "WeightsMethod")
68  else if (timeAlgoName == "crossCorrelationMethod") {
70  double startTime = ps.getParameter<double>("crossCorrelationStartTime");
71  double stopTime = ps.getParameter<double>("crossCorrelationStopTime");
72  double targetTimePrecision = ps.getParameter<double>("crossCorrelationTargetTimePrecision");
73  computeCC_ = std::make_unique<EcalUncalibRecHitTimingCCAlgo>(startTime, stopTime, targetTimePrecision);
74  } else if (timeAlgoName != "None")
75  edm::LogError("EcalUncalibRecHitError") << "No time estimation algorithm defined";
76 
77  // ratio method parameters
78  EBtimeFitParameters_ = ps.getParameter<std::vector<double>>("EBtimeFitParameters");
79  EEtimeFitParameters_ = ps.getParameter<std::vector<double>>("EEtimeFitParameters");
80  EBamplitudeFitParameters_ = ps.getParameter<std::vector<double>>("EBamplitudeFitParameters");
81  EEamplitudeFitParameters_ = ps.getParameter<std::vector<double>>("EEamplitudeFitParameters");
82  EBtimeFitLimits_.first = ps.getParameter<double>("EBtimeFitLimits_Lower");
83  EBtimeFitLimits_.second = ps.getParameter<double>("EBtimeFitLimits_Upper");
84  EEtimeFitLimits_.first = ps.getParameter<double>("EEtimeFitLimits_Lower");
85  EEtimeFitLimits_.second = ps.getParameter<double>("EEtimeFitLimits_Upper");
86  EBtimeConstantTerm_ = ps.getParameter<double>("EBtimeConstantTerm");
87  EEtimeConstantTerm_ = ps.getParameter<double>("EEtimeConstantTerm");
88  EBtimeNconst_ = ps.getParameter<double>("EBtimeNconst");
89  EEtimeNconst_ = ps.getParameter<double>("EEtimeNconst");
90  outOfTimeThreshG12pEB_ = ps.getParameter<double>("outOfTimeThresholdGain12pEB");
91  outOfTimeThreshG12mEB_ = ps.getParameter<double>("outOfTimeThresholdGain12mEB");
92  outOfTimeThreshG61pEB_ = ps.getParameter<double>("outOfTimeThresholdGain61pEB");
93  outOfTimeThreshG61mEB_ = ps.getParameter<double>("outOfTimeThresholdGain61mEB");
94  outOfTimeThreshG12pEE_ = ps.getParameter<double>("outOfTimeThresholdGain12pEE");
95  outOfTimeThreshG12mEE_ = ps.getParameter<double>("outOfTimeThresholdGain12mEE");
96  outOfTimeThreshG61pEE_ = ps.getParameter<double>("outOfTimeThresholdGain61pEE");
97  outOfTimeThreshG61mEE_ = ps.getParameter<double>("outOfTimeThresholdGain61mEE");
98  amplitudeThreshEB_ = ps.getParameter<double>("amplitudeThresholdEB");
99  amplitudeThreshEE_ = ps.getParameter<double>("amplitudeThresholdEE");
100 
101  // spike threshold
102  ebSpikeThresh_ = ps.getParameter<double>("ebSpikeThreshold");
103 
104  ebPulseShape_ = ps.getParameter<std::vector<double>>("ebPulseShape");
105  eePulseShape_ = ps.getParameter<std::vector<double>>("eePulseShape");
106 
107  // chi2 parameters for flags determination
108  kPoorRecoFlagEB_ = ps.getParameter<bool>("kPoorRecoFlagEB");
109  kPoorRecoFlagEE_ = ps.getParameter<bool>("kPoorRecoFlagEE");
110  chi2ThreshEB_ = ps.getParameter<double>("chi2ThreshEB_");
111  chi2ThreshEE_ = ps.getParameter<double>("chi2ThreshEE_");
112 }
113 
115  // common setup
117  peds = es.getHandle(pedsToken_);
118 
119  // for the multifit method
125 
126  // weights parameters for the time
127  grps = es.getHandle(grpsToken_);
128  wgts = es.getHandle(wgtsToken_);
129 
130  // which of the samples need be used
132 
133  // for the ratio method
136 
137  // for the time correction methods
139 
140  int nnoise = SampleVector::RowsAtCompileTime;
141  SampleMatrix& noisecorEBg12 = noisecors_[1][0];
142  SampleMatrix& noisecorEBg6 = noisecors_[1][1];
143  SampleMatrix& noisecorEBg1 = noisecors_[1][2];
144  SampleMatrix& noisecorEEg12 = noisecors_[0][0];
145  SampleMatrix& noisecorEEg6 = noisecors_[0][1];
146  SampleMatrix& noisecorEEg1 = noisecors_[0][2];
147 
148  for (int i = 0; i < nnoise; ++i) {
149  for (int j = 0; j < nnoise; ++j) {
150  int vidx = std::abs(j - i);
151  noisecorEBg12(i, j) = noisecovariances->EBG12SamplesCorrelation[vidx];
152  noisecorEEg12(i, j) = noisecovariances->EEG12SamplesCorrelation[vidx];
153  noisecorEBg6(i, j) = noisecovariances->EBG6SamplesCorrelation[vidx];
154  noisecorEEg6(i, j) = noisecovariances->EEG6SamplesCorrelation[vidx];
155  noisecorEBg1(i, j) = noisecovariances->EBG1SamplesCorrelation[vidx];
156  noisecorEEg1(i, j) = noisecovariances->EEG1SamplesCorrelation[vidx];
157  }
158  }
159 }
160 
162  unsigned int bunchspacing = 450;
163 
164  if (useLumiInfoRunHeader_) {
165  edm::Handle<unsigned int> bunchSpacingH;
166  evt.getByToken(bunchSpacing_, bunchSpacingH);
167  bunchspacing = *bunchSpacingH;
168  } else {
169  bunchspacing = bunchSpacingManual_;
170  }
171 
173  if (bunchspacing == 25) {
174  activeBX.resize(10);
175  activeBX << -5, -4, -3, -2, -1, 0, 1, 2, 3, 4;
176  } else {
177  //50ns configuration otherwise (also for no pileup)
178  activeBX.resize(5);
179  activeBX << -4, -2, 0, 2, 4;
180  }
181  }
182 }
183 
194  const std::vector<float>& amplitudeBins,
195  const std::vector<float>& shiftBins) {
196  // computed initially in ns. Than turned in the BX's, as
197  // EcalUncalibratedRecHit need be.
198  double theCorrection = 0;
199 
200  // sanity check for arrays
201  if (amplitudeBins.empty()) {
202  edm::LogError("EcalRecHitError") << "timeCorrAmplitudeBins is empty, forcing no time bias corrections.";
203 
204  return 0;
205  }
206 
207  if (amplitudeBins.size() != shiftBins.size()) {
208  edm::LogError("EcalRecHitError") << "Size of timeCorrAmplitudeBins different from "
209  "timeCorrShiftBins. Forcing no time bias corrections. ";
210 
211  return 0;
212  }
213 
214  // FIXME? what about a binary search?
215  int myBin = -1;
216  for (int bin = 0; bin < (int)amplitudeBins.size(); bin++) {
217  if (ampli > amplitudeBins[bin]) {
218  myBin = bin;
219  } else {
220  break;
221  }
222  }
223 
224  if (myBin == -1) {
225  theCorrection = shiftBins[0];
226  } else if (myBin == ((int)(amplitudeBins.size() - 1))) {
227  theCorrection = shiftBins[myBin];
228  } else {
229  // interpolate linearly between two assingned points
230  theCorrection = (shiftBins[myBin + 1] - shiftBins[myBin]);
231  theCorrection *= (((double)ampli) - amplitudeBins[myBin]) / (amplitudeBins[myBin + 1] - amplitudeBins[myBin]);
232  theCorrection += shiftBins[myBin];
233  }
234 
235  // convert ns into clocks
236  constexpr double inv25 = 1. / 25.;
237  return theCorrection * inv25;
238 }
239 
241  const EcalDigiCollection& digis,
243  if (digis.empty())
244  return;
245 
246  // assume all digis come from the same subdetector (either barrel or endcap)
247  DetId detid(digis.begin()->id());
248  bool barrel = (detid.subdetId() == EcalBarrel);
249 
251  if (barrel) {
259  } else {
267  }
268 
269  FullSampleVector fullpulse(FullSampleVector::Zero());
270  FullSampleMatrix fullpulsecov(FullSampleMatrix::Zero());
271 
272  result.reserve(result.size() + digis.size());
273  for (auto itdg = digis.begin(); itdg != digis.end(); ++itdg) {
274  DetId detid(itdg->id());
275 
276  const EcalSampleMask* sampleMask_ = sampleMaskHand_.product();
277 
278  // intelligence for recHit computation
279  float offsetTime = 0;
280 
281  const EcalPedestals::Item* aped = nullptr;
282  const EcalMGPAGainRatio* aGain = nullptr;
283  const EcalXtalGroupId* gid = nullptr;
284  const EcalPulseShapes::Item* aPulse = nullptr;
285  const EcalPulseCovariances::Item* aPulseCov = nullptr;
286 
287  if (barrel) {
288  unsigned int hashedIndex = EBDetId(detid).hashedIndex();
289  aped = &peds->barrel(hashedIndex);
290  aGain = &gains->barrel(hashedIndex);
291  gid = &grps->barrel(hashedIndex);
292  aPulse = &pulseshapes->barrel(hashedIndex);
293  aPulseCov = &pulsecovariances->barrel(hashedIndex);
294  offsetTime = offtime->getEBValue();
295  } else {
296  unsigned int hashedIndex = EEDetId(detid).hashedIndex();
297  aped = &peds->endcap(hashedIndex);
298  aGain = &gains->endcap(hashedIndex);
299  gid = &grps->endcap(hashedIndex);
300  aPulse = &pulseshapes->endcap(hashedIndex);
301  aPulseCov = &pulsecovariances->endcap(hashedIndex);
302  offsetTime = offtime->getEEValue();
303  }
304 
305  double pedVec[3] = {aped->mean_x12, aped->mean_x6, aped->mean_x1};
306  double pedRMSVec[3] = {aped->rms_x12, aped->rms_x6, aped->rms_x1};
307  double gainRatios[3] = {1., aGain->gain12Over6(), aGain->gain6Over1() * aGain->gain12Over6()};
308 
309  for (int i = 0; i < EcalPulseShape::TEMPLATESAMPLES; ++i)
310  fullpulse(i + 7) = aPulse->pdfval[i];
311 
312  for (int i = 0; i < EcalPulseShape::TEMPLATESAMPLES; i++)
313  for (int j = 0; j < EcalPulseShape::TEMPLATESAMPLES; j++)
314  fullpulsecov(i + 7, j + 7) = aPulseCov->covval[i][j];
315 
316  // compute the right bin of the pulse shape using time calibration constants
318  EcalTimeCalibConstant itimeconst = 0;
319  if (it != itime->end()) {
320  itimeconst = (*it);
321  } else {
322  edm::LogError("EcalRecHitError") << "No time intercalib const found for xtal " << detid.rawId()
323  << "! something wrong with EcalTimeCalibConstants in your DB? ";
324  }
325 
326  int lastSampleBeforeSaturation = -2;
327  for (unsigned int iSample = 0; iSample < EcalDataFrame::MAXSAMPLES; iSample++) {
328  if (((EcalDataFrame)(*itdg)).sample(iSample).gainId() == 0) {
329  lastSampleBeforeSaturation = iSample - 1;
330  break;
331  }
332  }
333 
334  // === amplitude computation ===
335 
336  if (lastSampleBeforeSaturation == 4) { // saturation on the expected max sample
337  result.emplace_back((*itdg).id(), 4095 * 12, 0, 0, 0);
338  auto& uncalibRecHit = result.back();
339  uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kSaturated);
340  // do not propagate the default chi2 = -1 value to the calib rechit (mapped to 64), set it to 0 when saturation
341  uncalibRecHit.setChi2(0);
342  } else if (lastSampleBeforeSaturation >=
343  -1) { // saturation on other samples: cannot extrapolate from the fourth one
344  int gainId = ((EcalDataFrame)(*itdg)).sample(5).gainId();
345  if (gainId == 0)
346  gainId = 3;
347  auto pedestal = pedVec[gainId - 1];
348  auto gainratio = gainRatios[gainId - 1];
349  double amplitude = ((double)(((EcalDataFrame)(*itdg)).sample(5).adc()) - pedestal) * gainratio;
350  result.emplace_back((*itdg).id(), amplitude, 0, 0, 0);
351  auto& uncalibRecHit = result.back();
352  uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kSaturated);
353  // do not propagate the default chi2 = -1 value to the calib rechit (mapped to 64), set it to 0 when saturation
354  uncalibRecHit.setChi2(0);
355  } else {
356  // multifit
357  const SampleMatrixGainArray& noisecors = noisecor(barrel);
358 
359  result.push_back(multiFitMethod_.makeRecHit(*itdg, aped, aGain, noisecors, fullpulse, fullpulsecov, activeBX));
360  auto& uncalibRecHit = result.back();
361 
362  // === time computation ===
363  if (timealgo_ == ratioMethod) {
364  // ratio method
365  constexpr float clockToNsConstant = 25.;
366  constexpr float invClockToNs = 1. / clockToNsConstant;
367  if (not barrel) {
368  ratioMethod_endcap_.init(*itdg, *sampleMask_, pedVec, pedRMSVec, gainRatios);
373  double theTimeCorrectionEE = timeCorrection(
374  uncalibRecHit.amplitude(), timeCorrBias_->EETimeCorrAmplitudeBins, timeCorrBias_->EETimeCorrShiftBins);
375 
376  uncalibRecHit.setJitter(crh.timeMax - 5 + theTimeCorrectionEE);
377  uncalibRecHit.setJitterError(
378  std::sqrt(std::pow(crh.timeError, 2) + std::pow(EEtimeConstantTerm_ * invClockToNs, 2)));
379 
380  // consider flagging as kOutOfTime only if above noise
381  if (uncalibRecHit.amplitude() > pedRMSVec[0] * amplitudeThreshEE_) {
382  float outOfTimeThreshP = outOfTimeThreshG12pEE_;
383  float outOfTimeThreshM = outOfTimeThreshG12mEE_;
384  // determine if gain has switched away from gainId==1 (x12 gain)
385  // and determine cuts (number of 'sigmas') to ose for kOutOfTime
386  // >3k ADC is necessasry condition for gain switch to occur
387  if (uncalibRecHit.amplitude() > 3000.) {
388  for (int iSample = 0; iSample < EEDataFrame::MAXSAMPLES; iSample++) {
389  int GainId = ((EcalDataFrame)(*itdg)).sample(iSample).gainId();
390  if (GainId != 1) {
391  outOfTimeThreshP = outOfTimeThreshG61pEE_;
392  outOfTimeThreshM = outOfTimeThreshG61mEE_;
393  break;
394  }
395  }
396  }
397  float correctedTime = (crh.timeMax - 5) * clockToNsConstant + itimeconst + offsetTime;
398  float cterm = EEtimeConstantTerm_;
399  float sigmaped = pedRMSVec[0]; // approx for lower gains
400  float nterm = EEtimeNconst_ * sigmaped / uncalibRecHit.amplitude();
401  float sigmat = std::sqrt(nterm * nterm + cterm * cterm);
402  if ((correctedTime > sigmat * outOfTimeThreshP) || (correctedTime < -sigmat * outOfTimeThreshM)) {
403  uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kOutOfTime);
404  }
405  }
406 
407  } else {
408  ratioMethod_barrel_.init(*itdg, *sampleMask_, pedVec, pedRMSVec, gainRatios);
414 
415  double theTimeCorrectionEB = timeCorrection(
416  uncalibRecHit.amplitude(), timeCorrBias_->EBTimeCorrAmplitudeBins, timeCorrBias_->EBTimeCorrShiftBins);
417 
418  uncalibRecHit.setJitter(crh.timeMax - 5 + theTimeCorrectionEB);
419  uncalibRecHit.setJitterError(std::hypot(crh.timeError, EBtimeConstantTerm_ / clockToNsConstant));
420 
421  // consider flagging as kOutOfTime only if above noise
422  if (uncalibRecHit.amplitude() > pedRMSVec[0] * amplitudeThreshEB_) {
423  float outOfTimeThreshP = outOfTimeThreshG12pEB_;
424  float outOfTimeThreshM = outOfTimeThreshG12mEB_;
425  // determine if gain has switched away from gainId==1 (x12 gain)
426  // and determine cuts (number of 'sigmas') to ose for kOutOfTime
427  // >3k ADC is necessasry condition for gain switch to occur
428  if (uncalibRecHit.amplitude() > 3000.) {
429  for (int iSample = 0; iSample < EBDataFrame::MAXSAMPLES; iSample++) {
430  int GainId = ((EcalDataFrame)(*itdg)).sample(iSample).gainId();
431  if (GainId != 1) {
432  outOfTimeThreshP = outOfTimeThreshG61pEB_;
433  outOfTimeThreshM = outOfTimeThreshG61mEB_;
434  break;
435  }
436  }
437  }
438  float correctedTime = (crh.timeMax - 5) * clockToNsConstant + itimeconst + offsetTime;
439  float cterm = EBtimeConstantTerm_;
440  float sigmaped = pedRMSVec[0]; // approx for lower gains
441  float nterm = EBtimeNconst_ * sigmaped / uncalibRecHit.amplitude();
442  float sigmat = std::sqrt(nterm * nterm + cterm * cterm);
443  if ((correctedTime > sigmat * outOfTimeThreshP) || (correctedTime < -sigmat * outOfTimeThreshM)) {
444  uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kOutOfTime);
445  }
446  }
447  }
448  } else if (timealgo_ == weightsMethod) {
449  // weights method on the PU subtracted pulse shape
450  std::vector<double> amplitudes;
451  for (unsigned int ibx = 0; ibx < activeBX.size(); ++ibx)
452  amplitudes.push_back(uncalibRecHit.outOfTimeAmplitude(ibx));
453 
454  EcalTBWeights::EcalTDCId tdcid(1);
455  EcalTBWeights::EcalTBWeightMap const& wgtsMap = wgts->getMap();
456  EcalTBWeights::EcalTBWeightMap::const_iterator wit;
457  wit = wgtsMap.find(std::make_pair(*gid, tdcid));
458  if (wit == wgtsMap.end()) {
459  edm::LogError("EcalUncalibRecHitError")
460  << "No weights found for EcalGroupId: " << gid->id() << " and EcalTDCId: " << tdcid
461  << "\n skipping digi with id: " << detid.rawId();
462  result.pop_back();
463  continue;
464  }
465  const EcalWeightSet& wset = wit->second; // this is the EcalWeightSet
466 
469 
470  weights[0] = &mat1;
471  weights[1] = &mat2;
472 
473  double timerh;
474  if (detid.subdetId() == EcalEndcap) {
475  timerh = weightsMethod_endcap_.time(*itdg, amplitudes, aped, aGain, fullpulse, weights);
476  } else {
477  timerh = weightsMethod_barrel_.time(*itdg, amplitudes, aped, aGain, fullpulse, weights);
478  }
479  uncalibRecHit.setJitter(timerh);
480  uncalibRecHit.setJitterError(0.); // not computed with weights
481 
482  } else if (timealgo_ == crossCorrelationMethod) {
483  std::vector<double> amplitudes(activeBX.size());
484  for (unsigned int ibx = 0; ibx < activeBX.size(); ++ibx)
485  amplitudes[ibx] = uncalibRecHit.outOfTimeAmplitude(ibx);
486 
487  float jitterError = 0.;
488  float jitter = computeCC_->computeTimeCC(*itdg, amplitudes, aped, aGain, fullpulse, uncalibRecHit, jitterError);
489 
490  uncalibRecHit.setJitter(jitter);
491  uncalibRecHit.setJitterError(jitterError);
492 
493  } else { // no time method;
494  uncalibRecHit.setJitter(0.);
495  uncalibRecHit.setJitterError(0.);
496  }
497  }
498 
499  // set flags if gain switch has occurred
500  auto& uncalibRecHit = result.back();
501  if (((EcalDataFrame)(*itdg)).hasSwitchToGain6())
502  uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kHasSwitchToGain6);
503  if (((EcalDataFrame)(*itdg)).hasSwitchToGain1())
504  uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kHasSwitchToGain1);
505  }
506 }
507 
510  psd0.addNode((edm::ParameterDescription<std::vector<double>>("EBPulseShapeTemplate",
511  {1.13979e-02,
512  7.58151e-01,
513  1.00000e+00,
514  8.87744e-01,
515  6.73548e-01,
516  4.74332e-01,
517  3.19561e-01,
518  2.15144e-01,
519  1.47464e-01,
520  1.01087e-01,
521  6.93181e-02,
522  4.75044e-02},
523  true) and
524  edm::ParameterDescription<std::vector<double>>("EEPulseShapeTemplate",
525  {1.16442e-01,
526  7.56246e-01,
527  1.00000e+00,
528  8.97182e-01,
529  6.86831e-01,
530  4.91506e-01,
531  3.44111e-01,
532  2.45731e-01,
533  1.74115e-01,
534  1.23361e-01,
535  8.74288e-02,
536  6.19570e-02},
537  true)));
538 
539  psd0.addNode((edm::ParameterDescription<std::string>("EEdigiCollection", "", true) and
540  edm::ParameterDescription<std::string>("EBdigiCollection", "", true) and
541  edm::ParameterDescription<std::string>("ESdigiCollection", "", true) and
542  edm::ParameterDescription<bool>("UseLCcorrection", true, false) and
543  edm::ParameterDescription<std::vector<double>>(
544  "EBCorrNoiseMatrixG12",
545  {1.00000, 0.71073, 0.55721, 0.46089, 0.40449, 0.35931, 0.33924, 0.32439, 0.31581, 0.30481},
546  true) and
547  edm::ParameterDescription<std::vector<double>>(
548  "EECorrNoiseMatrixG12",
549  {1.00000, 0.71373, 0.44825, 0.30152, 0.21609, 0.14786, 0.11772, 0.10165, 0.09465, 0.08098},
550  true) and
551  edm::ParameterDescription<std::vector<double>>(
552  "EBCorrNoiseMatrixG06",
553  {1.00000, 0.70946, 0.58021, 0.49846, 0.45006, 0.41366, 0.39699, 0.38478, 0.37847, 0.37055},
554  true) and
555  edm::ParameterDescription<std::vector<double>>(
556  "EECorrNoiseMatrixG06",
557  {1.00000, 0.71217, 0.47464, 0.34056, 0.26282, 0.20287, 0.17734, 0.16256, 0.15618, 0.14443},
558  true) and
559  edm::ParameterDescription<std::vector<double>>(
560  "EBCorrNoiseMatrixG01",
561  {1.00000, 0.73354, 0.64442, 0.58851, 0.55425, 0.53082, 0.51916, 0.51097, 0.50732, 0.50409},
562  true) and
563  edm::ParameterDescription<std::vector<double>>(
564  "EECorrNoiseMatrixG01",
565  {1.00000, 0.72698, 0.62048, 0.55691, 0.51848, 0.49147, 0.47813, 0.47007, 0.46621, 0.46265},
566  true) and
567  edm::ParameterDescription<bool>("EcalPreMixStage1", false, true) and
568  edm::ParameterDescription<bool>("EcalPreMixStage2", false, true)));
569 
570  psd0.addOptionalNode(
571  (edm::ParameterDescription<std::vector<double>>(
572  "EBPulseShapeCovariance",
573  {3.001e-06, 1.233e-05, 0.000e+00, -4.416e-06, -4.571e-06, -3.614e-06, -2.636e-06, -1.286e-06, -8.410e-07,
574  -5.296e-07, 0.000e+00, 0.000e+00, 1.233e-05, 6.154e-05, 0.000e+00, -2.200e-05, -2.309e-05, -1.838e-05,
575  -1.373e-05, -7.334e-06, -5.088e-06, -3.745e-06, -2.428e-06, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,
576  0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,
577  -4.416e-06, -2.200e-05, 0.000e+00, 8.319e-06, 8.545e-06, 6.792e-06, 5.059e-06, 2.678e-06, 1.816e-06,
578  1.223e-06, 8.245e-07, 5.589e-07, -4.571e-06, -2.309e-05, 0.000e+00, 8.545e-06, 9.182e-06, 7.219e-06,
579  5.388e-06, 2.853e-06, 1.944e-06, 1.324e-06, 9.083e-07, 6.335e-07, -3.614e-06, -1.838e-05, 0.000e+00,
580  6.792e-06, 7.219e-06, 6.016e-06, 4.437e-06, 2.385e-06, 1.636e-06, 1.118e-06, 7.754e-07, 5.556e-07,
581  -2.636e-06, -1.373e-05, 0.000e+00, 5.059e-06, 5.388e-06, 4.437e-06, 3.602e-06, 1.917e-06, 1.322e-06,
582  9.079e-07, 6.529e-07, 4.752e-07, -1.286e-06, -7.334e-06, 0.000e+00, 2.678e-06, 2.853e-06, 2.385e-06,
583  1.917e-06, 1.375e-06, 9.100e-07, 6.455e-07, 4.693e-07, 3.657e-07, -8.410e-07, -5.088e-06, 0.000e+00,
584  1.816e-06, 1.944e-06, 1.636e-06, 1.322e-06, 9.100e-07, 9.115e-07, 6.062e-07, 4.436e-07, 3.422e-07,
585  -5.296e-07, -3.745e-06, 0.000e+00, 1.223e-06, 1.324e-06, 1.118e-06, 9.079e-07, 6.455e-07, 6.062e-07,
586  7.217e-07, 4.862e-07, 3.768e-07, 0.000e+00, -2.428e-06, 0.000e+00, 8.245e-07, 9.083e-07, 7.754e-07,
587  6.529e-07, 4.693e-07, 4.436e-07, 4.862e-07, 6.509e-07, 4.418e-07, 0.000e+00, 0.000e+00, 0.000e+00,
588  5.589e-07, 6.335e-07, 5.556e-07, 4.752e-07, 3.657e-07, 3.422e-07, 3.768e-07, 4.418e-07, 6.142e-07},
589  true) and
590  edm::ParameterDescription<std::vector<double>>(
591  "EEPulseShapeCovariance",
592  {3.941e-05, 3.333e-05, 0.000e+00, -1.449e-05, -1.661e-05, -1.424e-05, -1.183e-05, -6.842e-06, -4.915e-06,
593  -3.411e-06, 0.000e+00, 0.000e+00, 3.333e-05, 2.862e-05, 0.000e+00, -1.244e-05, -1.431e-05, -1.233e-05,
594  -1.032e-05, -5.883e-06, -4.154e-06, -2.902e-06, -2.128e-06, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,
595  0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00, 0.000e+00,
596  -1.449e-05, -1.244e-05, 0.000e+00, 5.840e-06, 6.649e-06, 5.720e-06, 4.812e-06, 2.708e-06, 1.869e-06,
597  1.330e-06, 9.186e-07, 6.446e-07, -1.661e-05, -1.431e-05, 0.000e+00, 6.649e-06, 7.966e-06, 6.898e-06,
598  5.794e-06, 3.157e-06, 2.184e-06, 1.567e-06, 1.084e-06, 7.575e-07, -1.424e-05, -1.233e-05, 0.000e+00,
599  5.720e-06, 6.898e-06, 6.341e-06, 5.347e-06, 2.859e-06, 1.991e-06, 1.431e-06, 9.839e-07, 6.886e-07,
600  -1.183e-05, -1.032e-05, 0.000e+00, 4.812e-06, 5.794e-06, 5.347e-06, 4.854e-06, 2.628e-06, 1.809e-06,
601  1.289e-06, 9.020e-07, 6.146e-07, -6.842e-06, -5.883e-06, 0.000e+00, 2.708e-06, 3.157e-06, 2.859e-06,
602  2.628e-06, 1.863e-06, 1.296e-06, 8.882e-07, 6.108e-07, 4.283e-07, -4.915e-06, -4.154e-06, 0.000e+00,
603  1.869e-06, 2.184e-06, 1.991e-06, 1.809e-06, 1.296e-06, 1.217e-06, 8.669e-07, 5.751e-07, 3.882e-07,
604  -3.411e-06, -2.902e-06, 0.000e+00, 1.330e-06, 1.567e-06, 1.431e-06, 1.289e-06, 8.882e-07, 8.669e-07,
605  9.522e-07, 6.717e-07, 4.293e-07, 0.000e+00, -2.128e-06, 0.000e+00, 9.186e-07, 1.084e-06, 9.839e-07,
606  9.020e-07, 6.108e-07, 5.751e-07, 6.717e-07, 7.911e-07, 5.493e-07, 0.000e+00, 0.000e+00, 0.000e+00,
607  6.446e-07, 7.575e-07, 6.886e-07, 6.146e-07, 4.283e-07, 3.882e-07, 4.293e-07, 5.493e-07, 7.027e-07},
608  true)),
609  true);
610 
612  psd.addNode(
613  edm::ParameterDescription<std::vector<int>>("activeBXs", {-5, -4, -3, -2, -1, 0, 1, 2, 3, 4}, true) and
614  edm::ParameterDescription<bool>("ampErrorCalculation", true, true) and
615  edm::ParameterDescription<bool>("useLumiInfoRunHeader", true, true) and
616  edm::ParameterDescription<int>("bunchSpacing", 0, true) and
617  edm::ParameterDescription<bool>("doPrefitEB", false, true) and
618  edm::ParameterDescription<bool>("doPrefitEE", false, true) and
619  edm::ParameterDescription<double>("prefitMaxChiSqEB", 25., true) and
620  edm::ParameterDescription<double>("prefitMaxChiSqEE", 10., true) and
621  edm::ParameterDescription<bool>("dynamicPedestalsEB", false, true) and
622  edm::ParameterDescription<bool>("dynamicPedestalsEE", false, true) and
623  edm::ParameterDescription<bool>("mitigateBadSamplesEB", false, true) and
624  edm::ParameterDescription<bool>("mitigateBadSamplesEE", false, true) and
625  edm::ParameterDescription<bool>("gainSwitchUseMaxSampleEB", false, true) and
626  edm::ParameterDescription<bool>("gainSwitchUseMaxSampleEE", false, true) and
627  edm::ParameterDescription<bool>("selectiveBadSampleCriteriaEB", false, true) and
628  edm::ParameterDescription<bool>("selectiveBadSampleCriteriaEE", false, true) and
629  edm::ParameterDescription<double>("addPedestalUncertaintyEB", 0., true) and
630  edm::ParameterDescription<double>("addPedestalUncertaintyEE", 0., true) and
631  edm::ParameterDescription<bool>("simplifiedNoiseModelForGainSwitch", true, true) and
632  edm::ParameterDescription<std::string>("timealgo", "RatioMethod", true) and
633  edm::ParameterDescription<std::vector<double>>("EBtimeFitParameters",
634  {-2.015452e+00,
635  3.130702e+00,
636  -1.234730e+01,
637  4.188921e+01,
638  -8.283944e+01,
639  9.101147e+01,
640  -5.035761e+01,
641  1.105621e+01},
642  true) and
643  edm::ParameterDescription<std::vector<double>>("EEtimeFitParameters",
644  {-2.390548e+00,
645  3.553628e+00,
646  -1.762341e+01,
647  6.767538e+01,
648  -1.332130e+02,
649  1.407432e+02,
650  -7.541106e+01,
651  1.620277e+01},
652  true) and
653  edm::ParameterDescription<std::vector<double>>("EBamplitudeFitParameters", {1.138, 1.652}, true) and
654  edm::ParameterDescription<std::vector<double>>("EEamplitudeFitParameters", {1.890, 1.400}, true) and
655  edm::ParameterDescription<double>("EBtimeFitLimits_Lower", 0.2, true) and
656  edm::ParameterDescription<double>("EBtimeFitLimits_Upper", 1.4, true) and
657  edm::ParameterDescription<double>("EEtimeFitLimits_Lower", 0.2, true) and
658  edm::ParameterDescription<double>("EEtimeFitLimits_Upper", 1.4, true) and
659  edm::ParameterDescription<double>("EBtimeConstantTerm", .6, true) and
660  edm::ParameterDescription<double>("EEtimeConstantTerm", 1.0, true) and
661  edm::ParameterDescription<double>("EBtimeNconst", 28.5, true) and
662  edm::ParameterDescription<double>("EEtimeNconst", 31.8, true) and
663  edm::ParameterDescription<double>("outOfTimeThresholdGain12pEB", 5, true) and
664  edm::ParameterDescription<double>("outOfTimeThresholdGain12mEB", 5, true) and
665  edm::ParameterDescription<double>("outOfTimeThresholdGain61pEB", 5, true) and
666  edm::ParameterDescription<double>("outOfTimeThresholdGain61mEB", 5, true) and
667  edm::ParameterDescription<double>("outOfTimeThresholdGain12pEE", 1000, true) and
668  edm::ParameterDescription<double>("outOfTimeThresholdGain12mEE", 1000, true) and
669  edm::ParameterDescription<double>("outOfTimeThresholdGain61pEE", 1000, true) and
670  edm::ParameterDescription<double>("outOfTimeThresholdGain61mEE", 1000, true) and
671  edm::ParameterDescription<double>("amplitudeThresholdEB", 10, true) and
672  edm::ParameterDescription<double>("amplitudeThresholdEE", 10, true) and
673  edm::ParameterDescription<double>("ebSpikeThreshold", 1.042, true) and
674  edm::ParameterDescription<std::vector<double>>(
675  "ebPulseShape", {5.2e-05, -5.26e-05, 6.66e-05, 0.1168, 0.7575, 1., 0.8876, 0.6732, 0.4741, 0.3194}, true) and
676  edm::ParameterDescription<std::vector<double>>(
677  "eePulseShape", {5.2e-05, -5.26e-05, 6.66e-05, 0.1168, 0.7575, 1., 0.8876, 0.6732, 0.4741, 0.3194}, true) and
678  edm::ParameterDescription<bool>("kPoorRecoFlagEB", true, true) and
679  edm::ParameterDescription<bool>("kPoorRecoFlagEE", false, true) and
680  edm::ParameterDescription<double>("chi2ThreshEB_", 65.0, true) and
681  edm::ParameterDescription<double>("chi2ThreshEE_", 50.0, true) and
682  edm::ParameterDescription<double>("crossCorrelationStartTime", -25.0, true) and
683  edm::ParameterDescription<double>("crossCorrelationStopTime", 25.0, true) and
684  edm::ParameterDescription<double>("crossCorrelationTargetTimePrecision", 0.01, true) and
685  edm::ParameterDescription<edm::ParameterSetDescription>("EcalPulseShapeParameters", psd0, true));
686 
687  return psd;
688 }
689 
696  "EcalUncalibRecHitWorkerMultiFit");
edm::ESGetToken< EcalPulseCovariances, EcalPulseCovariancesRcd > pulseConvariancesToken_
edm::ParameterSetDescription getAlgoDescription() override
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
std::array< SampleMatrixGainArray, 2 > noisecors_
EcalPulseShapesMap EcalPulseShapes
int hashedIndex() const
get a compact index for arrays
Definition: EBDetId.h:82
const edm::EventSetup & c
void computeAmplitude(std::vector< double > &amplitudeFitParameters)
unsigned size(int bx) const
EcalUncalibRecHitMultiFitAlgo multiFitMethod_
std::array< SampleMatrix, NGains > SampleMatrixGainArray
EcalUncalibRecHitTimeWeightsAlgo< EEDataFrame > weightsMethod_endcap_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
edm::ESGetToken< EcalPulseShapes, EcalPulseShapesRcd > pulseShapesToken_
Eigen::Matrix< double, FullSampleVectorSize, FullSampleVectorSize > FullSampleMatrix
double time(const C &dataFrame, const std::vector< double > &amplitudes, const EcalPedestals::Item *aped, const EcalMGPAGainRatio *aGain, const FullSampleVector &fullpulse, const EcalWeightSet::EcalWeightMatrix **weights)
Compute time.
Eigen::Matrix< double, FullSampleVectorSize, 1 > FullSampleVector
std::map< std::pair< EcalXtalGroupId, EcalTDCId >, EcalWeightSet > EcalTBWeightMap
Definition: EcalTBWeights.h:18
void computeTime(std::vector< double > &timeFitParameters, std::pair< double, double > &timeFitLimits, std::vector< double > &amplitudeFitParameters)
edm::ESGetToken< EcalGainRatios, EcalGainRatiosRcd > gainsToken_
EcalCondObjectContainer< EcalXtalGroupId > EcalWeightXtalGroups
static const int TEMPLATESAMPLES
edm::ESHandle< EcalWeightXtalGroups > grps
void push_back(T const &t)
ParameterDescriptionNode * addNode(ParameterDescriptionNode const &node)
math::Matrix< 3, 10 >::type EcalWeightMatrix
Definition: EcalWeightSet.h:19
const_iterator begin() const
The iterator returned can not safely be used across threads.
Log< level::Error, false > LogError
edm::ESGetToken< EcalSamplesCorrelation, EcalSamplesCorrelationRcd > noiseConvariancesToken_
edm::EDGetTokenT< unsigned int > bunchSpacing_
EcalUncalibRecHitRatioMethodAlgo< EBDataFrame > ratioMethod_barrel_
tuple result
Definition: mps_fire.py:311
edm::ESGetToken< EcalPedestals, EcalPedestalsRcd > pedsToken_
int hashedIndex(int ieta, int iphi)
Definition: EcalPyUtils.cc:36
edm::ESGetToken< EcalTBWeights, EcalTBWeightsRcd > wgtsToken_
EcalGainRatioMap EcalGainRatios
EcalUncalibRecHitRatioMethodAlgo< EEDataFrame > ratioMethod_endcap_
T sqrt(T t)
Definition: SSEVec.h:19
edm::ESHandle< EcalPulseShapes > pulseshapes
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
EcalPulseCovariancesMap EcalPulseCovariances
ParameterDescriptionNode * addOptionalNode(ParameterDescriptionNode const &node, bool writeToCfi)
edm::ESHandle< EcalTimeCalibConstants > itime
double timeCorrection(float ampli, const std::vector< float > &amplitudeBins, const std::vector< float > &shiftBins)
float gain6Over1() const
std::unique_ptr< EcalUncalibRecHitTimingCCAlgo > computeCC_
EcalWeightMatrix & getWeightsAfterGainSwitch()
Definition: EcalWeightSet.h:27
EcalUncalibratedRecHit makeRecHit(const EcalDataFrame &dataFrame, const EcalPedestals::Item *aped, const EcalMGPAGainRatio *aGain, const SampleMatrixGainArray &noisecors, const FullSampleVector &fullpulse, const FullSampleMatrix &fullpulsecov, const BXVector &activeBX)
compute rechits
edm::ESHandle< EcalTimeBiasCorrections > timeCorrBias_
Definition: DetId.h:17
EcalPedestalsMap EcalPedestals
Definition: EcalPedestals.h:50
edm::ESGetToken< EcalTimeOffsetConstant, EcalTimeOffsetConstantRcd > offtimeToken_
void set(const edm::EventSetup &es) override
edm::ESHandle< EcalPulseCovariances > pulsecovariances
constexpr int gainId(sample_type sample)
get the gainId (2 bits)
int hashedIndex() const
Definition: EEDetId.h:183
EcalWeightMatrix & getWeightsBeforeGainSwitch()
Definition: EcalWeightSet.h:26
EcalTimeCalibConstantMap EcalTimeCalibConstants
std::vector< Item >::const_iterator const_iterator
void resize(int bx, unsigned size)
T const * product() const
Definition: ESHandle.h:86
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
unsigned int id() const
float gain12Over6() const
edm::ESGetToken< EcalTimeBiasCorrections, EcalTimeBiasCorrectionsRcd > timeCorrBiasToken_
Eigen::Matrix< double, SampleVectorSize, SampleVectorSize > SampleMatrix
EcalUncalibRecHitTimeWeightsAlgo< EBDataFrame > weightsMethod_barrel_
edm::ESHandle< EcalTimeOffsetConstant > offtime
const_iterator end() const
float EcalTimeCalibConstant
size_type size() const
edm::ESGetToken< EcalSampleMask, EcalSampleMaskRcd > sampleMaskToken_
edm::ESGetToken< EcalTimeCalibConstants, EcalTimeCalibConstantsRcd > itimeToken_
void reserve(size_type n)
#define DEFINE_EDM_PLUGIN(factory, type, name)
void run(const edm::Event &evt, const EcalDigiCollection &digis, EcalUncalibratedRecHitCollection &result) override
const SampleMatrix & noisecor(bool barrel, int gain) const
static constexpr int MAXSAMPLES
Definition: EcalDataFrame.h:48
edm::ESHandle< EcalSampleMask > sampleMaskHand_
const EcalWeightSet::EcalWeightMatrix * weights[2]
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:157
edm::ESHandle< EcalSamplesCorrelation > noisecovariances
void init(const C &dataFrame, const EcalSampleMask &sampleMask, const double *pedestals, const double *pedestalRMSes, const double *gainRatios)
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
edm::ESGetToken< EcalWeightXtalGroups, EcalWeightXtalGroupsRcd > grpsToken_
const_reference back() const