CMS 3D CMS Logo

HcalDigisClient.cc
Go to the documentation of this file.
2 // -*- C++ -*-
3 //
4 // Package: HcalDigisClient
5 // Class: HcalDigisClient
6 //
14 //
15 // Original Author: Ali Fahim,22 R-013,+41227672649,
16 // Created: Wed Mar 23 22:54:28 CET 2011
17 //
18 //
19 
20 // system include files
21 
23  outputFile_ = iConfig.getUntrackedParameter<std::string>("outputFile", "HcalDigisClient.root");
24  dirName_ = iConfig.getParameter<std::string>("DQMDirName");
25  msm_ = new std::map<std::string, MonitorElement*>();
26 }
27 
29 
32  std::vector<MonitorElement*> hcalMEs;
33  // Since out folders are fixed to three, we can just go over these three folders
34  // i.e., CaloTowersV/CaloTowersTask, HcalRecHitsV/HcalRecHitTask, NoiseRatesV/NoiseRatesTask.
35  std::vector<std::string> fullPathHLTFolders = ig.getSubdirs();
36  for (unsigned int i = 0; i < fullPathHLTFolders.size(); i++) {
37  ig.setCurrentFolder(fullPathHLTFolders[i]);
38  std::vector<std::string> fullSubPathHLTFolders = ig.getSubdirs();
39  for (unsigned int j = 0; j < fullSubPathHLTFolders.size(); j++) {
40  if (strcmp(fullSubPathHLTFolders[j].c_str(), "HcalDigisV/HcalDigiTask") == 0) {
41  hcalMEs = ig.getContents(fullSubPathHLTFolders[j]);
42  ig.setCurrentFolder("HcalDigisV/HcalDigiTask");
43  if (!HcalDigisEndjob(hcalMEs, "HB", ib))
44  edm::LogError("HcalDigisClient") << "Error in HcalDigisEndjob! HB";
45  if (!HcalDigisEndjob(hcalMEs, "HE", ib))
46  edm::LogError("HcalDigisClient") << "Error in HcalDigisEndjob! HE";
47  if (!HcalDigisEndjob(hcalMEs, "HO", ib))
48  edm::LogError("HcalDigisClient") << "Error in HcalDigisEndjob! HO";
49  if (!HcalDigisEndjob(hcalMEs, "HF", ib))
50  edm::LogError("HcalDigisClient") << "Error in HcalDigisEndjob! HF";
51  }
52  }
53  }
54 }
55 
56 int HcalDigisClient::HcalDigisEndjob(const std::vector<MonitorElement*>& hcalMEs,
57  std::string subdet_,
59  using namespace std;
60  string strtmp;
61 
62  MonitorElement* nevtot(nullptr);
63 
64  std::vector<MonitorElement*> ieta_iphi_occupancy_maps;
65  std::vector<std::string> depthID;
66 
67  edm::LogVerbatim("HcalDigisClient") << " Number of histos " << hcalMEs.size();
68 
69  for (unsigned int ih = 0; ih < hcalMEs.size(); ih++) {
70  if (hcalMEs[ih]->getName() == "nevtot") {
71  nevtot = hcalMEs[ih];
72  continue;
73  }
74 
75  //We search the occupancy maps corresponding to this subdetector
76  if ((hcalMEs[ih]->getName().find("HcalDigiTask_ieta_iphi_occupancy_map_depth") != std::string::npos) &&
77  (hcalMEs[ih]->getName().find(subdet_) != std::string::npos)) {
78  ieta_iphi_occupancy_maps.push_back(hcalMEs[ih]);
79 
80  std::string start = "depth";
81  std::string end = "_H";
82 
83  int position = hcalMEs[ih]->getName().find(start) + start.length();
84  int length = hcalMEs[ih]->getName().find(end) - position;
85 
86  depthID.push_back(hcalMEs[ih]->getName().substr(position, length));
87 
88  continue;
89  }
90  }
91 
92  if (hcalMEs.empty()) {
93  edm::LogError("HcalDigisClient") << "No nevtot or maps histo found...";
94  return 0;
95  }
96  if (!nevtot) {
97  edm::LogError("HcalDigisClient") << "No nevtot histoo found...";
98  return 0;
99  }
100  if (ieta_iphi_occupancy_maps.empty()) {
101  edm::LogError("HcalDigisClient") << "No maps histos found...";
102  return 0;
103  }
104 
105  int ev = nevtot->getEntries();
106 
107  if (ev <= 0) {
108  edm::LogError("HcalDigisClient") << "normalization factor <= 0!";
109  return 0;
110  }
111 
112  float fev = (float)nevtot->getEntries();
113 
114  int depths = ieta_iphi_occupancy_maps.size();
115 
116  HistLim ietaLim(85, -42.5, 42.5);
117 
118  for (int depth = 1; depth <= depths; depth++) {
119  strtmp = "HcalDigiTask_occupancy_vs_ieta_depth" + str(depth) + "_" + subdet_;
120  book1D(ib, strtmp, ietaLim);
121  }
122 
123  std::vector<float> sumphi(depths, 0);
124  std::vector<float> sumphie(depths, 0);
125 
126  float phi_factor;
127  float cnorm;
128  float enorm;
129 
130  for (int depth = 1; depth <= depths; depth++) {
131  int nx = ieta_iphi_occupancy_maps[depth - 1]->getNbinsX();
132  int ny = ieta_iphi_occupancy_maps[depth - 1]->getNbinsY();
133 
134  for (int i = 1; i <= nx; i++) {
135  for (int j = 1; j <= ny; j++) {
136  // occupancies
137  cnorm = ieta_iphi_occupancy_maps[depth - 1]->getBinContent(i, j) / fev;
138  enorm = ieta_iphi_occupancy_maps[depth - 1]->getBinError(i, j) / fev;
139  ieta_iphi_occupancy_maps[depth - 1]->setBinContent(i, j, cnorm);
140  ieta_iphi_occupancy_maps[depth - 1]->setBinError(i, j, enorm);
141 
142  } //for loop over NbinsYU
143  } //for loop over NbinsX
144  } //for loop over the occupancy maps
145 
146  for (int i = 1; i <= 82; i++) {
147  int ieta = i - 42; // -41 -1, 0 40
148  if (ieta >= 0)
149  ieta += 1; // -41 -1, 1 41 - to make it detector-like
150 
151  if (ieta >= -20 && ieta <= 20) {
152  phi_factor = 72.;
153  } else {
154  if (ieta >= 40 || ieta <= -40)
155  phi_factor = 18.;
156  else
157  phi_factor = 36.;
158  }
159 
160  //zero the sumphi and sumphie vector at the start of each ieta ring
161  sumphi.assign(depths, 0);
162  sumphie.assign(depths, 0);
163 
164  for (int iphi = 1; iphi <= 72; iphi++) {
165  for (int depth = 1; depth <= depths; depth++) {
166  int binIeta = ieta_iphi_occupancy_maps[depth - 1]->getTH2F()->GetXaxis()->FindBin(ieta);
167  int binIphi = ieta_iphi_occupancy_maps[depth - 1]->getTH2F()->GetYaxis()->FindBin(iphi);
168 
169  float content = ieta_iphi_occupancy_maps[depth - 1]->getBinContent(binIeta, binIphi);
170  float econtent = ieta_iphi_occupancy_maps[depth - 1]->getBinError(binIeta, binIphi);
171 
172  sumphi[depth - 1] += content;
173  sumphie[depth - 1] += econtent * econtent;
174 
175  } //for loop over depths
176  } //for loop over phi
177 
178  //double deta = double(ieta);
179 
180  // occupancies vs ieta
181  for (int depth = 1; depth <= depths; depth++) {
182  strtmp = "HcalDigiTask_occupancy_vs_ieta_depth" + depthID[depth - 1] + "_" + subdet_;
183  MonitorElement* ME = msm_->find(strtmp)->second;
184  int ietabin = ME->getTH1F()->GetXaxis()->FindBin(float(ieta));
185 
186  if (sumphi[depth - 1] > 1.e-30) {
187  cnorm = sumphi[depth - 1] / phi_factor;
188  enorm = sqrt(sumphie[depth - 1]) / phi_factor;
189  ME->setBinContent(ietabin, cnorm);
190  ME->setBinError(ietabin, enorm);
191  }
192  }
193  } // end of i-loop
194 
195  return 1;
196 }
197 
199  if (!msm_->count(name))
200  return nullptr;
201  else
202  return msm_->find(name)->second;
203 }
204 
206  std::stringstream out;
207  out << x;
208  return out.str();
209 }
210 
211 double HcalDigisClient::integralMETH2D(MonitorElement* ME, int i0, int i1, int j0, int j1) {
212  double sum(0);
213  for (int i = i0; i <= i1; i++) {
214  for (int j = j0; j <= j1; j++) {
215  sum += ME->getBinContent(i, j);
216  }
217  }
218 
219  return sum;
220 }
221 
223  int nx = ME->getNbinsX();
224  int ny = ME->getNbinsY();
225 
226  double content(0);
227  double error(0);
228  for (int i = 1; i <= nx; i++) {
229  for (int j = 1; j <= ny; j++) {
230  content = ME->getBinContent(i, j);
231  error = ME->getBinError(i, j);
232  content *= s;
233  error *= s;
234  ME->setBinContent(i, j, content);
235  ME->setBinError(i, j, error);
236  }
237  }
238 }
239 
240 //define this as a plug-in
Definition: start.py:1
Log< level::Info, true > LogVerbatim
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
std::map< std::string, MonitorElement * > * msm_
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
std::string dirName_
void scaleMETH2D(MonitorElement *ME, double s)
virtual void runClient(DQMStore::IBooker &ib, DQMStore::IGetter &ig)
Log< level::Error, false > LogError
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
HcalDigisClient(const edm::ParameterSet &)
T getUntrackedParameter(std::string const &, T const &) const
std::string outputFile_
Definition: ME.h:11
~HcalDigisClient() override
T sqrt(T t)
Definition: SSEVec.h:19
int HcalDigisEndjob(const std::vector< MonitorElement *> &hcalMEs, std::string subdet_, DQMStore::IBooker &ib)
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
virtual double getEntries() const
get # of entries
std::string str(int x)
void book1D(DQMStore::IBooker &ib, std::string name, int n, double min, double max)
MonitorElement * monitor(std::string name)
double integralMETH2D(MonitorElement *ME, int i0, int i1, int j0, int j1)
static int position[264][3]
Definition: ReadPGInfo.cc:289
std::string getName(const G4String &)
Definition: ForwardName.cc:3
ib
Definition: cuy.py:661
virtual std::vector< dqm::harvesting::MonitorElement * > getContents(std::string const &path) const
Definition: DQMStore.cc:625
virtual DQM_DEPRECATED std::vector< std::string > getSubdirs() const
Definition: DQMStore.cc:739