CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
EcalUncalibRecHitWorkerMultiFit.cc
Go to the documentation of this file.
2 
8 
20 
24 
27  noisecorEBg12(SampleMatrix::Zero()), noisecorEEg12(SampleMatrix::Zero()),
28  noisecorEBg6(SampleMatrix::Zero()), noisecorEEg6(SampleMatrix::Zero()),
29  noisecorEBg1(SampleMatrix::Zero()), noisecorEEg1(SampleMatrix::Zero()),
30  fullpulseEB(FullSampleVector::Zero()),fullpulseEE(FullSampleVector::Zero()),
31  fullpulsecovEB(FullSampleMatrix::Zero()),fullpulsecovEE(FullSampleMatrix::Zero()) {
32 
33  // get the BX for the pulses to be activated
34  std::vector<int32_t> activeBXs = ps.getParameter< std::vector<int32_t> >("activeBXs");
35  activeBX.resize(activeBXs.size());
36  for (unsigned int ibx=0; ibx<activeBXs.size(); ++ibx) {
37  activeBX.coeffRef(ibx) = activeBXs[ibx];
38  }
39 
40  // uncertainty calculation (CPU intensive)
41  ampErrorCalculation_ = ps.getParameter<bool>("ampErrorCalculation");
42  useLumiInfoRunHeader_ = ps.getParameter<bool>("useLumiInfoRunHeader");
43 
44  if (useLumiInfoRunHeader_) {
45  bunchSpacing_ = c.consumes<int>(edm::InputTag("addPileupInfo","bunchSpacing"));
46  }
47 
48  doPrefitEB_ = ps.getParameter<bool>("doPrefitEB");
49  doPrefitEE_ = ps.getParameter<bool>("doPrefitEE");
50 
51  prefitMaxChiSqEB_ = ps.getParameter<double>("prefitMaxChiSqEB");
52  prefitMaxChiSqEE_ = ps.getParameter<double>("prefitMaxChiSqEE");
53 
54  // algorithm to be used for timing
55  timealgo_ = ps.getParameter<std::string>("timealgo");
56 
57  // ratio method parameters
58  EBtimeFitParameters_ = ps.getParameter<std::vector<double> >("EBtimeFitParameters");
59  EEtimeFitParameters_ = ps.getParameter<std::vector<double> >("EEtimeFitParameters");
60  EBamplitudeFitParameters_ = ps.getParameter<std::vector<double> >("EBamplitudeFitParameters");
61  EEamplitudeFitParameters_ = ps.getParameter<std::vector<double> >("EEamplitudeFitParameters");
62  EBtimeFitLimits_.first = ps.getParameter<double>("EBtimeFitLimits_Lower");
63  EBtimeFitLimits_.second = ps.getParameter<double>("EBtimeFitLimits_Upper");
64  EEtimeFitLimits_.first = ps.getParameter<double>("EEtimeFitLimits_Lower");
65  EEtimeFitLimits_.second = ps.getParameter<double>("EEtimeFitLimits_Upper");
66  EBtimeConstantTerm_=ps.getParameter<double>("EBtimeConstantTerm");
67  EEtimeConstantTerm_=ps.getParameter<double>("EEtimeConstantTerm");
68  EBtimeNconst_=ps.getParameter<double>("EBtimeNconst");
69  EEtimeNconst_=ps.getParameter<double>("EEtimeNconst");
70  outOfTimeThreshG12pEB_ = ps.getParameter<double>("outOfTimeThresholdGain12pEB");
71  outOfTimeThreshG12mEB_ = ps.getParameter<double>("outOfTimeThresholdGain12mEB");
72  outOfTimeThreshG61pEB_ = ps.getParameter<double>("outOfTimeThresholdGain61pEB");
73  outOfTimeThreshG61mEB_ = ps.getParameter<double>("outOfTimeThresholdGain61mEB");
74  outOfTimeThreshG12pEE_ = ps.getParameter<double>("outOfTimeThresholdGain12pEE");
75  outOfTimeThreshG12mEE_ = ps.getParameter<double>("outOfTimeThresholdGain12mEE");
76  outOfTimeThreshG61pEE_ = ps.getParameter<double>("outOfTimeThresholdGain61pEE");
77  outOfTimeThreshG61mEE_ = ps.getParameter<double>("outOfTimeThresholdGain61mEE");
78  amplitudeThreshEB_ = ps.getParameter<double>("amplitudeThresholdEB");
79  amplitudeThreshEE_ = ps.getParameter<double>("amplitudeThresholdEE");
80 
81  // spike threshold
82  ebSpikeThresh_ = ps.getParameter<double>("ebSpikeThreshold");
83 
84  ebPulseShape_ = ps.getParameter<std::vector<double> >("ebPulseShape");
85  eePulseShape_ = ps.getParameter<std::vector<double> >("eePulseShape");
86 
87  // chi2 parameters for flags determination
88  kPoorRecoFlagEB_ = ps.getParameter<bool>("kPoorRecoFlagEB");
89  kPoorRecoFlagEE_ = ps.getParameter<bool>("kPoorRecoFlagEE");;
90  chi2ThreshEB_=ps.getParameter<double>("chi2ThreshEB_");
91  chi2ThreshEE_=ps.getParameter<double>("chi2ThreshEE_");
92 
93 }
94 
95 
96 
97 void
99 {
100 
101  // common setup
102  es.get<EcalGainRatiosRcd>().get(gains);
103  es.get<EcalPedestalsRcd>().get(peds);
104 
105  // for the multifit method
110 
111  // weights parameters for the time
112  es.get<EcalWeightXtalGroupsRcd>().get(grps);
113  es.get<EcalTBWeightsRcd>().get(wgts);
114 
115  // which of the samples need be used
117 
118  // for the ratio method
121 
122  // for the time correction methods
124 }
125 
126 void
128 {
129 
130  if (useLumiInfoRunHeader_) {
131 
132  int bunchspacing = 450;
133 
134  if (evt.isRealData()) {
135  edm::RunNumber_t run = evt.run();
136  if (run == 178003 ||
137  run == 178004 ||
138  run == 209089 ||
139  run == 209106 ||
140  run == 209109 ||
141  run == 209146 ||
142  run == 209148 ||
143  run == 209151) {
144  bunchspacing = 25;
145  }
146  else {
147  bunchspacing = 50;
148  }
149  }
150  else {
151  edm::Handle<int> bunchSpacingH;
152  evt.getByToken(bunchSpacing_,bunchSpacingH);
153  bunchspacing = *bunchSpacingH;
154  }
155 
156  if (bunchspacing == 25) {
157  activeBX.resize(10);
158  activeBX << -5,-4,-3,-2,-1,0,1,2,3,4;
159  }
160  else {
161  //50ns configuration otherwise (also for no pileup)
162  activeBX.resize(5);
163  activeBX << -4,-2,0,2,4;
164  }
165  }
166 
167 }
168 
179  float ampli,
180  const std::vector<float>& amplitudeBins,
181  const std::vector<float>& shiftBins) {
182 
183  // computed initially in ns. Than turned in the BX's, as
184  // EcalUncalibratedRecHit need be.
185  double theCorrection = 0;
186 
187  // sanity check for arrays
188  if (amplitudeBins.size() == 0) {
189  edm::LogError("EcalRecHitError")
190  << "timeCorrAmplitudeBins is empty, forcing no time bias corrections.";
191 
192  return 0;
193  }
194 
195  if (amplitudeBins.size() != shiftBins.size()) {
196  edm::LogError("EcalRecHitError")
197  << "Size of timeCorrAmplitudeBins different from "
198  "timeCorrShiftBins. Forcing no time bias corrections. ";
199 
200  return 0;
201  }
202 
203  int myBin = -1;
204  for (int bin = 0; bin < (int) amplitudeBins.size(); bin++) {
205  if (ampli > amplitudeBins.at(bin)) {
206  myBin = bin;
207  } else {
208  break;
209  }
210  }
211 
212  if (myBin == -1) {
213  theCorrection = shiftBins.at(0);
214  } else if (myBin == ((int)(amplitudeBins.size() - 1))) {
215  theCorrection = shiftBins.at(myBin);
216  } else if (-1 < myBin && myBin < ((int) amplitudeBins.size() - 1)) {
217  // interpolate linearly between two assingned points
218  theCorrection = (shiftBins.at(myBin + 1) - shiftBins.at(myBin));
219  theCorrection *= (((double) ampli) - amplitudeBins.at(myBin)) /
220  (amplitudeBins.at(myBin + 1) - amplitudeBins.at(myBin));
221  theCorrection += shiftBins.at(myBin);
222  } else {
223  edm::LogError("EcalRecHitError")
224  << "Assigning time correction impossible. Setting it to 0 ";
225  theCorrection = 0.;
226  }
227 
228  // convert ns into clocks
229  return theCorrection / 25.;
230 }
231 
232 
233 
234 bool
238 {
239  DetId detid(itdg->id());
240 
241  const EcalSampleMask *sampleMask_ = sampleMaskHand_.product();
242 
243  // intelligence for recHit computation
244  EcalUncalibratedRecHit uncalibRecHit;
245  float offsetTime = 0;
246 
247  const EcalPedestals::Item * aped = 0;
248  const EcalMGPAGainRatio * aGain = 0;
249  const EcalXtalGroupId * gid = 0;
250  const EcalPulseShapes::Item * aPulse = 0;
251  const EcalPulseCovariances::Item * aPulseCov = 0;
252 
253  if (detid.subdetId()==EcalEndcap) {
254  unsigned int hashedIndex = EEDetId(detid).hashedIndex();
255  aped = &peds->endcap(hashedIndex);
256  aGain = &gains->endcap(hashedIndex);
257  gid = &grps->endcap(hashedIndex);
258  aPulse = &pulseshapes->endcap(hashedIndex);
259  aPulseCov = &pulsecovariances->endcap(hashedIndex);
262  offsetTime = offtime->getEEValue();
263  } else {
264  unsigned int hashedIndex = EBDetId(detid).hashedIndex();
265  aped = &peds->barrel(hashedIndex);
266  aGain = &gains->barrel(hashedIndex);
267  gid = &grps->barrel(hashedIndex);
268  aPulse = &pulseshapes->barrel(hashedIndex);
269  aPulseCov = &pulsecovariances->barrel(hashedIndex);
272  offsetTime = offtime->getEBValue();
273  }
274 
275  pedVec[0] = aped->mean_x12;
276  pedVec[1] = aped->mean_x6;
277  pedVec[2] = aped->mean_x1;
278  pedRMSVec[0] = aped->rms_x12;
279  pedRMSVec[1] = aped->rms_x6;
280  pedRMSVec[2] = aped->rms_x1;
281  gainRatios[0] = 1.;
282  gainRatios[1] = aGain->gain12Over6();
283  gainRatios[2] = aGain->gain6Over1()*aGain->gain12Over6();
284 
285  int nnoise = noisecovariances->EBG12SamplesCorrelation.size();
286  for (int i=0; i<nnoise; ++i) {
287  for (int j=0; j<nnoise; ++j) {
288  int vidx = std::abs(j-i);
289  noisecorEBg12(i,j) = noisecovariances->EBG12SamplesCorrelation[vidx];
290  noisecorEEg12(i,j) = noisecovariances->EEG12SamplesCorrelation[vidx];
291  noisecorEBg6(i,j) = noisecovariances->EBG6SamplesCorrelation[vidx];
292  noisecorEEg6(i,j) = noisecovariances->EEG6SamplesCorrelation[vidx];
293  noisecorEBg1(i,j) = noisecovariances->EBG1SamplesCorrelation[vidx];
294  noisecorEEg1(i,j) = noisecovariances->EEG1SamplesCorrelation[vidx];
295  }
296  }
297 
298  for (int i=0; i<EcalPulseShape::TEMPLATESAMPLES; ++i) {
299  fullpulseEB(i+7) = aPulse->pdfval[i];
300  fullpulseEE(i+7) = aPulse->pdfval[i];
301  }
302 
303  for(int k=0; k<std::pow(EcalPulseShape::TEMPLATESAMPLES,2); ++k) {
306  fullpulsecovEB(i+7,j+7) = aPulseCov->covval[i][j];
307  fullpulsecovEE(i+7,j+7) = aPulseCov->covval[i][j];
308  }
309 
310  // compute the right bin of the pulse shape using time calibration constants
312  EcalTimeCalibConstant itimeconst = 0;
313  if( it != itime->end() ) {
314  itimeconst = (*it);
315  } else {
316  edm::LogError("EcalRecHitError") << "No time intercalib const found for xtal "
317  << detid.rawId()
318  << "! something wrong with EcalTimeCalibConstants in your DB? ";
319  }
320 
321  // === amplitude computation ===
322  int leadingSample = ((EcalDataFrame)(*itdg)).lastUnsaturatedSample();
323 
324  if ( leadingSample == 4 ) { // saturation on the expected max sample
325  uncalibRecHit = EcalUncalibratedRecHit( (*itdg).id(), 4095*12, 0, 0, 0);
327  // do not propagate the default chi2 = -1 value to the calib rechit (mapped to 64), set it to 0 when saturation
328  uncalibRecHit.setChi2(0);
329  } else if ( leadingSample >= 0 ) { // saturation on other samples: cannot extrapolate from the fourth one
330  double pedestal = 0.;
331  double gainratio = 1.;
332  int gainId = ((EcalDataFrame)(*itdg)).sample(5).gainId();
333 
334  if (gainId==0 || gainId==3) {
335  pedestal = aped->mean_x1;
336  gainratio = aGain->gain6Over1()*aGain->gain12Over6();
337  }
338  else if (gainId==1) {
339  pedestal = aped->mean_x12;
340  gainratio = 1.;
341  }
342  else if (gainId==2) {
343  pedestal = aped->mean_x6;
344  gainratio = aGain->gain12Over6();
345  }
346  double amplitude = ((double)(((EcalDataFrame)(*itdg)).sample(5).adc()) - pedestal) * gainratio;
347  uncalibRecHit = EcalUncalibratedRecHit( (*itdg).id(), amplitude, 0, 0, 0);
349  // do not propagate the default chi2 = -1 value to the calib rechit (mapped to 64), set it to 0 when saturation
350  uncalibRecHit.setChi2(0);
351  } else {
352  // multifit
353  bool barrel = detid.subdetId()==EcalBarrel;
354  int gain = 12;
355  if (((EcalDataFrame)(*itdg)).hasSwitchToGain6()) {
356  gain = 6;
357  }
358  if (((EcalDataFrame)(*itdg)).hasSwitchToGain1()) {
359  gain = 1;
360  }
361  const SampleMatrix &noisecormat = noisecor(barrel,gain);
362  const FullSampleVector &fullpulse = barrel ? fullpulseEB : fullpulseEE;
363  const FullSampleMatrix &fullpulsecov = barrel ? fullpulsecovEB : fullpulsecovEE;
364 
365  uncalibRecHit = multiFitMethod_.makeRecHit(*itdg, aped, aGain, noisecormat,fullpulse,fullpulsecov,activeBX);
366 
367  // === time computation ===
368  if(timealgo_.compare("RatioMethod")==0) {
369  // ratio method
370  float const clockToNsConstant = 25.;
371  if (detid.subdetId()==EcalEndcap) {
372  ratioMethod_endcap_.init( *itdg, *sampleMask_, pedVec, pedRMSVec, gainRatios );
376  double theTimeCorrectionEE = timeCorrection(uncalibRecHit.amplitude(),
377  timeCorrBias_->EETimeCorrAmplitudeBins, timeCorrBias_->EETimeCorrShiftBins);
378 
379  uncalibRecHit.setJitter( crh.timeMax - 5 + theTimeCorrectionEE);
380  uncalibRecHit.setJitterError( std::sqrt(pow(crh.timeError,2) + std::pow(EEtimeConstantTerm_,2)/std::pow(clockToNsConstant,2)) );
381 
382  // consider flagging as kOutOfTime only if above noise
383  if (uncalibRecHit.amplitude() > pedRMSVec[0] * amplitudeThreshEE_){
384  float outOfTimeThreshP = outOfTimeThreshG12pEE_;
385  float outOfTimeThreshM = outOfTimeThreshG12mEE_;
386  // determine if gain has switched away from gainId==1 (x12 gain)
387  // and determine cuts (number of 'sigmas') to ose for kOutOfTime
388  // >3k ADC is necessasry condition for gain switch to occur
389  if (uncalibRecHit.amplitude() > 3000.){
390  for (int iSample = 0; iSample < EEDataFrame::MAXSAMPLES; iSample++) {
391  int GainId = ((EcalDataFrame)(*itdg)).sample(iSample).gainId();
392  if (GainId!=1) {
393  outOfTimeThreshP = outOfTimeThreshG61pEE_;
394  outOfTimeThreshM = outOfTimeThreshG61mEE_;
395  break;
396  }
397  }}
398  float correctedTime = (crh.timeMax-5) * clockToNsConstant + itimeconst + offsetTime;
399  float cterm = EEtimeConstantTerm_;
400  float sigmaped = pedRMSVec[0]; // approx for lower gains
401  float nterm = EEtimeNconst_*sigmaped/uncalibRecHit.amplitude();
402  float sigmat = std::sqrt( nterm*nterm + cterm*cterm );
403  if ( ( correctedTime > sigmat*outOfTimeThreshP ) ||
404  ( correctedTime < (-1.*sigmat*outOfTimeThreshM) ))
405  { uncalibRecHit.setFlagBit( EcalUncalibratedRecHit::kOutOfTime ); }
406  }
407 
408  } else {
409  ratioMethod_barrel_.init( *itdg, *sampleMask_, pedVec, pedRMSVec, gainRatios );
414 
415  double theTimeCorrectionEB = timeCorrection(uncalibRecHit.amplitude(),
416  timeCorrBias_->EBTimeCorrAmplitudeBins, timeCorrBias_->EBTimeCorrShiftBins);
417 
418  uncalibRecHit.setJitter( crh.timeMax - 5 + theTimeCorrectionEB);
419  uncalibRecHit.setJitterError( std::sqrt(std::pow(crh.timeError,2) + std::pow(EBtimeConstantTerm_,2)/std::pow(clockToNsConstant,2)) );
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  float correctedTime = (crh.timeMax-5) * clockToNsConstant + itimeconst + offsetTime;
437  float cterm = EBtimeConstantTerm_;
438  float sigmaped = pedRMSVec[0]; // approx for lower gains
439  float nterm = EBtimeNconst_*sigmaped/uncalibRecHit.amplitude();
440  float sigmat = std::sqrt( nterm*nterm + cterm*cterm );
441  if ( ( correctedTime > sigmat*outOfTimeThreshP ) ||
442  ( correctedTime < (-1.*sigmat*outOfTimeThreshM) ))
443  { uncalibRecHit.setFlagBit( EcalUncalibratedRecHit::kOutOfTime ); }
444  }
445 
446  }
447  } else if(timealgo_.compare("WeightsMethod")==0) {
448  // weights method on the PU subtracted pulse shape
449  std::vector<double> amplitudes;
450  for(unsigned int ibx=0; ibx<activeBX.size(); ++ibx) amplitudes.push_back(uncalibRecHit.outOfTimeAmplitude(ibx));
451 
452  EcalTBWeights::EcalTDCId tdcid(1);
453  EcalTBWeights::EcalTBWeightMap const & wgtsMap = wgts->getMap();
454  EcalTBWeights::EcalTBWeightMap::const_iterator wit;
455  wit = wgtsMap.find( std::make_pair(*gid,tdcid) );
456  if( wit == wgtsMap.end() ) {
457  edm::LogError("EcalUncalibRecHitError") << "No weights found for EcalGroupId: "
458  << gid->id() << " and EcalTDCId: " << tdcid
459  << "\n skipping digi with id: " << detid.rawId();
460 
461  return false;
462  }
463  const EcalWeightSet& wset = wit->second; // this is the EcalWeightSet
464 
467 
468  weights[0] = &mat1;
469  weights[1] = &mat2;
470 
471  double timerh;
472  if (detid.subdetId()==EcalEndcap) {
473  timerh = weightsMethod_endcap_.time( *itdg, amplitudes, aped, aGain, fullpulse, weights);
474  } else {
475  timerh = weightsMethod_barrel_.time( *itdg, amplitudes, aped, aGain, fullpulse, weights);
476  }
477  uncalibRecHit.setJitter( timerh );
478  uncalibRecHit.setJitterError( 0. ); // not computed with weights
479  } else if(timealgo_.compare("None")==0) {
480  uncalibRecHit.setJitter( 0. );
481  uncalibRecHit.setJitterError( 0. );
482  } else {
483  edm::LogError("EcalUncalibRecHitError") << "No time estimation algorithm called "
484  << timealgo_
485  << "\n setting jitter to 0. and jitter uncertainty to 10000. ";
486 
487  uncalibRecHit.setJitter( 0. );
488  uncalibRecHit.setJitterError( 10000. );
489  }
490  }
491 
492  // set flags if gain switch has occurred
493  if( ((EcalDataFrame)(*itdg)).hasSwitchToGain6() ) uncalibRecHit.setFlagBit( EcalUncalibratedRecHit::kHasSwitchToGain6 );
494  if( ((EcalDataFrame)(*itdg)).hasSwitchToGain1() ) uncalibRecHit.setFlagBit( EcalUncalibratedRecHit::kHasSwitchToGain1 );
495 
496  // set quality flags based on chi2 of the the fit
497  /*
498  if(detid.subdetId()==EcalEndcap) {
499  if(kPoorRecoFlagEE_ && uncalibRecHit.chi2()>chi2ThreshEE_) {
500  bool samplesok = true;
501  for (int sample =0; sample < EcalDataFrame::MAXSAMPLES; ++sample) {
502  if (!sampleMask_->useSampleEE(sample)) {
503  samplesok = false;
504  break;
505  }
506  }
507  if (samplesok) uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kPoorReco);
508  }
509  } else {
510  if(kPoorRecoFlagEB_ && uncalibRecHit.chi2()>chi2ThreshEB_) {
511  bool samplesok = true;
512  for (int sample =0; sample < EcalDataFrame::MAXSAMPLES; ++sample) {
513  if (!sampleMask_->useSampleEB(sample)) {
514  samplesok = false;
515  break;
516  }
517  }
518  if (samplesok) uncalibRecHit.setFlagBit(EcalUncalibratedRecHit::kPoorReco);
519  }
520  }
521  */
522 
523  // put the recHit in the collection
524  if (detid.subdetId()==EcalEndcap) {
525  result.push_back( uncalibRecHit );
526  } else {
527  result.push_back( uncalibRecHit );
528  }
529 
530  return true;
531 }
532 
533 
535  if (barrel) {
536  if (gain==6) {
537  return noisecorEBg6;
538  }
539  else if (gain==1) {
540  return noisecorEBg1;
541  }
542  else {
543  return noisecorEBg12;
544  }
545  }
546  else {
547  if (gain==6) {
548  return noisecorEEg6;
549  }
550  else if (gain==1) {
551  return noisecorEEg1;
552  }
553  else {
554  return noisecorEEg12;
555  }
556  }
557 
558  return noisecorEBg12;
559 
560 }
561 
564 
566  psd0.addNode((edm::ParameterDescription<std::vector<double>>("EBPulseShapeTemplate", {1.13979e-02, 7.58151e-01, 1.00000e+00, 8.87744e-01, 6.73548e-01, 4.74332e-01, 3.19561e-01, 2.15144e-01, 1.47464e-01, 1.01087e-01, 6.93181e-02, 4.75044e-02}, true) and
567  edm::ParameterDescription<std::vector<double>>("EEPulseShapeTemplate", {1.16442e-01, 7.56246e-01, 1.00000e+00, 8.97182e-01, 6.86831e-01, 4.91506e-01, 3.44111e-01, 2.45731e-01, 1.74115e-01, 1.23361e-01, 8.74288e-02, 6.19570e-02}, true)));
568 
569  psd0.addNode((edm::ParameterDescription<std::string>("EEdigiCollection", "", true) and
570  edm::ParameterDescription<std::string>("EBdigiCollection", "", true) and
571  edm::ParameterDescription<std::string>("ESdigiCollection", "", true) and
572  edm::ParameterDescription<bool>("UseLCcorrection", true, false) and
573  edm::ParameterDescription<std::vector<double>>("EBCorrNoiseMatrixG12", {1.00000, 0.71073, 0.55721, 0.46089, 0.40449, 0.35931, 0.33924, 0.32439, 0.31581, 0.30481 }, true) and
574  edm::ParameterDescription<std::vector<double>>("EECorrNoiseMatrixG12", {1.00000, 0.71373, 0.44825, 0.30152, 0.21609, 0.14786, 0.11772, 0.10165, 0.09465, 0.08098 }, true) and
575  edm::ParameterDescription<std::vector<double>>("EBCorrNoiseMatrixG06", {1.00000, 0.70946, 0.58021, 0.49846, 0.45006, 0.41366, 0.39699, 0.38478, 0.37847, 0.37055 }, true) and
576  edm::ParameterDescription<std::vector<double>>("EECorrNoiseMatrixG06", {1.00000, 0.71217, 0.47464, 0.34056, 0.26282, 0.20287, 0.17734, 0.16256, 0.15618, 0.14443 }, true) and
577  edm::ParameterDescription<std::vector<double>>("EBCorrNoiseMatrixG01", {1.00000, 0.73354, 0.64442, 0.58851, 0.55425, 0.53082, 0.51916, 0.51097, 0.50732, 0.50409 }, true) and
578  edm::ParameterDescription<std::vector<double>>("EECorrNoiseMatrixG01", {1.00000, 0.72698, 0.62048, 0.55691, 0.51848, 0.49147, 0.47813, 0.47007, 0.46621, 0.46265 }, true) and
579  edm::ParameterDescription<bool>("EcalPreMixStage1", false, true) and
580  edm::ParameterDescription<bool>("EcalPreMixStage2", false, true)));
581 
582  psd0.addOptionalNode((edm::ParameterDescription<std::vector<double>>("EBPulseShapeCovariance", {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, -5.296e-07, 0.000e+00, 0.000e+00,
583  1.233e-05, 6.154e-05, 0.000e+00, -2.200e-05, -2.309e-05, -1.838e-05, -1.373e-05, -7.334e-06, -5.088e-06, -3.745e-06, -2.428e-06, 0.000e+00,
584  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, 0.000e+00, 0.000e+00, 0.000e+00,
585  -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, 1.223e-06, 8.245e-07, 5.589e-07,
586  -4.571e-06, -2.309e-05, 0.000e+00, 8.545e-06, 9.182e-06, 7.219e-06, 5.388e-06, 2.853e-06, 1.944e-06, 1.324e-06, 9.083e-07, 6.335e-07,
587  -3.614e-06, -1.838e-05, 0.000e+00, 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,
588  -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, 9.079e-07, 6.529e-07, 4.752e-07,
589  -1.286e-06, -7.334e-06, 0.000e+00, 2.678e-06, 2.853e-06, 2.385e-06, 1.917e-06, 1.375e-06, 9.100e-07, 6.455e-07, 4.693e-07, 3.657e-07,
590  -8.410e-07, -5.088e-06, 0.000e+00, 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,
591  -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, 7.217e-07, 4.862e-07, 3.768e-07,
592  0.000e+00, -2.428e-06, 0.000e+00, 8.245e-07, 9.083e-07, 7.754e-07, 6.529e-07, 4.693e-07, 4.436e-07, 4.862e-07, 6.509e-07, 4.418e-07,
593  0.000e+00, 0.000e+00, 0.000e+00, 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}, true) and
594  edm::ParameterDescription<std::vector<double>>("EEPulseShapeCovariance", {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, -3.411e-06, 0.000e+00, 0.000e+00,
595  3.333e-05, 2.862e-05, 0.000e+00, -1.244e-05, -1.431e-05, -1.233e-05, -1.032e-05, -5.883e-06, -4.154e-06, -2.902e-06, -2.128e-06, 0.000e+00,
596  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, 0.000e+00, 0.000e+00, 0.000e+00,
597  -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, 1.330e-06, 9.186e-07, 6.446e-07,
598  -1.661e-05, -1.431e-05, 0.000e+00, 6.649e-06, 7.966e-06, 6.898e-06, 5.794e-06, 3.157e-06, 2.184e-06, 1.567e-06, 1.084e-06, 7.575e-07,
599  -1.424e-05, -1.233e-05, 0.000e+00, 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, 1.289e-06, 9.020e-07, 6.146e-07,
601  -6.842e-06, -5.883e-06, 0.000e+00, 2.708e-06, 3.157e-06, 2.859e-06, 2.628e-06, 1.863e-06, 1.296e-06, 8.882e-07, 6.108e-07, 4.283e-07,
602  -4.915e-06, -4.154e-06, 0.000e+00, 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,
603  -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, 9.522e-07, 6.717e-07, 4.293e-07,
604  0.000e+00, -2.128e-06, 0.000e+00, 9.186e-07, 1.084e-06, 9.839e-07, 9.020e-07, 6.108e-07, 5.751e-07, 6.717e-07, 7.911e-07, 5.493e-07,
605  0.000e+00, 0.000e+00, 0.000e+00, 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}, true)), true);
606 
608  psd.addNode(edm::ParameterDescription<std::vector<int>>("activeBXs", {-5,-4,-3,-2,-1,0,1,2,3,4}, true) and
609  edm::ParameterDescription<bool>("ampErrorCalculation", true, true) and
610  edm::ParameterDescription<bool>("useLumiInfoRunHeader", true, true) and
611  edm::ParameterDescription<bool>("doPrefitEB", false, true) and
612  edm::ParameterDescription<bool>("doPrefitEE", false, true) and
613  edm::ParameterDescription<double>("prefitMaxChiSqEB", 25., true) and
614  edm::ParameterDescription<double>("prefitMaxChiSqEE", 10., true) and
615  edm::ParameterDescription<std::string>("timealgo", "RatioMethod", true) and
616  edm::ParameterDescription<std::vector<double>>("EBtimeFitParameters", {-2.015452e+00, 3.130702e+00, -1.234730e+01, 4.188921e+01, -8.283944e+01, 9.101147e+01, -5.035761e+01, 1.105621e+01}, true) and
617  edm::ParameterDescription<std::vector<double>>("EEtimeFitParameters", {-2.390548e+00, 3.553628e+00, -1.762341e+01, 6.767538e+01, -1.332130e+02, 1.407432e+02, -7.541106e+01, 1.620277e+01}, true) and
618  edm::ParameterDescription<std::vector<double>>("EBamplitudeFitParameters", {1.138,1.652}, true) and
619  edm::ParameterDescription<std::vector<double>>("EEamplitudeFitParameters", {1.890,1.400}, true) and
620  edm::ParameterDescription<double>("EBtimeFitLimits_Lower", 0.2, true) and
621  edm::ParameterDescription<double>("EBtimeFitLimits_Upper", 1.4, true) and
622  edm::ParameterDescription<double>("EEtimeFitLimits_Lower", 0.2, true) and
623  edm::ParameterDescription<double>("EEtimeFitLimits_Upper", 1.4, true) and
624  edm::ParameterDescription<double>("EBtimeConstantTerm", .6, true) and
625  edm::ParameterDescription<double>("EEtimeConstantTerm", 1.0, true) and
626  edm::ParameterDescription<double>("EBtimeNconst", 28.5, true) and
627  edm::ParameterDescription<double>("EEtimeNconst", 31.8, true) and
628  edm::ParameterDescription<double>("outOfTimeThresholdGain12pEB", 5, true) and
629  edm::ParameterDescription<double>("outOfTimeThresholdGain12mEB", 5, true) and
630  edm::ParameterDescription<double>("outOfTimeThresholdGain61pEB", 5, true) and
631  edm::ParameterDescription<double>("outOfTimeThresholdGain61mEB", 5, true) and
632  edm::ParameterDescription<double>("outOfTimeThresholdGain12pEE", 1000, true) and
633  edm::ParameterDescription<double>("outOfTimeThresholdGain12mEE", 1000, true) and
634  edm::ParameterDescription<double>("outOfTimeThresholdGain61pEE", 1000, true) and
635  edm::ParameterDescription<double>("outOfTimeThresholdGain61mEE", 1000, true) and
636  edm::ParameterDescription<double>("amplitudeThresholdEB", 10, true) and
637  edm::ParameterDescription<double>("amplitudeThresholdEE", 10, true) and
638  edm::ParameterDescription<double>("ebSpikeThreshold", 1.042, true) and
639  edm::ParameterDescription<std::vector<double>>("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
640  edm::ParameterDescription<std::vector<double>>("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
641  edm::ParameterDescription<bool>("kPoorRecoFlagEB", true, true) and
642  edm::ParameterDescription<bool>("kPoorRecoFlagEE", false, true) and
643  edm::ParameterDescription<double>("chi2ThreshEB_", 65.0, true) and
644  edm::ParameterDescription<double>("chi2ThreshEE_", 50.0, true) and
645  edm::ParameterDescription<edm::ParameterSetDescription>("EcalPulseShapeParameters", psd0, true));
646 
647  return psd;
648 }
649 
650 
651 
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
std::map< std::pair< EcalXtalGroupId, EcalTDCId >, EcalWeightSet > EcalTBWeightMap
Definition: EcalTBWeights.h:21
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
int hashedIndex() const
get a compact index for arrays
Definition: EBDetId.h:86
void computeAmplitude(std::vector< double > &amplitudeFitParameters)
unsigned size(int bx) const
EcalUncalibRecHitMultiFitAlgo multiFitMethod_
boost::transform_iterator< IterHelp, boost::counting_iterator< int > > const_iterator
const SampleMatrix & noisecor(bool barrel, int gain) const
void setJitter(float jitter)
int gainId(sample_type sample)
get the gainId (2 bits)
EcalUncalibRecHitTimeWeightsAlgo< EEDataFrame > weightsMethod_endcap_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:449
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.
EcalUncalibratedRecHit makeRecHit(const EcalDataFrame &dataFrame, const EcalPedestals::Item *aped, const EcalMGPAGainRatio *aGain, const SampleMatrix &noisecor, const FullSampleVector &fullpulse, const FullSampleMatrix &fullpulsecov, const BXVector &activeBX)
compute rechits
void computeTime(std::vector< double > &timeFitParameters, std::pair< double, double > &timeFitLimits, std::vector< double > &amplitudeFitParameters)
static const int TEMPLATESAMPLES
edm::ESHandle< EcalWeightXtalGroups > grps
void push_back(T const &t)
ParameterDescriptionNode * addNode(ParameterDescriptionNode const &node)
void setJitterError(float jitterErr)
Eigen::Matrix< double, 19, 19 > FullSampleMatrix
bool isRealData() const
Definition: EventBase.h:64
EcalUncalibRecHitRatioMethodAlgo< EBDataFrame > ratioMethod_barrel_
const unsigned int id() const
bool run(const edm::Event &evt, const EcalDigiCollection::const_iterator &digi, EcalUncalibratedRecHitCollection &result) override
int hashedIndex(int ieta, int iphi)
Definition: EcalPyUtils.cc:42
float outOfTimeAmplitude(int bx) const
EcalUncalibRecHitRatioMethodAlgo< EEDataFrame > ratioMethod_endcap_
T sqrt(T t)
Definition: SSEVec.h:48
edm::ESHandle< EcalPulseShapes > pulseshapes
tuple result
Definition: query.py:137
RunNumber_t run() const
Definition: Event.h:85
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int j
Definition: DBlmapReader.cc:9
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
Eigen::Matrix< double, 19, 1 > FullSampleVector
EcalWeightMatrix & getWeightsAfterGainSwitch()
Definition: EcalWeightSet.h:30
float covval[EcalPulseShape::TEMPLATESAMPLES][EcalPulseShape::TEMPLATESAMPLES]
edm::ESHandle< EcalTimeBiasCorrections > timeCorrBias_
edm::ParameterSetDescription getAlgoDescription()
Definition: DetId.h:18
void set(const edm::EventSetup &es) override
edm::ESHandle< EcalPulseCovariances > pulsecovariances
int hashedIndex() const
Definition: EEDetId.h:182
EcalWeightMatrix & getWeightsBeforeGainSwitch()
Definition: EcalWeightSet.h:29
const T & get() const
Definition: EventSetup.h:55
std::vector< Item >::const_iterator const_iterator
void resize(int bx, unsigned size)
T const * product() const
Definition: ESHandle.h:86
float gain12Over6() const
EcalUncalibRecHitTimeWeightsAlgo< EBDataFrame > weightsMethod_barrel_
edm::ESHandle< EcalTimeOffsetConstant > offtime
float pdfval[TEMPLATESAMPLES]
float EcalTimeCalibConstant
math::Matrix< 3, 10 >::type EcalWeightMatrix
Definition: EcalWeightSet.h:22
#define DEFINE_EDM_PLUGIN(factory, type, name)
unsigned int RunNumber_t
edm::ESHandle< EcalSampleMask > sampleMaskHand_
const EcalWeightSet::EcalWeightMatrix * weights[2]
edm::ESHandle< EcalSamplesCorrelation > noisecovariances
Eigen::Matrix< double, 10, 10 > SampleMatrix
void init(const C &dataFrame, const EcalSampleMask &sampleMask, const double *pedestals, const double *pedestalRMSes, const double *gainRatios)
static const int MAXSAMPLES
Definition: EcalDataFrame.h:48
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40