File indexing completed on 2023-03-17 10:44:02
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022 #include "FWCore/Framework/interface/EventSetup.h"
0023
0024 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0025
0026 #include "DataFormats/SiPixelDigi/interface/SiPixelCalibDigiError.h"
0027 #include "CondFormats/SiPixelObjects/interface/SiPixelFrameConverter.h"
0028 #include "CondFormats/SiPixelObjects/interface/SiPixelFedCablingMap.h"
0029 #include "CondFormats/SiPixelObjects/interface/ElectronicIndex.h"
0030 #include "CondFormats/SiPixelObjects/interface/DetectorIndex.h"
0031 #include "CondFormats/SiPixelObjects/interface/LocalPixel.h"
0032
0033 #include "SiPixelCalibDigiProducer.h"
0034 #include <sstream>
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047 SiPixelCalibDigiProducer::SiPixelCalibDigiProducer(const edm::ParameterSet& iConfig)
0048 : src_(iConfig.getParameter<edm::InputTag>("src")),
0049 iEventCounter_(0),
0050 ignore_non_pattern_(iConfig.getParameter<bool>("ignoreNonPattern")),
0051 control_pattern_size_(iConfig.getParameter<bool>("checkPatternEachEvent")),
0052 includeErrors_(iConfig.getUntrackedParameter<bool>("includeErrors", false)),
0053 errorType(iConfig.getUntrackedParameter<int>("errorTypeNumber", 1)),
0054 conf_(iConfig),
0055 number_of_pixels_per_pattern_(0),
0056 use_realeventnumber_(iConfig.getParameter<bool>("useRealEventNumber"))
0057
0058 {
0059 tPixelDigi = consumes<edm::DetSetVector<PixelDigi>>(src_);
0060
0061 produces<edm::DetSetVector<SiPixelCalibDigi>>();
0062 if (includeErrors_)
0063 produces<edm::DetSetVector<SiPixelCalibDigiError>>();
0064
0065 calibToken_ = esConsumes<SiPixelCalibConfiguration, SiPixelCalibConfigurationRcd>();
0066 trackerGeomToken_ = esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>();
0067 cablingMapToken_ = esConsumes<SiPixelFedCablingMap, SiPixelFedCablingMapRcd>();
0068 }
0069
0070 SiPixelCalibDigiProducer::~SiPixelCalibDigiProducer() {
0071
0072
0073 }
0074
0075
0076
0077
0078
0079
0080
0081
0082 bool SiPixelCalibDigiProducer::store() {
0083
0084 if (iEventCounter_ % pattern_repeat_ == 0) {
0085
0086 return true;
0087 } else if (iEventCounter_ == calib_->expectedTotalEvents())
0088 return true;
0089 else
0090 return false;
0091 return true;
0092 }
0093
0094
0095
0096 void SiPixelCalibDigiProducer::fill(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0097
0098 short icalibpoint = calib_->vcalIndexForEvent(iEventCounter_);
0099 edm::Handle<edm::DetSetVector<PixelDigi>> pixelDigis;
0100 iEvent.getByToken(tPixelDigi, pixelDigis);
0101
0102 edm::LogInfo("SiPixelCalibProducer") << "in fill(), calibpoint " << icalibpoint << " ndigis " << pixelDigis->size()
0103 << std::endl;
0104
0105 edm::DetSetVector<PixelDigi>::const_iterator digiIter;
0106 for (digiIter = pixelDigis->begin(); digiIter != pixelDigis->end(); ++digiIter) {
0107 uint32_t detid = digiIter->id;
0108 edm::DetSet<PixelDigi>::const_iterator ipix;
0109
0110 for (ipix = digiIter->data.begin(); ipix != digiIter->end(); ++ipix) {
0111
0112 fillPixel(detid, ipix->row(), ipix->column(), icalibpoint, ipix->adc());
0113 }
0114 }
0115 }
0116
0117
0118
0119
0120 bool SiPixelCalibDigiProducer::checkFED(uint32_t detid) {
0121
0122
0123 if (detid_to_fedid_[detid])
0124 return true;
0125 for (int fedid = 0; fedid <= 40; ++fedid) {
0126
0127 SiPixelFrameConverter converter(theCablingMap_.product(), fedid);
0128 if (converter.hasDetUnit(detid)) {
0129 detid_to_fedid_[detid] = fedid;
0130 edm::LogInfo("SiPixelCalibDigiProducer")
0131 << "matched detid " << detid << " to fed " << detid_to_fedid_[detid] << std::endl;
0132 return true;
0133 }
0134 }
0135 return false;
0136 }
0137
0138
0139
0140
0141 void SiPixelCalibDigiProducer::fillPixel(uint32_t detid, short row, short col, short ipoint, short adc) {
0142
0143
0144
0145 if (!checkFED(detid)) {
0146 edm::LogError("SiPixelCalibDigiProducer") << " was unable to match detid " << detid << " to a FED!" << std::endl;
0147 return;
0148 }
0149 if (!checkPixel(detid, row, col)) {
0150 return;
0151 }
0152
0153
0154 pixelstruct temppixelworker;
0155 temppixelworker.first = detid;
0156 temppixelworker.second.first = row;
0157 temppixelworker.second.second = col;
0158 std::map<pixelstruct, SiPixelCalibDigi>::const_iterator ipix = intermediate_data_.find(temppixelworker);
0159
0160 if (ipix == intermediate_data_.end()) {
0161 SiPixelCalibDigi tempdigi(calib_->nVCal());
0162 tempdigi.setrowcol(row, col);
0163 intermediate_data_[temppixelworker] = tempdigi;
0164 }
0165
0166 intermediate_data_[temppixelworker].fill(ipoint, adc);
0167 return;
0168 }
0169
0170
0171
0172 void SiPixelCalibDigiProducer::clear() {
0173
0174
0175
0176
0177
0178
0179
0180
0181 uint32_t tempsize = intermediate_data_.size();
0182 if (tempsize > number_of_pixels_per_pattern_) {
0183 edm::LogError("SiPixelCalibDigiProducer") << "Number of pixels in pattern is now: " << tempsize << ", size is was "
0184 << number_of_pixels_per_pattern_ << std::endl;
0185 number_of_pixels_per_pattern_ = tempsize;
0186 }
0187
0188 intermediate_data_.erase(intermediate_data_.begin(), intermediate_data_.end());
0189 intermediate_data_.clear();
0190
0191
0192 error_data_.erase(error_data_.begin(), error_data_.end());
0193 error_data_.clear();
0194 }
0195
0196
0197
0198
0199 void SiPixelCalibDigiProducer::setPattern() {
0200
0201 uint32_t patternnumber = (iEventCounter_ - 1) / pattern_repeat_;
0202 uint32_t rowpatternnumber = patternnumber / calib_->nColumnPatterns();
0203 uint32_t colpatternnumber = patternnumber % calib_->nColumnPatterns();
0204 edm::LogInfo("SiPixelCalibDigiProducer")
0205 << " rowpatternnumbers = " << rowpatternnumber << " " << colpatternnumber << " " << patternnumber << std::endl;
0206
0207 std::vector<short> calibcols = calib_->getColumnPattern();
0208 std::vector<short> calibrows = calib_->getRowPattern();
0209 std::vector<short> temprowvals(0);
0210 std::vector<short> tempcolvals(0);
0211 uint32_t nminuscol = 0;
0212 uint32_t nminusrow = 0;
0213 uint32_t npatterns = 0;
0214 for (uint32_t icol = 0; icol < calibcols.size(); icol++) {
0215 if (calibcols[icol] == -1) {
0216 nminuscol++;
0217 } else if (nminuscol == colpatternnumber) {
0218
0219 short val = calibcols[icol];
0220 tempcolvals.push_back(val);
0221 } else if (nminuscol > colpatternnumber)
0222 break;
0223 }
0224 for (uint32_t irow = 0; irow < calibrows.size(); irow++) {
0225
0226 if (calibrows[irow] == -1)
0227 nminusrow++;
0228 else if (nminusrow == rowpatternnumber) {
0229 short val = calibrows[irow];
0230 temprowvals.push_back(val);
0231 } else if (nminusrow > rowpatternnumber)
0232 break;
0233 }
0234
0235 while (currentpattern_.size() > temprowvals.size() * tempcolvals.size()) {
0236 currentpattern_.erase(currentpattern_.end());
0237 }
0238 for (uint32_t irow = 0; irow < temprowvals.size(); irow++) {
0239 for (uint32_t icol = 0; icol < tempcolvals.size(); icol++) {
0240 std::pair<short, short> pattern(temprowvals[irow], tempcolvals[icol]);
0241 npatterns++;
0242 if (npatterns > currentpattern_.size())
0243 currentpattern_.push_back(pattern);
0244 else
0245 currentpattern_[npatterns - 1] = pattern;
0246 }
0247 }
0248 }
0249
0250
0251
0252
0253 void SiPixelCalibDigiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0254
0255 using namespace edm;
0256 calib_ = iSetup.getHandle(calibToken_);
0257 theGeometry_ = iSetup.getHandle(trackerGeomToken_);
0258 theCablingMap_ = iSetup.getHandle(cablingMapToken_);
0259 pattern_repeat_ = calib_->getNTriggers() * calib_->nVCal();
0260 if (use_realeventnumber_) {
0261 iEventCounter_ = iEvent.id().event() - 1;
0262 } else
0263 iEventCounter_++;
0264 if (iEventCounter_ % pattern_repeat_ == 1)
0265 setPattern();
0266
0267
0268 fill(iEvent, iSetup);
0269
0270 auto pOut = std::make_unique<edm::DetSetVector<SiPixelCalibDigi>>();
0271 auto pErr = std::make_unique<edm::DetSetVector<SiPixelCalibDigiError>>();
0272
0273
0274 if (store()) {
0275
0276 for (std::map<pixelstruct, SiPixelCalibDigi>::const_iterator idet = intermediate_data_.begin();
0277 idet != intermediate_data_.end();
0278 ++idet) {
0279 uint32_t detid = idet->first.first;
0280 if (!control_pattern_size_) {
0281 if (!checkPixel(idet->first.first, idet->first.second.first, idet->first.second.second))
0282 continue;
0283 }
0284
0285 SiPixelCalibDigi tempdigi = idet->second;
0286 edm::DetSet<SiPixelCalibDigi>& detSet = pOut->find_or_insert(detid);
0287 detSet.data.push_back(tempdigi);
0288 }
0289 if (includeErrors_) {
0290 for (std::map<pixelstruct, SiPixelCalibDigiError>::const_iterator ierr = error_data_.begin();
0291 ierr != error_data_.end();
0292 ++ierr) {
0293 uint32_t detid = ierr->first.first;
0294 SiPixelCalibDigiError temperror = ierr->second;
0295 edm::DetSet<SiPixelCalibDigiError>& errSet = pErr->find_or_insert(detid);
0296 errSet.data.push_back(temperror);
0297 }
0298 }
0299 edm::LogInfo("INFO") << "now filling event " << iEventCounter_ << " as pixel pattern changes every "
0300 << pattern_repeat_ << " events..." << std::endl;
0301 clear();
0302 }
0303 iEvent.put(std::move(pOut));
0304 if (includeErrors_)
0305 iEvent.put(std::move(pErr));
0306 }
0307
0308
0309 bool SiPixelCalibDigiProducer::checkPixel(uint32_t detid, short row, short col) {
0310 if (!control_pattern_size_ && !store())
0311 return true;
0312
0313 if (!ignore_non_pattern_)
0314 return true;
0315
0316 edm::LogInfo("SiPixelCalibDigiProducer") << "Event" << iEventCounter_ << ",now in checkpixel() " << std::endl;
0317 if (currentpattern_.empty())
0318 setPattern();
0319
0320 uint32_t fedid = detid_to_fedid_[detid];
0321
0322 SiPixelFrameConverter formatter(theCablingMap_.product(), fedid);
0323 sipixelobjects::DetectorIndex detector = {detid, row, col};
0324 sipixelobjects::ElectronicIndex cabling;
0325
0326 formatter.toCabling(cabling, detector);
0327
0328
0329
0330 sipixelobjects::LocalPixel::DcolPxid loc;
0331 loc.dcol = cabling.dcol;
0332 loc.pxid = cabling.pxid;
0333 sipixelobjects::LocalPixel locpixel(loc);
0334 currentpair_.first = locpixel.rocRow();
0335 currentpair_.second = locpixel.rocCol();
0336
0337 for (uint32_t i = 0; i < currentpattern_.size(); ++i) {
0338
0339 if (currentpair_ == currentpattern_[i]) {
0340 return true;
0341 }
0342 }
0343 std::ostringstream errorlog;
0344 errorlog << "DETID " << detid << ", row, col (offline)=" << row << "," << col
0345 << " row, col (ROC) =" << currentpair_.first << "," << currentpair_.second
0346 << " found no match in list of patterns: ";
0347 for (uint32_t i = 0; i < currentpattern_.size(); ++i) {
0348 if (i != 0 && i != currentpattern_.size() - 1)
0349 errorlog << " ";
0350 errorlog << "(";
0351 errorlog << currentpattern_[i].first;
0352 errorlog << ",";
0353 errorlog << currentpattern_[i].second;
0354 errorlog << ")";
0355 }
0356 edm::LogError("ERROR") << errorlog.str() << std::endl;
0357 if (includeErrors_) {
0358
0359 pixelstruct temppixelworker;
0360 temppixelworker.first = detid;
0361 temppixelworker.second.first = row;
0362 temppixelworker.second.second = col;
0363 std::map<pixelstruct, SiPixelCalibDigiError>::const_iterator ierr = error_data_.find(temppixelworker);
0364 if (ierr == error_data_.end()) {
0365 SiPixelCalibDigiError temperr(row, col, 1);
0366 error_data_[temppixelworker] = temperr;
0367 }
0368 }
0369
0370 return false;
0371 }
0372
0373 DEFINE_FWK_MODULE(SiPixelCalibDigiProducer);