File indexing completed on 2025-06-12 23:30:06
0001 #include <iostream>
0002 #include "FWCore/Utilities/interface/Exception.h"
0003
0004 #include "RecoPPS/Local/interface/RPixDetClusterizer.h"
0005 #include "DataFormats/CTPPSDetId/interface/CTPPSPixelDetId.h"
0006
0007 namespace {
0008 constexpr int maxCol = CTPPSPixelCluster::MAXCOL;
0009 constexpr int maxRow = CTPPSPixelCluster::MAXROW;
0010 constexpr double highRangeCal = 1800.;
0011 constexpr double lowRangeCal = 260.;
0012 constexpr int rocMask = 0xE000;
0013 constexpr int rocOffset = 13;
0014 }
0015
0016 RPixDetClusterizer::RPixDetClusterizer(edm::ParameterSet const &conf) : SeedVector_(0) {
0017 verbosity_ = conf.getUntrackedParameter<int>("RPixVerbosity");
0018 SeedADCThreshold_ = conf.getParameter<int>("SeedADCThreshold");
0019 ADCThreshold_ = conf.getParameter<int>("ADCThreshold");
0020 ElectronADCGain_ = conf.getParameter<double>("ElectronADCGain");
0021 VcaltoElectronGain_ = conf.getParameter<int>("VCaltoElectronGain");
0022 VcaltoElectronOffset_ = conf.getParameter<int>("VCaltoElectronOffset");
0023 doSingleCalibration_ = conf.getParameter<bool>("doSingleCalibration");
0024 }
0025
0026 RPixDetClusterizer::~RPixDetClusterizer() {}
0027
0028 void RPixDetClusterizer::buildClusters(unsigned int detId,
0029 const std::vector<CTPPSPixelDigi> &digi,
0030 std::vector<CTPPSPixelCluster> &clusters,
0031 const CTPPSPixelGainCalibrations *pcalibrations,
0032 const CTPPSPixelAnalysisMask *maskera) {
0033 std::map<uint32_t, CTPPSPixelROCAnalysisMask> const &mask = maskera->analysisMask;
0034 std::set<std::pair<unsigned char, unsigned char> > maskedPixels;
0035
0036 bool fullyMaskedRoc[6][6] = {{false}};
0037 CTPPSPixelDetId currentId(detId);
0038
0039
0040 CTPPSPixelIndices pI;
0041 for (auto const &det : mask) {
0042 uint32_t planeId = det.first & ~rocMask;
0043
0044 if (planeId == detId) {
0045 unsigned int rocNum = (det.first & rocMask) >> rocOffset;
0046 if (rocNum > 5)
0047 throw cms::Exception("InvalidRocNumber") << "roc number from mask: " << rocNum;
0048 if (det.second.fullMask)
0049 fullyMaskedRoc[currentId.plane()][rocNum] = true;
0050 for (auto const &paio : det.second.maskedPixels) {
0051 std::pair<unsigned char, unsigned char> pa = paio;
0052 int modCol, modRow;
0053 int status = pI.transformToModule(pa.second, pa.first, rocNum, modCol, modRow);
0054
0055 if (status == 0) {
0056 std::pair<int, int> modPix(modRow, modCol);
0057 maskedPixels.insert(modPix);
0058 }
0059 }
0060 }
0061 }
0062
0063 if (verbosity_)
0064 edm::LogInfo("RPixDetClusterizer") << detId << " received digi.size()=" << digi.size();
0065
0066
0067
0068 rpix_digi_set_.clear();
0069 rpix_digi_set_.insert(digi.begin(), digi.end());
0070 SeedVector_.clear();
0071
0072
0073 calib_rpix_digi_map_.clear();
0074
0075 for (auto const &RPdit : rpix_digi_set_) {
0076 uint8_t row = RPdit.row();
0077 uint8_t column = RPdit.column();
0078 if (row > maxRow || column > maxCol)
0079 throw cms::Exception("CTPPSDigiOutofRange") << " row = " << row << " column = " << column;
0080
0081 std::pair<unsigned char, unsigned char> pixel = std::make_pair(row, column);
0082 unsigned short adc = RPdit.adc();
0083 unsigned short electrons = 0;
0084
0085
0086 const bool is_in = maskedPixels.find(pixel) != maskedPixels.end();
0087
0088
0089 int piano = currentId.plane();
0090 int rocId = pI.getROCId(column, row);
0091
0092 if (!is_in && !fullyMaskedRoc[piano][rocId]) {
0093
0094 electrons = calibrate(detId, adc, row, column, pcalibrations);
0095 if (electrons < SeedADCThreshold_ * ElectronADCGain_)
0096 electrons = SeedADCThreshold_ * ElectronADCGain_;
0097 RPixCalibDigi calibDigi(row, column, adc, electrons);
0098 unsigned int index = column * maxRow + row;
0099 calib_rpix_digi_map_.insert(std::pair<unsigned int, RPixCalibDigi>(index, calibDigi));
0100 SeedVector_.push_back(calibDigi);
0101 }
0102 }
0103 if (verbosity_)
0104 edm::LogInfo("RPixDetClusterizer") << " RPix set size = " << calib_rpix_digi_map_.size();
0105
0106 for (auto SeedIt : SeedVector_) {
0107 make_cluster(SeedIt, clusters);
0108 }
0109 }
0110
0111 void RPixDetClusterizer::make_cluster(RPixCalibDigi const &aSeed, std::vector<CTPPSPixelCluster> &clusters) {
0112
0113 unsigned int seedIndex = aSeed.column() * maxRow + aSeed.row();
0114 if (calib_rpix_digi_map_.find(seedIndex) == calib_rpix_digi_map_.end()) {
0115 return;
0116 }
0117
0118 RPixTempCluster atempCluster;
0119
0120
0121 atempCluster.addPixel(aSeed.row(), aSeed.column(), aSeed.electrons());
0122 calib_rpix_digi_map_.erase(seedIndex);
0123
0124 while (!atempCluster.empty()) {
0125
0126 auto curInd = atempCluster.top();
0127 atempCluster.pop();
0128 for (auto c = std::max(0, int(atempCluster.col[curInd]) - 1);
0129 c < std::min(int(atempCluster.col[curInd]) + 2, maxCol);
0130 ++c) {
0131 for (auto r = std::max(0, int(atempCluster.row[curInd]) - 1);
0132 r < std::min(int(atempCluster.row[curInd]) + 2, maxRow);
0133 ++r) {
0134 unsigned int currIndex = c * maxRow + r;
0135 if (calib_rpix_digi_map_.find(currIndex) != calib_rpix_digi_map_.end()) {
0136 if (!atempCluster.addPixel(r, c, calib_rpix_digi_map_[currIndex].electrons())) {
0137 clusters.emplace_back(atempCluster.isize, atempCluster.adc, atempCluster.row, atempCluster.col);
0138 return;
0139 }
0140 calib_rpix_digi_map_.erase(currIndex);
0141 }
0142 }
0143 }
0144
0145 }
0146
0147 clusters.emplace_back(atempCluster.isize, atempCluster.adc, atempCluster.row, atempCluster.col);
0148 }
0149
0150 int RPixDetClusterizer::calibrate(
0151 unsigned int detId, int adc, int row, int col, const CTPPSPixelGainCalibrations *pcalibrations) {
0152 float gain = 0;
0153 float pedestal = 0;
0154 int electrons = 0;
0155
0156 if (!doSingleCalibration_) {
0157 const CTPPSPixelGainCalibration &DetCalibs = pcalibrations->getGainCalibration(detId);
0158
0159 if (DetCalibs.getDetId() != 0) {
0160 gain = DetCalibs.getGain(col, row) * highRangeCal / lowRangeCal;
0161 pedestal = DetCalibs.getPed(col, row);
0162 float vcal = (adc - pedestal) * gain;
0163 electrons = int(vcal * VcaltoElectronGain_ + VcaltoElectronOffset_);
0164 } else {
0165 gain = ElectronADCGain_;
0166 pedestal = 0;
0167 electrons = int(adc * gain - pedestal);
0168 }
0169
0170 if (verbosity_)
0171 edm::LogInfo("RPixCalibration") << "calibrate: adc = " << adc << " electrons = " << electrons
0172 << " gain = " << gain << " pedestal = " << pedestal;
0173 } else {
0174 gain = ElectronADCGain_;
0175 pedestal = 0;
0176 electrons = int(adc * gain - pedestal);
0177 if (verbosity_)
0178 edm::LogInfo("RPixCalibration") << "calibrate: adc = " << adc << " electrons = " << electrons
0179 << " ElectronADCGain = " << gain << " pedestal = " << pedestal;
0180 }
0181 if (electrons < 0) {
0182 edm::LogInfo("RPixCalibration") << "RPixDetClusterizer::calibrate: *** electrons < 0 *** --> " << electrons
0183 << " --> electrons = 0";
0184 electrons = 0;
0185 }
0186
0187 return electrons;
0188 }