CMS 3D CMS Logo

SiStripCorrelateBadStripAndNoise.cc
Go to the documentation of this file.
2 
5 
7  : cacheID_quality(0xFFFFFFFF), cacheID_noise(0xFFFFFFFF) {
8  // now do what ever initialization is needed
9  if (!edm::Service<SiStripDetInfoFileReader>().isAvailable()) {
10  edm::LogError("TkLayerMap") << "\n------------------------------------------"
11  "\nUnAvailable Service SiStripDetInfoFileReader: please insert in "
12  "the configuration file an instance like"
13  "\n\tprocess.SiStripDetInfoFileReader = "
14  "cms.Service(\"SiStripDetInfoFileReader\")"
15  "\n------------------------------------------";
16  }
17 
19  file = new TFile("correlTest.root", "RECREATE");
20  tkmap = new TrackerMap();
21 }
22 
24 
25 //
26 
29  return;
32 
33  edm::LogInfo("") << "[SiStripCorrelateBadStripAndNoise::beginRun] cacheID_quality " << cacheID_quality
34  << " cacheID_noise " << cacheID_noise << std::endl;
35 
38 
39  DoAnalysis(es);
40 }
41 
43  // Retrieve tracker topology from geometry
45  es.get<TrackerTopologyRcd>().get(tTopoHandle);
46  const TrackerTopology *const tTopo = tTopoHandle.product();
47 
48  // Loop on quality bad stirps
49  // for each strip, look at the noise
50  // evalaute the mean apv noise and the ratio among strip noise and
51  // meanApvNoise put the value in the histo in terms of ratio Vs percentage of
52  // badStrips per APV
53 
54  // Fill an histo per subdet and layer (and plus && minus for TEC/TID)
55  edm::LogInfo("") << "[Doanalysis]";
56  iterateOnDets(tTopo);
57 }
58 
62 
63  for (SiStripBadStrip::RegistryIterator rp = rbegin; rp != rend; ++rp) {
64  const uint32_t detid = rp->detid;
65 
67  qualityHandle_->getDataVectorBegin() + rp->iend);
68  iterateOnBadStrips(detid, tTopo, sqrange);
69  }
70 }
71 
73  const TrackerTopology *tTopo,
74  SiStripQuality::Range &sqrange) {
75  float percentage = 0;
76  for (int it = 0; it < sqrange.second - sqrange.first; it++) {
77  unsigned int firstStrip = qualityHandle_->decode(*(sqrange.first + it)).firstStrip;
78  unsigned int range = qualityHandle_->decode(*(sqrange.first + it)).range;
79 
80  correlateWithNoise(detid, tTopo, firstStrip, range);
81 
82  edm::LogInfo("range") << range;
83  percentage += range;
84  }
85  if (percentage != 0)
86  percentage /= 128. * fr->getNumberOfApvsAndStripLength(detid).first;
87  if (percentage > 1)
88  edm::LogError("SiStripQualityStatistics") << "PROBLEM detid " << detid << " value " << percentage << std::endl;
89 
90  //------- Global Statistics on percentage of bad components along the IOVs
91  //------//
92  if (percentage != 0)
93  edm::LogInfo("") << "percentage " << detid << " " << percentage;
94 }
95 
97  const TrackerTopology *tTopo,
98  const uint32_t &firstStrip,
99  const uint32_t &range) {
100  std::vector<TH2F *> histos;
101 
102  SiStripNoises::Range noiseRange = noiseHandle_->getRange(detid);
103  edm::LogInfo("Domenico") << "detid " << detid << " first " << firstStrip << " range " << range;
104  float meanAPVNoise = getMeanNoise(noiseRange, firstStrip / 128, 128);
105 
106  // float meanNoiseHotStrips=getMeanNoise(noiseRange,firstStrip,range);
107  for (size_t theStrip = firstStrip; theStrip < firstStrip + range; theStrip++) {
108  float meanNoiseHotStrips = getMeanNoise(noiseRange, theStrip, 1);
109 
110  // Get the histogram for this detid
111  getHistos(detid, tTopo, histos);
112  float yvalue = range < 21 ? 1. * range : 21;
113 
114  for (size_t i = 0; i < histos.size(); ++i)
115  histos[i]->Fill(meanNoiseHotStrips / meanAPVNoise - 1., yvalue);
116 
117  if (meanNoiseHotStrips / meanAPVNoise - 1. < -0.3)
118  tkmap->fillc(detid, 0xFF0000);
119  else
120  tkmap->fillc(detid, 0x0000FF);
121  }
122 }
123 
125  const uint32_t &firstStrip,
126  const uint32_t &range) {
127  float mean = 0;
128  for (size_t istrip = firstStrip; istrip < firstStrip + range; istrip++) {
129  mean += noiseHandle_->getNoise(istrip, noiseRange);
130  }
131  return mean / (1. * range);
132 }
133 
135  const TrackerTopology *tTopo,
136  std::vector<TH2F *> &histos) {
137  histos.clear();
138 
139  int subdet = -999;
140  int component = -999;
141  SiStripDetId a(detid);
142  if (a.subdetId() == 3) {
143  subdet = 0;
144  component = tTopo->tibLayer(detid);
145  } else if (a.subdetId() == 4) {
146  subdet = 1;
147  component = tTopo->tidSide(detid) == 2 ? tTopo->tidWheel(detid) : tTopo->tidWheel(detid) + 3;
148  } else if (a.subdetId() == 5) {
149  subdet = 2;
150  component = tTopo->tobLayer(detid);
151  } else if (a.subdetId() == 6) {
152  subdet = 3;
153  component = tTopo->tecSide(detid) == 2 ? tTopo->tecWheel(detid) : tTopo->tecWheel(detid) + 9;
154  }
155 
156  int index = 100 + subdet * 100 + component;
157 
158  histos.push_back(getHisto(subdet));
159  histos.push_back(getHisto(index));
160 }
161 
162 TH2F *SiStripCorrelateBadStripAndNoise::getHisto(const long unsigned int &index) {
163  if (vTH2.size() < index + 1)
164  vTH2.resize(index + 1, nullptr);
165 
166  if (vTH2[index] == nullptr) {
167  char name[128];
168  sprintf(name, "%lu", index);
169  edm::LogInfo("") << "[getHisto] creating index " << index << std::endl;
170  vTH2[index] = new TH2F(name, name, 50, -2., 2., 21, 0.5, 21.5);
171  }
172 
173  return vTH2[index];
174 }
175 
177  for (size_t i = 0; i < vTH2.size(); i++)
178  if (vTH2[i] != nullptr)
179  vTH2[i]->Write();
180 
181  file->Write();
182  file->Close();
183 
184  tkmap->save(true, 0, 0, "testTkMap.png");
185 }
unsigned short range
void beginRun(const edm::Run &run, const edm::EventSetup &es) override
unsigned int tibLayer(const DetId &id) const
TH2F * getHisto(const long unsigned int &index)
const std::pair< unsigned short, double > getNumberOfApvsAndStripLength(uint32_t detId) const
void iterateOnBadStrips(const uint32_t &detid, const TrackerTopology *tTopo, SiStripQuality::Range &sqrange)
unsigned int tidWheel(const DetId &id) const
void getHistos(const uint32_t &detid, const TrackerTopology *tTopo, std::vector< TH2F * > &histos)
Registry::const_iterator RegistryIterator
static float getNoise(uint16_t strip, const Range &range)
Definition: SiStripNoises.h:74
unsigned long long getQualityCache(const edm::EventSetup &eSetup)
RegistryIterator getRegistryVectorEnd() const
unsigned int tidSide(const DetId &id) const
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
void iterateOnDets(const TrackerTopology *tTopo)
void save(bool print_total=true, float minval=0., float maxval=0., std::string s="svgmap.svg", int width=1500, int height=800)
Definition: TrackerMap.cc:699
unsigned long long getNoiseCache(const edm::EventSetup &eSetup)
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:41
edm::ESHandle< SiStripQuality > qualityHandle_
SiStripCorrelateBadStripAndNoise(const edm::ParameterSet &)
void fillc(int idmod, int RGBcode)
Definition: TrackerMap.h:109
float getMeanNoise(const SiStripNoises::Range &noiseRange, const uint32_t &first, const uint32_t &range)
ContainerIterator getDataVectorBegin() const
Detector identifier class for the strip tracker.
Definition: SiStripDetId.h:17
unsigned short firstStrip
RegistryIterator getRegistryVectorBegin() const
double a
Definition: hdecay.h:121
const Range getRange(const uint32_t detID) const
std::pair< ContainerIterator, ContainerIterator > Range
T get() const
Definition: EventSetup.h:71
std::pair< ContainerIterator, ContainerIterator > Range
Definition: SiStripNoises.h:50
unsigned int tecWheel(const DetId &id) const
T const * product() const
Definition: ESHandle.h:86
data decode(const unsigned int &value) const
unsigned int tobLayer(const DetId &id) const
Definition: Run.h:45
void correlateWithNoise(const uint32_t &detid, const TrackerTopology *tTopo, const uint32_t &firstStrip, const uint32_t &range)
unsigned int tecSide(const DetId &id) const