11 #include <vxd/dataobjects/VxdID.h>
12 #include <pxd/calibration/PXDHotPixelMaskCalibrationAlgorithm.h>
13 #include <pxd/dbobjects/PXDMaskedPixelPar.h>
14 #include <pxd/dbobjects/PXDDeadPixelPar.h>
15 #include <pxd/dbobjects/PXDOccupancyInfoPar.h>
17 #include <boost/format.hpp>
29 PXDHotPixelMaskCalibrationAlgorithm::PXDHotPixelMaskCalibrationAlgorithm():
CalibrationAlgorithm(
"PXDHotPixelMaskCollector"),
30 forceContinueMasking(false), minEvents(10000), minHits(20), pixelMultiplier(7), maskDrains(true),
31 drainMultiplier(7), maskRows(true), rowMultiplier(7)
34 " -------------------------- PXDHotPixelMak Calibration Algorithm ------------------------\n"
36 " Algorithm which masks hot pixels with too large occupancy and dead pixels w/o no hits. \n"
37 " ----------------------------------------------------------------------------------------\n"
44 auto collector_pxdhits = getObjectPtr<TH1I>(
"PXDHits");
45 auto collector_pxdhitcounts = getObjectPtr<TH1I>(
"PXDHitCounts");
48 auto nevents = collector_pxdhits->GetEntries();
51 B2INFO(
"Not enough data: Only " << nevents <<
" events were collected!");
54 B2WARNING(
"Not enough data: Only " << nevents <<
" events were collected! The masking continous but the mask may be empty.");
58 B2RESULT(
"Found total of " << nevents <<
" events in collected data.");
61 unsigned long long int nPXDHits = 0;
63 for (
auto sensBin = 1; sensBin <= collector_pxdhitcounts->GetXaxis()->GetNbins(); sensBin++) {
65 string sensorDescr = collector_pxdhitcounts->GetXaxis()->GetBinLabel(sensBin);
66 VxdID id(sensorDescr);
70 unsigned long long int nSensorHits = collector_pxdhitcounts->GetBinContent(sensBin);
71 nPXDHits += nSensorHits;
73 B2RESULT(
"Number of hits for sensor sensor=" <<
id <<
" is " << nSensorHits);
78 if (nPXDHits < minPXDHits) {
80 B2INFO(
"Not enough data: Only " << nPXDHits <<
" raw hits were collected!");
83 B2WARNING(
"Not enough data: Only " << nPXDHits <<
" raw hits were collected! The masking continous but the mask may be empty.");
87 B2RESULT(
"Found total of " << nPXDHits <<
" raw hits in collected data.");
90 map<VxdID, double> medianOfHitsMap;
91 for (
auto sensBin = 1; sensBin <= collector_pxdhitcounts->GetXaxis()->GetNbins(); sensBin++) {
93 string sensorDescr = collector_pxdhitcounts->GetXaxis()->GetBinLabel(sensBin);
94 VxdID id(sensorDescr);
96 medianOfHitsMap[id] = 0.0;
99 string name = str(format(
"PXD_%1%_PixelHitmap") %
id.
getID());
100 auto collector_pxdhitmap = getObjectPtr<TH1I>(name.c_str());
103 if (collector_pxdhitmap ==
nullptr)
continue;
106 int nBins = collector_pxdhitmap->GetXaxis()->GetNbins();
108 vector<double> hitVec(nBins);
110 for (
auto bin = 1; bin <= nBins; bin++) {
111 hitVec[bin - 1] = (double) collector_pxdhitmap->GetBinContent(bin);
113 double medianNumberOfHits;
114 TMath::Quantiles(nBins, 1, &hitVec[0], &medianNumberOfHits, &prob, kFALSE);
115 if (medianNumberOfHits <= 0) {
116 B2WARNING(
"Median number of hits per senor is smaller <1. Raise median to 1 instead.");
117 medianNumberOfHits = 1;
119 B2RESULT(
"Median of hits is " << medianNumberOfHits <<
" for sensor " <<
id);
123 medianOfHitsMap[id] = medianNumberOfHits;
126 if (medianNumberOfHits <
minHits) {
128 B2INFO(
"Not enough data: Median number of his is only " << medianNumberOfHits <<
"!");
131 B2WARNING(
"Not enough data: Median number of hits is only " << medianNumberOfHits <<
132 "! The masking continous but the mask may be empty.");
150 for (
auto sensBin = 1; sensBin <= collector_pxdhitcounts->GetXaxis()->GetNbins(); sensBin++) {
152 string sensorDescr = collector_pxdhitcounts->GetXaxis()->GetBinLabel(sensBin);
153 VxdID id(sensorDescr);
156 int nSensorHits = collector_pxdhitcounts->GetBinContent(sensBin);
157 if (nSensorHits == 0) {
163 string name = str(format(
"PXD_%1%_PixelHitmap") %
id.
getID());
164 auto collector_pxdhitmap = getObjectPtr<TH1I>(name.c_str());
165 if (collector_pxdhitmap ==
nullptr) {
166 B2WARNING(
"Cannot find PixelHitmap although there should be hits. This is strange!");
170 double medianNumberOfHits = medianOfHitsMap[id];
171 int nBins = collector_pxdhitmap->GetXaxis()->GetNbins();
174 if (medianNumberOfHits >=
minHits) {
178 vector<int> hitsAlongDrain(
c_nDrains, 0);
181 for (
auto bin = 1; bin <= nBins; bin++) {
186 int drainID = uCell * 4 + vCell % 4;
187 int nhits = collector_pxdhitmap->GetBinContent(bin);
188 hitsAlongDrain[drainID] += nhits;
189 hitsAlongRow[vCell] += nhits;
193 B2INFO(
"Entering masking of dead rows ...");
194 for (
auto vCell = 0; vCell <
c_nVCells; vCell++) {
196 int nhits = hitsAlongRow[vCell];
200 B2RESULT(
"Dead row with vCell=" << vCell <<
" on sensor " <<
id);
205 B2INFO(
"Entering masking of dead drains ...");
206 for (
auto drainID = 0; drainID <
c_nDrains; drainID++) {
208 int nhits = hitsAlongDrain[drainID];
212 B2RESULT(
"Dead drain line at drainID=" << drainID <<
" on sensor " <<
id);
217 B2INFO(
"Entering masking of single dead pixels ...");
218 for (
auto bin = 1; bin <= nBins; bin++) {
220 int nhits = collector_pxdhitmap->GetBinContent(bin);
224 int drainID = uCell * 4 + vCell % 4;
226 if (nhits == 0 && !deadPixelsPar->
isDeadRow(
id.getID(), vCell) && !deadPixelsPar->
isDeadDrain(
id.getID(), drainID)) {
229 B2RESULT(
"Dead single pixel with ucell=" << uCell <<
", vcell=" << vCell <<
" on sensor " <<
id);
237 vector<float> unmaskedHitsAlongDrain(
c_nDrains, 0);
238 vector<int> unmaskedCellsAlongDrain(
c_nDrains, 0);
241 vector<float> unmaskedHitsAlongRow(
c_nVCells, 0);
242 vector<int> unmaskedCellsAlongRow(
c_nVCells, 0);
246 B2RESULT(
"Pixel hit threshold is " << pixelHitThr <<
" for sensor " <<
id);
249 for (
auto bin = 1; bin <= nBins; bin++) {
254 int drainID = uCell * 4 + vCell % 4;
257 float nhits = collector_pxdhitmap->GetBinContent(bin);
260 if (nhits > pixelHitThr) {
264 B2RESULT(
"Masking single pixel with ucell=" << uCell <<
", vcell=" << vCell <<
" on sensor " <<
id);
270 ++unmaskedCellsAlongDrain[drainID];
271 unmaskedHitsAlongDrain[drainID] += nhits;
272 ++unmaskedCellsAlongRow[vCell];
273 unmaskedHitsAlongRow[vCell] += nhits;
279 B2RESULT(
"Drain hit threshold is " << drainHitThr <<
" for sensor " <<
id);
281 for (
auto drainID = 0; drainID <
c_nDrains; drainID++) {
282 if (unmaskedCellsAlongDrain[drainID] > 0) {
284 float nhits = unmaskedHitsAlongDrain[drainID] / unmaskedCellsAlongDrain[drainID];
286 if (nhits > drainHitThr) {
287 for (
auto iGate = 0; iGate < 192; iGate++) {
288 int uCell = drainID / 4;
289 int vCell = drainID % 4 + iGate * 4;
292 B2RESULT(
"Masking drain line with drainID=" << drainID <<
" on sensor " <<
id);
300 B2RESULT(
"Row hit threshold is " << rowHitThr <<
" for sensor " <<
id);
302 for (
auto vCell = 0; vCell <
c_nVCells; vCell++) {
303 if (unmaskedCellsAlongRow[vCell] > 0) {
305 float nhits = unmaskedHitsAlongRow[vCell] / unmaskedCellsAlongRow[vCell];
307 if (nhits > rowHitThr) {
308 for (
auto uCell = 0; uCell <
c_nUCells; uCell++)
311 B2RESULT(
"Masking complete row with vCell=" << vCell <<
" on sensor " <<
id);
321 int numberOfUnmaskedHits = 0;
323 for (
auto bin = 1; bin <= nBins; bin++) {
327 if (maskedPixelsPar->
pixelOK(
id.getID(), pixID)) {
328 numberOfUnmaskedHits += collector_pxdhitmap->GetBinContent(bin);
333 float meanOccupancyAfterMasking = (float)numberOfUnmaskedHits / nevents /
c_nVCells /
c_nUCells;
334 B2RESULT(
"Hotpixel filtered occupancy sensor=" <<
id <<
" is " << meanOccupancyAfterMasking);
340 auto id = elem.first;
341 auto singles = elem.second;
342 B2RESULT(
"SensorID " <<
VxdID(
id) <<
" has filtered occupancy of " << occupancyInfoPar->
getOccupancy(
id));
343 B2RESULT(
"SensorID " <<
VxdID(
id) <<
" has fraction of masked pixels of " << (
float)singles.size() / (
c_nVCells *
c_nUCells));
352 B2INFO(
"PXDHotPixelMask Calibration Successful");