File indexing completed on 2023-03-17 11:19:02
0001 #include "RecoLocalFastTime/FTLCommonAlgos/interface/MTDUncalibratedRecHitAlgoBase.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003
0004 #include "CommonTools/Utils/interface/FormulaEvaluator.h"
0005
0006 class ETLUncalibRecHitAlgo : public ETLUncalibratedRecHitAlgoBase {
0007 public:
0008
0009 ETLUncalibRecHitAlgo(const edm::ParameterSet& conf, edm::ConsumesCollector& sumes)
0010 : MTDUncalibratedRecHitAlgoBase<ETLDataFrame>(conf, sumes),
0011 adcNBits_(conf.getParameter<uint32_t>("adcNbits")),
0012 adcSaturation_(conf.getParameter<double>("adcSaturation")),
0013 adcLSB_(adcSaturation_ / (1 << adcNBits_)),
0014 toaLSBToNS_(conf.getParameter<double>("toaLSB_ns")),
0015 tofDelay_(conf.getParameter<double>("tofDelay")),
0016 timeError_(conf.getParameter<std::string>("timeResolutionInNs")) {}
0017
0018
0019 ~ETLUncalibRecHitAlgo() override {}
0020
0021
0022 void getEvent(const edm::Event&) final {}
0023 void getEventSetup(const edm::EventSetup&) final {}
0024
0025
0026 FTLUncalibratedRecHit makeRecHit(const ETLDataFrame& dataFrame) const final;
0027
0028 private:
0029 const uint32_t adcNBits_;
0030 const double adcSaturation_;
0031 const double adcLSB_;
0032 const double toaLSBToNS_;
0033 const double tofDelay_;
0034 const reco::FormulaEvaluator timeError_;
0035 };
0036
0037 FTLUncalibratedRecHit ETLUncalibRecHitAlgo::makeRecHit(const ETLDataFrame& dataFrame) const {
0038 constexpr int iSample = 2;
0039 const auto& sample = dataFrame.sample(iSample);
0040
0041 const std::array<double, 1> amplitudeV = {{double(sample.data()) * adcLSB_}};
0042
0043
0044 double time = double(sample.toa()) * toaLSBToNS_ - tofDelay_;
0045 unsigned char flag = 0;
0046
0047 LogDebug("ETLUncalibRecHit") << "ADC+: set the charge to: " << amplitudeV[0] << ' ' << sample.data() << ' ' << adcLSB_
0048 << ' ' << std::endl;
0049 LogDebug("ETLUncalibRecHit") << "ADC+: set the time to: " << time << ' ' << sample.toa() << ' ' << toaLSBToNS_ << ' '
0050 << std::endl;
0051 LogDebug("ETLUncalibRecHit") << "Final uncalibrated amplitude : " << amplitudeV[0] << std::endl;
0052
0053 const std::array<double, 1> emptyV = {{0.}};
0054 double timeError = timeError_.evaluate(amplitudeV, emptyV);
0055
0056 return FTLUncalibratedRecHit(dataFrame.id(),
0057 dataFrame.row(),
0058 dataFrame.column(),
0059 {amplitudeV[0], 0.f},
0060 {time, 0.f},
0061 timeError,
0062 -1.f,
0063 -1.f,
0064 flag);
0065 }
0066
0067 #include "FWCore/Framework/interface/MakerMacros.h"
0068 DEFINE_EDM_PLUGIN(ETLUncalibratedRecHitAlgoFactory, ETLUncalibRecHitAlgo, "ETLUncalibRecHitAlgo");