File indexing completed on 2022-11-29 00:56:52
0001 #include "RecoLocalFastTime/FTLCommonAlgos/interface/MTDUncalibratedRecHitAlgoBase.h"
0002 #include "RecoLocalFastTime/FTLClusterizer/interface/BTLRecHitsErrorEstimatorIM.h"
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004
0005 #include "CommonTools/Utils/interface/FormulaEvaluator.h"
0006
0007 class BTLUncalibRecHitAlgo : public BTLUncalibratedRecHitAlgoBase {
0008 public:
0009
0010 BTLUncalibRecHitAlgo(const edm::ParameterSet& conf, edm::ConsumesCollector& sumes)
0011 : MTDUncalibratedRecHitAlgoBase<BTLDataFrame>(conf, sumes),
0012 adcNBits_(conf.getParameter<uint32_t>("adcNbits")),
0013 adcSaturation_(conf.getParameter<double>("adcSaturation")),
0014 adcLSB_(adcSaturation_ / (1 << adcNBits_)),
0015 toaLSBToNS_(conf.getParameter<double>("toaLSB_ns")),
0016 timeError_(conf.getParameter<std::string>("timeResolutionInNs")),
0017 timeCorr_p0_(conf.getParameter<double>("timeCorr_p0")),
0018 timeCorr_p1_(conf.getParameter<double>("timeCorr_p1")),
0019 timeCorr_p2_(conf.getParameter<double>("timeCorr_p2")),
0020 c_LYSO_(conf.getParameter<double>("c_LYSO")) {}
0021
0022
0023 ~BTLUncalibRecHitAlgo() override {}
0024
0025
0026 void getEvent(const edm::Event&) final {}
0027 void getEventSetup(const edm::EventSetup&) final {}
0028
0029
0030 FTLUncalibratedRecHit makeRecHit(const BTLDataFrame& dataFrame) const final;
0031
0032 private:
0033 const uint32_t adcNBits_;
0034 const double adcSaturation_;
0035 const double adcLSB_;
0036 const double toaLSBToNS_;
0037 const reco::FormulaEvaluator timeError_;
0038 const double timeCorr_p0_;
0039 const double timeCorr_p1_;
0040 const double timeCorr_p2_;
0041 const double c_LYSO_;
0042 };
0043
0044 FTLUncalibratedRecHit BTLUncalibRecHitAlgo::makeRecHit(const BTLDataFrame& dataFrame) const {
0045
0046
0047
0048
0049
0050
0051 std::pair<float, float> amplitude(0., 0.);
0052 std::pair<float, float> time(0., 0.);
0053
0054 unsigned char flag = 0;
0055
0056 const auto& sampleRight = dataFrame.sample(0);
0057 const auto& sampleLeft = dataFrame.sample(1);
0058
0059 double nHits = 0.;
0060
0061 if (sampleRight.data() > 0) {
0062 amplitude.first = float(sampleRight.data()) * adcLSB_;
0063 time.first = float(sampleRight.toa()) * toaLSBToNS_;
0064
0065 nHits += 1.;
0066
0067
0068 time.first -= timeCorr_p0_ * pow(amplitude.first, timeCorr_p1_) + timeCorr_p2_;
0069 flag |= 0x1;
0070 }
0071
0072
0073 if (sampleLeft.data() > 0) {
0074 amplitude.second = float(sampleLeft.data()) * adcLSB_;
0075 time.second = float(sampleLeft.toa()) * toaLSBToNS_;
0076
0077 nHits += 1.;
0078
0079
0080 time.second -= timeCorr_p0_ * pow(amplitude.second, timeCorr_p1_) + timeCorr_p2_;
0081 flag |= (0x1 << 1);
0082 }
0083
0084
0085
0086 const std::array<double, 1> amplitudeV = {{(amplitude.first + amplitude.second) / nHits}};
0087 const std::array<double, 1> emptyV = {{0.}};
0088
0089 double timeError = (nHits > 0. ? timeError_.evaluate(amplitudeV, emptyV) : -1.);
0090
0091
0092
0093 float position = 0.5f * (c_LYSO_ * (time.second - time.first));
0094 float positionError = BTLRecHitsErrorEstimatorIM::positionError();
0095
0096 LogDebug("BTLUncalibRecHit") << "ADC+: set the charge to: (" << amplitude.first << ", " << amplitude.second << ") ("
0097 << sampleRight.data() << ", " << sampleLeft.data() << " " << adcLSB_ << ' '
0098 << std::endl;
0099 LogDebug("BTLUncalibRecHit") << "TDC+: set the time to: (" << time.first << ", " << time.second << ") ("
0100 << sampleRight.toa() << ", " << sampleLeft.toa() << " " << toaLSBToNS_ << ' '
0101 << std::endl;
0102
0103 return FTLUncalibratedRecHit(
0104 dataFrame.id(), dataFrame.row(), dataFrame.column(), amplitude, time, timeError, position, positionError, flag);
0105 }
0106
0107 #include "FWCore/Framework/interface/MakerMacros.h"
0108 DEFINE_EDM_PLUGIN(BTLUncalibratedRecHitAlgoFactory, BTLUncalibRecHitAlgo, "BTLUncalibRecHitAlgo");