File indexing completed on 2024-04-06 12:30:33
0001
0002
0003
0004
0005
0006
0007 #include "DataFormats/Common/interface/Handle.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include <map>
0010 #include <memory>
0011
0012
0013 #include "DataMixingSiPixelWorker.h"
0014
0015 using namespace std;
0016
0017 namespace edm {
0018
0019
0020
0021 DataMixingSiPixelWorker::DataMixingSiPixelWorker() {}
0022
0023
0024 DataMixingSiPixelWorker::DataMixingSiPixelWorker(const edm::ParameterSet &ps, edm::ConsumesCollector &&iC)
0025 : label_(ps.getParameter<std::string>("Label"))
0026
0027 {
0028
0029
0030
0031
0032
0033
0034 pixeldigi_collectionSig_ = ps.getParameter<edm::InputTag>("pixeldigiCollectionSig");
0035 pixeldigi_collectionPile_ = ps.getParameter<edm::InputTag>("pixeldigiCollectionPile");
0036 PixelDigiCollectionDM_ = ps.getParameter<std::string>("PixelDigiCollectionDM");
0037
0038 PixelDigiToken_ = iC.consumes<edm::DetSetVector<PixelDigi>>(pixeldigi_collectionSig_);
0039 PixelDigiPToken_ = iC.consumes<edm::DetSetVector<PixelDigi>>(pixeldigi_collectionPile_);
0040
0041
0042 SiHitStorage_.clear();
0043 }
0044
0045
0046 DataMixingSiPixelWorker::~DataMixingSiPixelWorker() {}
0047
0048 void DataMixingSiPixelWorker::addSiPixelSignals(const edm::Event &e) {
0049
0050
0051 LogDebug("DataMixingSiPixelWorker") << "===============> adding MC signals for " << e.id();
0052
0053 Handle<edm::DetSetVector<PixelDigi>> input;
0054
0055 if (e.getByToken(PixelDigiToken_, input)) {
0056
0057 edm::DetSetVector<PixelDigi>::const_iterator DSViter = input->begin();
0058 for (; DSViter != input->end(); DSViter++) {
0059 #ifdef DEBUG
0060 LogDebug("DataMixingSiPixelWorker") << "Processing DetID " << DSViter->id;
0061 #endif
0062
0063 uint32_t detID = DSViter->id;
0064 edm::DetSet<PixelDigi>::const_iterator begin = (DSViter->data).begin();
0065 edm::DetSet<PixelDigi>::const_iterator end = (DSViter->data).end();
0066 edm::DetSet<PixelDigi>::const_iterator icopy;
0067
0068 OneDetectorMap LocalMap;
0069
0070 for (icopy = begin; icopy != end; icopy++) {
0071 LocalMap.insert(OneDetectorMap::value_type((icopy->channel()), *icopy));
0072 }
0073
0074 SiHitStorage_.insert(SiGlobalIndex::value_type(detID, LocalMap));
0075 }
0076 }
0077 }
0078
0079 void DataMixingSiPixelWorker::addSiPixelPileups(const int bcr,
0080 const EventPrincipal *ep,
0081 unsigned int eventNr,
0082 ModuleCallingContext const *mcc) {
0083 LogDebug("DataMixingSiPixelWorker") << "\n===============> adding pileups from event " << ep->id()
0084 << " for bunchcrossing " << bcr;
0085
0086
0087
0088
0089 std::shared_ptr<Wrapper<edm::DetSetVector<PixelDigi>> const> inputPTR =
0090 getProductByTag<edm::DetSetVector<PixelDigi>>(*ep, pixeldigi_collectionPile_, mcc);
0091
0092 if (inputPTR) {
0093 const edm::DetSetVector<PixelDigi> *input = const_cast<edm::DetSetVector<PixelDigi> *>(inputPTR->product());
0094
0095
0096
0097
0098
0099
0100 edm::DetSetVector<PixelDigi>::const_iterator DSViter = input->begin();
0101 for (; DSViter != input->end(); DSViter++) {
0102 #ifdef DEBUG
0103 LogDebug("DataMixingSiPixelWorker") << "Pileups: Processing DetID " << DSViter->id;
0104 #endif
0105
0106 uint32_t detID = DSViter->id;
0107 edm::DetSet<PixelDigi>::const_iterator begin = (DSViter->data).begin();
0108 edm::DetSet<PixelDigi>::const_iterator end = (DSViter->data).end();
0109 edm::DetSet<PixelDigi>::const_iterator icopy;
0110
0111
0112
0113 SiGlobalIndex::const_iterator itest;
0114
0115 itest = SiHitStorage_.find(detID);
0116
0117 if (itest != SiHitStorage_.end()) {
0118
0119 OneDetectorMap LocalMap = itest->second;
0120
0121
0122 for (icopy = begin; icopy != end; icopy++) {
0123 LocalMap.insert(OneDetectorMap::value_type((icopy->channel()), *icopy));
0124 }
0125
0126 SiHitStorage_[detID] = LocalMap;
0127
0128 } else {
0129
0130
0131 OneDetectorMap LocalMap;
0132
0133 for (icopy = begin; icopy != end; icopy++) {
0134 LocalMap.insert(OneDetectorMap::value_type((icopy->channel()), *icopy));
0135 }
0136
0137 SiHitStorage_.insert(SiGlobalIndex::value_type(detID, LocalMap));
0138 }
0139 }
0140 }
0141 }
0142
0143 void DataMixingSiPixelWorker::putSiPixel(edm::Event &e) {
0144
0145
0146 std::vector<edm::DetSet<PixelDigi>> vPixelDigi;
0147
0148
0149
0150
0151
0152
0153 for (SiGlobalIndex::const_iterator IDet = SiHitStorage_.begin(); IDet != SiHitStorage_.end(); IDet++) {
0154 edm::DetSet<PixelDigi> SPD(IDet->first);
0155
0156 OneDetectorMap LocalMap = IDet->second;
0157
0158
0159 int formerPixel = -1;
0160 int currentPixel;
0161 int ADCSum = 0;
0162
0163 OneDetectorMap::const_iterator iLocalchk;
0164
0165 for (OneDetectorMap::const_iterator iLocal = LocalMap.begin(); iLocal != LocalMap.end(); ++iLocal) {
0166 currentPixel = iLocal->first;
0167
0168 if (currentPixel == formerPixel) {
0169 ADCSum += (iLocal->second).adc();
0170 } else {
0171 if (formerPixel != -1) {
0172 if (ADCSum > 511)
0173 ADCSum = 255;
0174 else if (ADCSum > 253 && ADCSum < 512)
0175 ADCSum = 254;
0176 PixelDigi aHit(formerPixel, ADCSum);
0177 SPD.push_back(aHit);
0178 }
0179
0180 formerPixel = currentPixel;
0181 ADCSum = (iLocal->second).adc();
0182 }
0183
0184 iLocalchk = iLocal;
0185 if ((++iLocalchk) == LocalMap.end()) {
0186 if (ADCSum > 511)
0187 ADCSum = 255;
0188 else if (ADCSum > 253 && ADCSum < 512)
0189 ADCSum = 254;
0190 SPD.push_back(PixelDigi(formerPixel, ADCSum));
0191 }
0192
0193 }
0194
0195
0196 vPixelDigi.push_back(SPD);
0197
0198 }
0199
0200
0201 LogInfo("DataMixingSiPixelWorker") << "total # Merged Pixels: " << vPixelDigi.size();
0202
0203
0204
0205 std::unique_ptr<edm::DetSetVector<PixelDigi>> MyPixelDigis(new edm::DetSetVector<PixelDigi>(vPixelDigi));
0206
0207
0208
0209 e.put(std::move(MyPixelDigis), PixelDigiCollectionDM_);
0210
0211
0212 SiHitStorage_.clear();
0213 }
0214
0215 }