CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SiPixelStatusProducer.cc
Go to the documentation of this file.
1 
7 // C++ standard
8 #include <string>
9 // ROOT
10 #include "TMath.h"
11 // CMSSW FW
20 // CMSSW DataFormats
28 // "FED error 25"
30 // CMSSW CondFormats
39 // EDProducer related dataformat
42 // header file
44 
45 using namespace std;
46 
47 //--------------------------------------------------------------------------------------------------
49  // get parameter
50 
51  // badPixelFEDChannelCollections
52  std::vector<edm::InputTag> badPixelFEDChannelCollectionLabels_ =
53  iConfig.getParameter<edm::ParameterSet>("SiPixelStatusProducerParameters")
54  .getParameter<std::vector<edm::InputTag>>("badPixelFEDChannelCollections");
55  for (auto& t : badPixelFEDChannelCollectionLabels_)
56  theBadPixelFEDChannelsTokens_.push_back(consumes<PixelFEDChannelCollection>(t));
57  // badPixelFEDChannelCollections = cms.VInputTag(cms.InputTag('siPixelDigis'))
58 
59  fPixelClusterLabel_ = iConfig.getParameter<edm::ParameterSet>("SiPixelStatusProducerParameters")
60  .getUntrackedParameter<edm::InputTag>("pixelClusterLabel");
61  fSiPixelClusterToken_ = consumes<edmNew::DetSetVector<SiPixelCluster>>(fPixelClusterLabel_);
62  resetNLumi_ = iConfig.getParameter<edm::ParameterSet>("SiPixelStatusProducerParameters")
63  .getUntrackedParameter<int>("resetEveryNLumi", 1);
64 
65  ftotalevents = 0;
66  countLumi_ = 0;
67 
68  beginLumi_ = endLumi_ = -1;
69  endLumi_ = endRun_ = -1;
70 
71  produces<SiPixelDetectorStatus, edm::Transition::EndLuminosityBlock>("siPixelStatus");
72 }
73 
74 //--------------------------------------------------------------------------------------------------
76 
77 //--------------------------------------------------------------------------------------------------
79  edm::LogInfo("SiPixelStatusProducer") << "beginlumi setup " << endl;
80 
81  if (countLumi_ == 0 && resetNLumi_ > 0) {
82  beginLumi_ = lumiSeg.luminosityBlock();
83  beginRun_ = lumiSeg.run();
84  ftotalevents = 0;
85  }
86 
87  // The es watcher is acutally not needed if run parallel jobs for each lumi section
88  if (siPixelFedCablingMapWatcher_.check(iSetup) || trackerDIGIGeoWatcher_.check(iSetup) ||
89  trackerTopoWatcher_.check(iSetup)) {
90  coord_.init(iSetup);
91  iSetup.get<SiPixelFedCablingMapRcd>().get(fCablingMap);
92  fCablingMap_ = fCablingMap.product();
93  iSetup.get<TrackerDigiGeometryRecord>().get(fTG);
94 
95  fFedIds = fCablingMap_->det2fedMap();
96 
97  } // if conditionWatcher_.check(iSetup)
98 
99  // init the SiPixelDetectorStatus fDet and sensor size fSensors in the begining (when countLumi is zero)
100  if (countLumi_ == 0) {
101  for (TrackerGeometry::DetContainer::const_iterator it = fTG->dets().begin(); it != fTG->dets().end(); it++) {
102  const PixelGeomDetUnit* pgdu = dynamic_cast<const PixelGeomDetUnit*>((*it));
103  if (pgdu == nullptr)
104  continue;
105  DetId detId = (*it)->geographicalId();
106  int detid = detId.rawId();
107 
108  // don't want to use magic number row 80 column 52
109  const PixelTopology* topo = static_cast<const PixelTopology*>(&pgdu->specificTopology());
110  int rowsperroc = topo->rowsperroc();
111  int colsperroc = topo->colsperroc();
112 
113  int nROCrows = pgdu->specificTopology().nrows() / rowsperroc;
114  int nROCcolumns = pgdu->specificTopology().ncolumns() / colsperroc;
115  int nrocs = nROCrows * nROCcolumns;
116 
117  fDet.addModule(detid, nrocs);
118 
119  fSensors[detid] = std::make_pair(rowsperroc, colsperroc);
120  fSensorLayout[detid] = std::make_pair(nROCrows, nROCcolumns);
121 
122  std::map<int, int> rocIdMap;
123  for (int irow = 0; irow < nROCrows; irow++) {
124  for (int icol = 0; icol < nROCcolumns; icol++) {
125  int dummyOfflineRow = (rowsperroc / 2 - 1) + irow * rowsperroc;
126  int dummeOfflineColumn = (colsperroc / 2 - 1) + icol * colsperroc;
127  // encode irow, icol
128  int key = indexROC(irow, icol, nROCcolumns);
129 
130  int roc(-1), rocR(-1), rocC(-1);
131  onlineRocColRow(detId, dummyOfflineRow, dummeOfflineColumn, roc, rocR, rocC);
132 
133  int value = roc;
134  rocIdMap[key] = value;
135  }
136  }
137 
138  fRocIds[detid] = rocIdMap;
139  }
140 
141  } // init when countLumi = 0
142 
143  FEDerror25_.clear();
144  countLumi_++;
145 }
146 
147 //--------------------------------------------------------------------------------------------------
149  ftotalevents++;
150 
151  edm::LogInfo("SiPixelStatusProducer") << "start cluster analyzer " << endl;
152 
153  // ----------------------------------------------------------------------
154  // -- Pixel cluster analysis
155  // ----------------------------------------------------------------------
156 
158  if (!iEvent.getByToken(fSiPixelClusterToken_, hClusterColl)) {
159  edm::LogWarning("SiPixelStatusProducer")
160  << " edmNew::DetSetVector<SiPixelCluster> " << fPixelClusterLabel_ << " does not exist!" << endl;
161  return;
162  }
163 
164  iEvent.getByToken(fSiPixelClusterToken_, hClusterColl);
165 
166  if (hClusterColl.isValid()) {
167  for (const auto& clusters : *hClusterColl) { //loop over different clusters in a clusters vector (module)
168  for (const auto& clu : clusters) { // loop over cluster in a given detId (module)
169  int detid = clusters.detId();
170  int rowsperroc = fSensors[detid].first;
171  int colsperroc = fSensors[detid].second;
172 
173  int nROCcolumns = fSensorLayout[detid].second;
174 
175  int roc(-1);
176  std::map<int, int> fRocIds_detid;
177  if (fRocIds.find(detid) != fRocIds.end()) {
178  fRocIds_detid = fRocIds[detid];
179  }
180 
181  const vector<SiPixelCluster::Pixel>& pixvector = clu.pixels();
182  for (unsigned int i = 0; i < pixvector.size(); ++i) {
183  int mr0 = pixvector[i].x; // constant column direction is along x-axis,
184  int mc0 = pixvector[i].y; // constant row direction is along y-axis
185 
186  int irow = mr0 / rowsperroc;
187  int icol = mc0 / colsperroc;
188 
189  int key = indexROC(irow, icol, nROCcolumns);
190  if (fRocIds_detid.find(key) != fRocIds_detid.end()) {
191  roc = fRocIds_detid[key];
192  }
193 
194  fDet.fillDIGI(detid, roc);
195 
196  } // loop over pixels in a given cluster
197 
198  } // loop over cluster in a given detId (module)
199 
200  } // loop over detId-grouped clusters in cluster detId-grouped clusters-vector
201 
202  } // hClusterColl.isValid()
203  else {
204  edm::LogWarning("SiPixelStatusProducer")
205  << " edmNew::DetSetVector<SiPixelCluster> " << fPixelClusterLabel_ << " is NOT Valid!" << endl;
206  }
208 
209  // FEDerror25 per-ROC per-event
210  edm::Handle<PixelFEDChannelCollection> pixelFEDChannelCollectionHandle;
211 
212  // look over different resouces of takens
213  for (const edm::EDGetTokenT<PixelFEDChannelCollection>& tk : theBadPixelFEDChannelsTokens_) {
214  // collection has to exist
215  if (!iEvent.getByToken(tk, pixelFEDChannelCollectionHandle)) {
216  edm::LogWarning("SiPixelStatusProducer")
217  << " PixelFEDChannelCollection with index " << tk.index() << " does NOT exist!" << std::endl;
218  continue;
219  }
220  iEvent.getByToken(tk, pixelFEDChannelCollectionHandle);
221  // collection has to be valid
222  if (!pixelFEDChannelCollectionHandle.isValid()) {
223  edm::LogWarning("SiPixelStatusProducer")
224  << " PixelFEDChannelCollection with index " << tk.index() << " is NOT valid!" << endl;
225  continue;
226  }
227  // FEDerror channels for the current events
228  std::map<int, std::vector<PixelFEDChannel>> tmpFEDerror25;
229  for (const auto& disabledChannels : *pixelFEDChannelCollectionHandle) {
230  //loop over different PixelFED in a PixelFED vector (module)
231  for (const auto& ch : disabledChannels) {
232  DetId detId = disabledChannels.detId();
233  int detid = detId.rawId();
234 
235  if (ftotalevents == 1) {
236  // FEDerror25 channels for the first event in the lumi section
237  FEDerror25_[detid].push_back(ch);
238  } else
239  tmpFEDerror25[detid].push_back(ch);
240 
241  } // loop over different PixelFED in a PixelFED vector (different channel for a given module)
242 
243  } // loop over different (different DetId) PixelFED vectors in PixelFEDChannelCollection
244 
245  // Compare the current FEDerror list with the first event's FEDerror list
246  // and save the common channels
247  if (!tmpFEDerror25.empty() && !FEDerror25_.empty()) { // non-empty FEDerror lists
248 
249  std::map<int, std::vector<PixelFEDChannel>>::iterator itFEDerror25;
250  for (itFEDerror25 = FEDerror25_.begin(); itFEDerror25 != FEDerror25_.end(); itFEDerror25++) {
251  int detid = itFEDerror25->first;
252  if (tmpFEDerror25.find(detid) != tmpFEDerror25.end()) {
253  std::vector<PixelFEDChannel> chs = itFEDerror25->second;
254  std::vector<PixelFEDChannel> chs_tmp = tmpFEDerror25[detid];
255 
256  std::vector<PixelFEDChannel> chs_common;
257  for (unsigned int ich = 0; ich < chs.size(); ich++) {
258  PixelFEDChannel ch = chs[ich];
259  // look over the current FEDerror25 channels, save the common FED channels
260  for (unsigned int ich_tmp = 0; ich_tmp < chs_tmp.size(); ich_tmp++) {
261  PixelFEDChannel ch_tmp = chs_tmp[ich_tmp];
262  if ((ch.fed == ch_tmp.fed) && (ch.link == ch_tmp.link)) { // the same FED channel
263  chs_common.push_back(ch);
264  break;
265  }
266  }
267  }
268  // remove the full module from FEDerror25 list if no common channels are left
269  if (chs_common.empty())
270  FEDerror25_.erase(itFEDerror25);
271  // otherwise replace with the common channels
272  else {
273  FEDerror25_[detid].clear();
274  FEDerror25_[detid] = chs_common;
275  }
276  } else { // remove the full module from FEDerror25 list if the module doesn't appear in the current event's FEDerror25 list
277  FEDerror25_.erase(itFEDerror25);
278  }
279 
280  } // loop over modules that have FEDerror25 in the first event in the lumi section
281 
282  } // non-empty FEDerror lists
283 
284  } // look over different resouces of takens
285 
286  // no per-event collection put into iEvent
287  // If use produce() function and no collection is put into iEvent, produce() will not run in unScheduled mode
288  // Now since CMSSW_10_1_X, the accumulate() function will run whatsoever in the unScheduled mode
289  // But accumulate() is NOT available and will NOT be available for releases before CMSSW_10_1_X
290 }
291 
292 //--------------------------------------------------------------------------------------------------
294 
295 //--------------------------------------------------------------------------------------------------
297  edm::LogInfo("SiPixelStatusProducer") << "endlumi producer " << endl;
298 
299  endLumi_ = lumiSeg.luminosityBlock();
300  endRun_ = lumiSeg.run();
301 
302  // check if countLumi_ is large enough to read out/save data and reset for the next round
303  if (resetNLumi_ == -1)
304  return;
305  if (countLumi_ < resetNLumi_)
306  return;
307 
308  // set the FEDerror25 flag to be true for ROCs send out FEDerror25 for all events in the lumi section
309  if (!FEDerror25_.empty()) {
310  std::map<int, std::vector<PixelFEDChannel>>::iterator itFEDerror25;
311  for (itFEDerror25 = FEDerror25_.begin(); itFEDerror25 != FEDerror25_.end(); itFEDerror25++) {
312  int detid = itFEDerror25->first;
313  std::vector<PixelFEDChannel> chs = itFEDerror25->second;
314  for (unsigned int ich = 0; ich < chs.size(); ich++) {
315  fDet.fillFEDerror25(detid, chs[ich]);
316  }
317  }
318  }
319 
320  fDet.setRunRange(beginRun_, endRun_);
321  fDet.setLSRange(beginLumi_, endLumi_);
322  fDet.setNevents(ftotalevents);
323 
324  // save result
325  auto result = std::make_unique<SiPixelDetectorStatus>();
326  *result = fDet;
327 
328  // only save for the lumi sections with NON-ZERO events
329  lumiSeg.put(std::move(result), std::string("siPixelStatus"));
330  edm::LogInfo("SiPixelStatusProducer") << "new lumi-based data stored for run " << beginRun_ << " lumi from "
331  << beginLumi_ << " to " << endLumi_ << std::endl;
332 
333  // reset detector status and lumi-counter
334  fDet.resetDetectorStatus();
335  countLumi_ = 0;
336  ftotalevents = 0;
337  FEDerror25_.clear();
338 }
339 
340 //--------------------------------------------------------------------------------------------------
342  const DetId& detId, int offlineRow, int offlineCol, int& roc, int& row, int& col) {
343  int fedId = fFedIds[detId.rawId()];
344 
345  // from detector to cabling
347  sipixelobjects::DetectorIndex detector; //{detId.rawId(), offlineRow, offlineCol};
348  detector.rawId = detId.rawId();
349  detector.row = offlineRow;
350  detector.col = offlineCol;
351 
352  SiPixelFrameConverter converter(fCablingMap_, fedId);
353  converter.toCabling(cabling, detector);
354 
355  // then one can construct local pixel
357  loc.dcol = cabling.dcol;
358  loc.pxid = cabling.pxid;
359  // and get local(online) row/column
360  sipixelobjects::LocalPixel locpixel(loc);
361  col = locpixel.rocCol();
362  row = locpixel.rocRow();
363  //sipixelobjects::CablingPathToDetUnit path = {(unsigned int) fedId, (unsigned int)cabling.link, (unsigned int)cabling.roc};
364  //const sipixelobjects::PixelROC *theRoc = fCablingMap->findItem(path);
365  const sipixelobjects::PixelROC* theRoc = converter.toRoc(cabling.link, cabling.roc);
366  roc = theRoc->idInDetUnit();
367 
368  // has to be BPIX; has to be minus side; has to be half module
369  // for phase-I, there is no half module
370  if (detId.subdetId() == PixelSubdetector::PixelBarrel && coord_.side(detId) == 1 && coord_.half(detId)) {
371  roc += 8;
372  }
373 }
374 
375 //--------------------------------------------------------------------------------------------------
376 int SiPixelStatusProducer::indexROC(int irow, int icol, int nROCcolumns) {
377  return int(icol + irow * nROCcolumns);
378 
379  // generate the folling roc index that is going to map with ROC id as
380  // 8 9 10 11 12 13 14 15
381  // 0 1 2 3 4 5 6 7
382 }
383 
384 //--------------------------------------------------------------------------------------------------
385 //edmPythonConfigToCppValidation CalibTracker/SiPixelQuality/python/SiPixelStatusProducer_cfi.py
387  // siPixelStatusProducer
389  {
391  psd0.addUntracked<int>("resetEveryNLumi", 1);
392  psd0.addUntracked<edm::InputTag>("pixelClusterLabel", edm::InputTag("siPixelClusters", "", "RECO"));
393  psd0.add<std::vector<edm::InputTag>>("badPixelFEDChannelCollections",
394  {
395  edm::InputTag("siPixelDigis"),
396  });
397  desc.add<edm::ParameterSetDescription>("SiPixelStatusProducerParameters", psd0);
398  }
399  descriptions.add("siPixelStatusProducer", desc);
400 }
401 
T getParameter(std::string const &) const
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
virtual int nrows() const =0
unsigned int fed
virtual int rowsperroc() const =0
ParameterDescriptionBase * addUntracked(U const &iLabel, T const &value)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
void accumulate(edm::Event const &, const edm::EventSetup &iSetup) final
SiPixelStatusProducer(const edm::ParameterSet &)
unsigned int link
identify pixel inside single ROC
Definition: LocalPixel.h:7
LuminosityBlockNumber_t luminosityBlock() const
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
virtual int colsperroc() const =0
void put(std::unique_ptr< PROD > product)
Put a new product.
virtual void onlineRocColRow(const DetId &detId, int offlineRow, int offlineCol, int &roc, int &row, int &col) final
unsigned int idInDetUnit() const
id of this ROC in DetUnit etermined by token path
Definition: PixelROC.h:37
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
RunNumber_t run() const
Definition: value.py:1
ParameterDescriptionBase * add(U const &iLabel, T const &value)
bool isValid() const
Definition: HandleBase.h:70
double collumn and pixel ID in double collumn representation
Definition: LocalPixel.h:19
Definition: DetId.h:17
int toCabling(sipixelobjects::ElectronicIndex &cabling, const sipixelobjects::DetectorIndex &detector) const
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
void add(std::string const &label, ParameterSetDescription const &psetDescription)
void beginLuminosityBlock(edm::LuminosityBlock const &lumiSeg, const edm::EventSetup &iSetup) final
virtual int ncolumns() const =0
sipixelobjects::PixelROC const * toRoc(int link, int roc) const
T get() const
Definition: EventSetup.h:73
col
Definition: cuy.py:1010
virtual int indexROC(int irow, int icol, int nROCcolumns) final
void endLuminosityBlockProduce(edm::LuminosityBlock &lumiSeg, const edm::EventSetup &iSetup) final
def move(src, dest)
Definition: eostools.py:511
void endLuminosityBlock(edm::LuminosityBlock const &lumiSeg, const edm::EventSetup &iSetup) final