Local/src/RPixDetClusterizer.cc

0001 #include <iostream>
0002 #include "FWCore/Utilities/interface/Exception.h"
0003
0004 #include "RecoPPS/Local/interface/RPixDetClusterizer.h"
0005
0006 namespace {
0007   constexpr int maxCol = CTPPSPixelCluster::MAXCOL;
0008   constexpr int maxRow = CTPPSPixelCluster::MAXROW;
0009   constexpr double highRangeCal = 1800.;
0010   constexpr double lowRangeCal = 260.;
0011 }  // namespace
0012
0013 RPixDetClusterizer::RPixDetClusterizer(edm::ParameterSet const &conf) : params_(conf), SeedVector_(0) {
0014   verbosity_ = conf.getUntrackedParameter<int>("RPixVerbosity");
0015   SeedADCThreshold_ = conf.getParameter<int>("SeedADCThreshold");
0016   ADCThreshold_ = conf.getParameter<int>("ADCThreshold");
0017   ElectronADCGain_ = conf.getParameter<double>("ElectronADCGain");
0018   VcaltoElectronGain_ = conf.getParameter<int>("VCaltoElectronGain");
0019   VcaltoElectronOffset_ = conf.getParameter<int>("VCaltoElectronOffset");
0020   doSingleCalibration_ = conf.getParameter<bool>("doSingleCalibration");
0021 }
0022
0023 RPixDetClusterizer::~RPixDetClusterizer() {}
0024
0025 void RPixDetClusterizer::buildClusters(unsigned int detId,
0026                                        const std::vector<CTPPSPixelDigi> &digi,
0027                                        std::vector<CTPPSPixelCluster> &clusters,
0028                                        const CTPPSPixelGainCalibrations *pcalibrations,
0029                                        const CTPPSPixelAnalysisMask *maskera) {
0030   std::map<uint32_t, CTPPSPixelROCAnalysisMask> const &mask = maskera->analysisMask;
0031   std::map<uint32_t, CTPPSPixelROCAnalysisMask>::const_iterator mask_it = mask.find(detId);
0032
0033   std::set<std::pair<unsigned char, unsigned char> > maskedPixels;
0034   if (mask_it != mask.end())
0035     maskedPixels = mask_it->second.maskedPixels;
0036
0037   if (verbosity_)
0038     edm::LogInfo("RPixDetClusterizer") << detId << " received digi.size()=" << digi.size();
0039
0040   // creating a set of CTPPSPixelDigi's and fill it
0041
0042   rpix_digi_set_.clear();
0043   rpix_digi_set_.insert(digi.begin(), digi.end());
0044   SeedVector_.clear();
0045
0046   // calibrate digis here
0047   calib_rpix_digi_map_.clear();
0048
0049   for (auto const &RPdit : rpix_digi_set_) {
0050     uint8_t row = RPdit.row();
0051     uint8_t column = RPdit.column();
0052     if (row > maxRow || column > maxCol)
0053       throw cms::Exception("CTPPSDigiOutofRange") << " row = " << row << "  column = " << column;
0054
0055     std::pair<unsigned char, unsigned char> pixel = std::make_pair(row, column);
0056     unsigned short adc = RPdit.adc();
0057     unsigned short electrons = 0;
0058
0059     // check if pixel is noisy or dead (i.e. in the mask)
0060     const bool is_in = maskedPixels.find(pixel) != maskedPixels.end();
0061     if (!is_in) {
0062       //calibrate digi and store the new ones
0063       electrons = calibrate(detId, adc, row, column, pcalibrations);
0064       if (electrons < SeedADCThreshold_ * ElectronADCGain_)
0065         electrons = SeedADCThreshold_ * ElectronADCGain_;
0066       RPixCalibDigi calibDigi(row, column, adc, electrons);
0067       unsigned int index = column * maxRow + row;
0068       calib_rpix_digi_map_.insert(std::pair<unsigned int, RPixCalibDigi>(index, calibDigi));
0069       SeedVector_.push_back(calibDigi);
0070     }
0071   }
0072   if (verbosity_)
0073     edm::LogInfo("RPixDetClusterizer") << " RPix set size = " << calib_rpix_digi_map_.size();
0074
0075   for (auto SeedIt : SeedVector_) {
0076     make_cluster(SeedIt, clusters);
0077   }
0078 }
0079
0080 void RPixDetClusterizer::make_cluster(RPixCalibDigi const &aSeed, std::vector<CTPPSPixelCluster> &clusters) {
0081   /// check if seed already used
0082   unsigned int seedIndex = aSeed.column() * maxRow + aSeed.row();
0083   if (calib_rpix_digi_map_.find(seedIndex) == calib_rpix_digi_map_.end()) {
0084     return;
0085   }
0086   // creating a temporary cluster
0087   RPixTempCluster atempCluster;
0088
0089   // filling the cluster with the seed
0090   atempCluster.addPixel(aSeed.row(), aSeed.column(), aSeed.electrons());
0091   calib_rpix_digi_map_.erase(seedIndex);
0092
0093   while (!atempCluster.empty()) {
0094     //This is the standard algorithm to find and add a pixel
0095     auto curInd = atempCluster.top();
0096     atempCluster.pop();
0097     for (auto c = std::max(0, int(atempCluster.col[curInd]) - 1);
0098          c < std::min(int(atempCluster.col[curInd]) + 2, maxCol);
0099          ++c) {
0100       for (auto r = std::max(0, int(atempCluster.row[curInd]) - 1);
0101            r < std::min(int(atempCluster.row[curInd]) + 2, maxRow);
0102            ++r) {
0103         unsigned int currIndex = c * maxRow + r;
0104         if (calib_rpix_digi_map_.find(currIndex) != calib_rpix_digi_map_.end()) {
0105           if (!atempCluster.addPixel(r, c, calib_rpix_digi_map_[currIndex].electrons())) {
0106             clusters.emplace_back(atempCluster.isize, atempCluster.adc, atempCluster.row, atempCluster.col);
0107             return;
0108           }
0109           calib_rpix_digi_map_.erase(currIndex);
0110         }
0111       }
0112     }
0113
0114   }  // while accretion
0115
0116   clusters.emplace_back(atempCluster.isize, atempCluster.adc, atempCluster.row, atempCluster.col);
0117 }
0118
0119 int RPixDetClusterizer::calibrate(
0120     unsigned int detId, int adc, int row, int col, const CTPPSPixelGainCalibrations *pcalibrations) {
0121   float gain = 0;
0122   float pedestal = 0;
0123   int electrons = 0;
0124
0125   if (!doSingleCalibration_) {
0126     const CTPPSPixelGainCalibration &DetCalibs = pcalibrations->getGainCalibration(detId);
0127
0128     if (DetCalibs.getDetId() != 0) {
0129       gain = DetCalibs.getGain(col, row) * highRangeCal / lowRangeCal;
0130       pedestal = DetCalibs.getPed(col, row);
0131       float vcal = (adc - pedestal) * gain;
0132       electrons = int(vcal * VcaltoElectronGain_ + VcaltoElectronOffset_);
0133     } else {
0134       gain = ElectronADCGain_;
0135       pedestal = 0;
0136       electrons = int(adc * gain - pedestal);
0137     }
0138
0139     if (verbosity_)
0140       edm::LogInfo("RPixCalibration") << "calibrate:  adc = " << adc << " electrons = " << electrons
0141                                       << " gain = " << gain << " pedestal = " << pedestal;
0142   } else {
0143     gain = ElectronADCGain_;
0144     pedestal = 0;
0145     electrons = int(adc * gain - pedestal);
0146     if (verbosity_)
0147       edm::LogInfo("RPixCalibration") << "calibrate:  adc = " << adc << " electrons = " << electrons
0148                                       << " ElectronADCGain = " << gain << " pedestal = " << pedestal;
0149   }
0150   if (electrons < 0) {
0151     edm::LogInfo("RPixCalibration") << "RPixDetClusterizer::calibrate: *** electrons < 0 *** --> " << electrons
0152                                     << "  --> electrons = 0";
0153     electrons = 0;
0154   }
0155
0156   return electrons;
0157 }