File indexing completed on 2024-04-06 12:29:22
0001 #include <string>
0002 #include <algorithm>
0003 #include <numeric>
0004 #include <functional>
0005
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "Geometry/EcalMapping/interface/EcalElectronicsMapping.h"
0008 #include "Geometry/EcalMapping/interface/EcalMappingRcd.h"
0009 #include "Geometry/EcalAlgo/interface/EcalBarrelGeometry.h"
0010
0011 #include "SimCalorimetry/EcalEBTrigPrimAlgos/interface/EcalEBPhase2TrigPrimAlgo.h"
0012
0013 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0014 #include "DataFormats/EcalDigi/interface/EBDataFrame_Ph2.h"
0015 #include "DataFormats/EcalDetId/interface/EcalTrigTowerDetId.h"
0016 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0017 #include "DataFormats/EcalDetId/interface/EcalTriggerElectronicsId.h"
0018
0019 #include "CondFormats/EcalObjects/interface/EcalTPGPedestals.h"
0020 #include "CondFormats/DataRecord/interface/EcalTPGPedestalsRcd.h"
0021
0022 #include "DataFormats/EcalDigi/interface/EcalConstants.h"
0023
0024 #include <TTree.h>
0025 #include <TMath.h>
0026
0027
0028
0029 const unsigned int EcalEBPhase2TrigPrimAlgo::nrSamples_ =
0030 ecalPh2::sampleSize;
0031 const unsigned int EcalEBPhase2TrigPrimAlgo::maxNrTowers_ = 2448;
0032
0033 EcalEBPhase2TrigPrimAlgo::EcalEBPhase2TrigPrimAlgo(const EcalTrigTowerConstituentsMap *eTTmap,
0034 const CaloGeometry *theGeometry,
0035 int binofmax,
0036 bool debug)
0037 : eTTmap_(eTTmap),
0038 theGeometry_(theGeometry),
0039 binOfMaximum_(binofmax),
0040 debug_(debug)
0041
0042 {
0043 maxNrSamples_ = ecalPh2::sampleSize;
0044 this->init();
0045 }
0046
0047 void EcalEBPhase2TrigPrimAlgo::init() {
0048 theMapping_ = new EcalElectronicsMapping();
0049
0050 initStructures(towerMapEB_);
0051 hitTowers_.resize(maxNrTowers_);
0052
0053 linearizer_ = new EcalEBPhase2Linearizer(debug_);
0054 lin_out_.resize(maxNrSamples_);
0055
0056 amplitude_reconstructor_ = new EcalEBPhase2AmplitudeReconstructor(debug_);
0057 filt_out_.resize(maxNrSamples_);
0058
0059 tpFormatter_ = new EcalEBPhase2TPFormatter(debug_);
0060 outEt_.resize(maxNrSamples_);
0061 outTime_.resize(maxNrSamples_);
0062
0063
0064
0065 time_reconstructor_ = new EcalEBPhase2TimeReconstructor(debug_);
0066 time_out_.resize(maxNrSamples_);
0067 spike_tagger_ = new EcalEBPhase2SpikeTagger(debug_);
0068 }
0069
0070
0071 EcalEBPhase2TrigPrimAlgo::~EcalEBPhase2TrigPrimAlgo() {
0072 delete linearizer_;
0073 delete amplitude_reconstructor_;
0074 delete time_reconstructor_;
0075 delete spike_tagger_;
0076 delete tpFormatter_;
0077 delete theMapping_;
0078 }
0079
0080 void EcalEBPhase2TrigPrimAlgo::run(EBDigiCollectionPh2 const *digi, EcalEBPhase2TrigPrimDigiCollection &result) {
0081 if (debug_)
0082 LogDebug("") << " EcalEBPhase2TrigPrimAlgo: digi size " << digi->size() << std::endl;
0083
0084 EcalEBPhase2TriggerPrimitiveDigi tp;
0085 int firstSample = binOfMaximum_ - 1 - nrSamples_ / 2;
0086 int lastSample = binOfMaximum_ - 1 + nrSamples_ / 2;
0087
0088 if (debug_) {
0089 LogDebug("") << " binOfMaximum_ " << binOfMaximum_ << " nrSamples_" << nrSamples_ << std::endl;
0090 LogDebug("") << " first sample " << firstSample << " last " << lastSample << std::endl;
0091 }
0092
0093 clean(towerMapEB_);
0094 fillMap(digi, towerMapEB_);
0095
0096 int iChannel = 0;
0097 int nXinBCP = 0;
0098 for (int itow = 0; itow < nrTowers_; ++itow) {
0099 int index = hitTowers_[itow].first;
0100 const EcalTrigTowerDetId &thisTower = hitTowers_[itow].second;
0101 if (debug_)
0102 LogDebug("") << " Data for TOWER num " << itow << " index " << index << " TowerId " << thisTower << " zside "
0103 << thisTower.zside() << " ieta " << thisTower.ieta() << " iphi " << thisTower.iphi() << " size "
0104 << towerMapEB_[itow].size() << std::endl;
0105
0106
0107 int nxstals = 0;
0108 for (unsigned int iStrip = 0; iStrip < towerMapEB_[itow].size(); ++iStrip) {
0109 if (debug_)
0110 LogDebug("") << " Data for STRIP num " << iStrip << std::endl;
0111 std::vector<EBDataFrame_Ph2> &dataFrames =
0112 (towerMapEB_[index])[iStrip].second;
0113
0114 nxstals = (towerMapEB_[index])[iStrip].first;
0115 if (nxstals <= 0)
0116 continue;
0117 if (debug_)
0118 LogDebug("") << " Number of xTals " << nxstals << std::endl;
0119
0120
0121
0122
0123
0124 for (int iXstal = 0; iXstal < nxstals; iXstal++) {
0125 const EBDetId &myid = dataFrames[iXstal].id();
0126
0127 nXinBCP++;
0128 if (debug_) {
0129 LogDebug("") << " Data for TOWER num " << itow << " index " << index << " TowerId " << thisTower << " size "
0130 << towerMapEB_[itow].size() << std::endl;
0131 LogDebug("") << "nXinBCP " << nXinBCP << " myid rawId " << myid.rawId() << " xTal iEta " << myid.ieta()
0132 << " iPhi " << myid.iphi() << std::endl;
0133 }
0134
0135 tp = EcalEBPhase2TriggerPrimitiveDigi(myid);
0136 tp.setSize(nrSamples_);
0137
0138 iChannel++;
0139 if (debug_) {
0140 LogDebug("") << " " << std::endl;
0141 LogDebug("") << " ****** iChannel " << iChannel << std::endl;
0142 for (int i = 0; i < dataFrames[iXstal].size(); i++) {
0143 LogDebug("") << " " << dataFrames[iXstal][i].adc();
0144 }
0145 LogDebug("") << " " << std::endl;
0146 }
0147
0148 if (debug_) {
0149 LogDebug("") << std::endl;
0150 EBDetId id = dataFrames[iXstal].id();
0151 LogDebug("") << "iXstal= " << iXstal << std::endl;
0152 LogDebug("") << "iXstal= " << iXstal << " id " << id << " EcalDataFrame_Ph2 is: " << std::endl;
0153 for (int i = 0; i < dataFrames[iXstal].size(); i++) {
0154 LogDebug("") << " " << std::dec << dataFrames[iXstal][i].adc();
0155 }
0156 LogDebug("") << std::endl;
0157 }
0158
0159
0160 this->getLinearizer()->setParameters(dataFrames[iXstal].id(), ecaltpPed_, ecaltpLin_, ecaltpgBadX_);
0161 this->getLinearizer()->process(dataFrames[iXstal], lin_out_);
0162
0163 for (unsigned int i = 0; i < lin_out_.size(); i++) {
0164 if (lin_out_[i] > 0X3FFF)
0165 lin_out_[i] = 0X3FFF;
0166 }
0167
0168 if (debug_) {
0169 LogDebug("") << "EcalEBPhase2TrigPrimAlgo output of linearize for channel " << iXstal << std::endl;
0170 for (unsigned int i = 0; i < lin_out_.size(); i++) {
0171 LogDebug("") << " " << std::dec << lin_out_[i];
0172 }
0173 LogDebug("") << std::endl;
0174 }
0175
0176
0177 this->getSpikeTagger()->setParameters(dataFrames[iXstal].id(), ecaltpPed_, ecaltpLin_, ecaltpgBadX_);
0178 bool isASpike = this->getSpikeTagger()->process(lin_out_);
0179
0180
0181
0182
0183 this->getAmplitudeFinder()->setParameters(myid.rawId(), ecaltpgAmplWeightMap_, ecaltpgWeightGroup_);
0184 this->getAmplitudeFinder()->process(lin_out_, filt_out_);
0185
0186 if (debug_) {
0187 LogDebug("") << "EcalEBPhase2TrigPrimAlgo output of amp finder is a vector of size: " << std::dec
0188 << time_out_.size() << std::endl;
0189 for (unsigned int ix = 0; ix < filt_out_.size(); ix++) {
0190 LogDebug("") << std::dec << filt_out_[ix] << " ";
0191 }
0192 LogDebug("") << std::endl;
0193 }
0194
0195 if (debug_) {
0196 LogDebug("") << " Ampl "
0197 << " ";
0198 for (unsigned int ix = 0; ix < filt_out_.size(); ix++) {
0199 LogDebug("") << std::dec << filt_out_[ix] << " ";
0200 }
0201 LogDebug("") << std::endl;
0202 }
0203
0204
0205 this->getTimeFinder()->setParameters(myid.rawId(), ecaltpgTimeWeightMap_, ecaltpgWeightGroup_);
0206 this->getTimeFinder()->process(lin_out_, filt_out_, time_out_);
0207
0208 if (debug_) {
0209 LogDebug("") << " Time "
0210 << " ";
0211 for (unsigned int ix = 0; ix < time_out_.size(); ix++) {
0212 LogDebug("") << std::dec << time_out_[ix] << " ";
0213 }
0214 LogDebug("") << std::endl;
0215 }
0216
0217 if (debug_) {
0218 LogDebug("") << "EcalEBPhase2TrigPrimAlgo output of timefinder is a vector of size: " << std::dec
0219 << time_out_.size() << std::endl;
0220 for (unsigned int ix = 0; ix < time_out_.size(); ix++) {
0221 LogDebug("") << std::dec << time_out_[ix] << " ";
0222 }
0223 LogDebug("") << std::endl;
0224 }
0225
0226 this->getTPFormatter()->process(filt_out_, time_out_, outEt_, outTime_);
0227
0228 if (debug_) {
0229 LogDebug("") << " compressed Et "
0230 << " ";
0231 for (unsigned int iSample = 0; iSample < outEt_.size(); ++iSample) {
0232 LogDebug("") << outEt_[iSample] << " ";
0233 }
0234 LogDebug("") << std::endl;
0235
0236 LogDebug("") << " compressed time "
0237 << " ";
0238 for (unsigned int iSample = 0; iSample < outEt_.size(); ++iSample) {
0239 LogDebug("") << outTime_[iSample] << " ";
0240 }
0241 LogDebug("") << std::endl;
0242 }
0243
0244 if (debug_) {
0245 LogDebug("") << " EcalEBPhase2TrigPrimAlgo after getting the formatter " << std::endl;
0246 for (unsigned int iSample = 0; iSample < outEt_.size(); ++iSample) {
0247 LogDebug("") << " outEt " << outEt_[iSample] << " outTime " << outTime_[iSample] << " ";
0248 }
0249 LogDebug("") << std::endl;
0250 }
0251
0252
0253
0254
0255 int etInADC = 0;
0256 ;
0257 int64_t time = -999;
0258 int nSam = 0;
0259 for (int iSample = 0; iSample < 16; ++iSample) {
0260 etInADC = outEt_[iSample];
0261 time = outTime_[iSample];
0262 if (debug_) {
0263 LogDebug("") << "TrigPrimAlgo outEt " << outEt_[iSample] << " outTime " << outTime_[iSample] << std::endl;
0264 LogDebug("") << "TrigPrimAlgo etInADCt " << outEt_[iSample] << " outTime " << time << std::endl;
0265 }
0266
0267 tp.setSample(nSam, EcalEBPhase2TriggerPrimitiveSample(etInADC, isASpike, time));
0268 nSam++;
0269 }
0270
0271 result.push_back(tp);
0272
0273 }
0274
0275 }
0276
0277 if (debug_) {
0278 if (nXinBCP > 0)
0279 LogDebug("") << " Accepted xTals " << nXinBCP << std::endl;
0280 }
0281 }
0282 }
0283
0284
0285
0286 int EcalEBPhase2TrigPrimAlgo::findStripNr(const EBDetId &id) {
0287 int stripnr;
0288 int n = ((id.ic() - 1) % 100) / 20;
0289 if (id.ieta() < 0)
0290 stripnr = n + 1;
0291 else
0292 stripnr = nbMaxStrips_ - n;
0293 return stripnr;
0294 }