CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
RPCEfficiencySecond.cc
Go to the documentation of this file.
1 /***************************************
2 Original Author:Camilo Carrillo
3 ****************************************/
4 #include <sstream>
5 
19 
21 
23  SaveFile = iConfig.getUntrackedParameter<bool>("SaveFile", false);
24  NameFile = iConfig.getUntrackedParameter<std::string>("NameFile","RPCEfficiency.root");
25  folderPath = iConfig.getUntrackedParameter<std::string>("folderPath","RPC/RPCEfficiency/");
26  debug = iConfig.getUntrackedParameter<bool>("debug",false);
27  numberOfDisks_ = iConfig.getUntrackedParameter<int>("NumberOfEndcapDisks", 4);
28  innermostRings_ = iConfig.getUntrackedParameter<int>("NumberOfInnermostEndcapRings", 2);
29 
30 }
31 
33 
35 
37  rpcdqm::utils rpcUtils;
38  //Barrel
39  std::stringstream folderName;
40  std::stringstream MeName; std::stringstream MeTitle;
41  for (int w =-2; w<=2; w++){
42 
43  folderName.str("");
44  folderName<<folderPath<<"Wheel_"<<w;
45 
46  dbe->setCurrentFolder(folderName.str());
47 
48  MeName.str("");
49  MeName<<"EffDistroWheel_"<<w;
50  MeTitle.str("");
51  MeTitle<<"Efficiency Distribution for Wheel "<<w;
52  EffDistroW[w+2]=dbe->book1D(MeName.str(),MeTitle.str(),20,0.5,100.5);
53 
54  MeName.str("");
55  MeName<<"GlobEfficiencyWheel_"<<w;
56  MeTitle.str("");
57  MeTitle<<"Efficiency Wheel "<<w;
58  EffGlobW[w+2]=dbe->book1D(MeName.str(),MeTitle.str(),206,0.5,206.5);
59 
61  MeName.str("");
62  MeName<<"Efficiency_Roll_vs_Sector_Wheel_"<<w;
63  WheelSummary[w+2] = dbe->book2D(MeName.str(),MeName.str(),12, 0.5,12.5, 21, 0.5, 21.5);
64  rpcUtils.labelXAxisSector( WheelSummary[w+2] );
65  rpcUtils.labelYAxisRoll( WheelSummary[w+2], 0, w, true);
66 
67  dbe->setCurrentFolder(folderPath+"Azimutal/");
68  MeName.str("");
69  MeName<<"AzimutalDistroW"<<w;
70  MeTitle.str("");
71  MeTitle<<"Efficiency per Sector Wheel "<<w;
72  sectorEffW[w+2]= dbe->book1D(MeName.str(),MeTitle.str(),12,0.5,12.5);
73 
74  MeName.str("");
75  MeName<<"AzimutalDistroW"<<w<<"Ex";
76  MeTitle.str("");
77  MeTitle<<"Expected per Sector Wheel "<<w;
78  ExsectorEffW[w+2]= dbe->book1D(MeName.str(),MeTitle.str(),12,0.5,12.5);
79 
80  MeName.str("");
81  MeName<<"AzimutalDistroW"<<w<<"Oc";
82  MeTitle.str("");
83  MeTitle<<"Occupancy per Sector Wheel "<<w;
84  OcsectorEffW[w+2]= dbe->book1D(MeName.str(),MeTitle.str(),12,0.5,12.5);
85 
86  dbe->setCurrentFolder(folderPath+"BarrelPerLayer/");
87  MeName.str("");
88  MeName<<"ExpLayerW"<<w;
89  MeTitle.str("");
90  MeTitle<<"Expected Wheel "<<w;
91  ExpLayerW[w+2]= dbe->book1D(MeName.str(),MeTitle.str(),6,0.5,6.5);
92 
93  MeName.str("");
94  MeName<<"ObsLayerW"<<w;
95  MeTitle.str("");
96  MeTitle<<"Observed Wheel "<<w;
97  ObsLayerW[w+2]= dbe->book1D(MeName.str(),MeTitle.str(),6,0.5,6.5);
98  }
99 
100  //EndCap
101  int index = 0;
102  for (int d = (-1 *numberOfDisks_); d<=numberOfDisks_; d++){
103 
104  if (d==0) {continue;}
105 
106  folderName.str("");
107  folderName<<folderPath<<"Disk_"<<d;
108  dbe->setCurrentFolder(folderName.str());
109 
110  MeName.str("");
111  MeName<<"EffDistroDisk_"<<d;
112  MeTitle.str("");
113  MeTitle<<"Efficiency Distribution Disk "<<d;
114  EffDistroD[index]=dbe->book1D(MeName.str(),MeTitle.str(),20,0.5,100.5);
115 
116  MeName.str("");
117  MeName<<"GlobEfficiencyDisk_"<<d;
118  MeTitle.str("");
119  MeTitle<<"Efficiency Disk "<<d;
120  EffGlobD[index] = dbe->book1D(MeName.str(),MeTitle.str(),218,0.5,218.5);
121 
122 
124  MeName.str("");
125  MeName<<"Efficiency_Roll_vs_Segment_Disk_"<<d;
126  DiskSummary[index] = dbe->book2D(MeName.str(),MeName.str(),36,0.5,36.5,6,0.5,6.5);
127 
128  rpcUtils.labelXAxisSegment(DiskSummary[index]);
129  rpcUtils.labelYAxisRing(DiskSummary[index], innermostRings_ , true);
130 
131 
132 
133  dbe->setCurrentFolder(folderPath+"Azimutal/");
134  MeName.str("");
135  MeName<<"GregDistroR2D"<<d;
136  MeTitle.str("");
137  MeTitle<<"Efficiency for Station "<<d<<" Ring 2";
138  GregR2D[index]= dbe->book1D(MeName.str(),MeTitle.str(),36,0.5,36.5);
139 
140  MeName.str("");
141  MeName<<"GregDistroR3D"<<d;
142  MeTitle.str("");
143  MeTitle<<"Efficiency for Station "<<d<<" Ring 3";
144  GregR3D[index]= dbe->book1D(MeName.str(),MeTitle.str(),36,0.5,36.5);
145 
146  MeName.str("");
147  MeName<<"OcGregDistroR2D"<<d;
148  MeTitle.str("");
149  MeTitle<<"Occupancy Distribution for Station "<<d<<" Ring 2";
150  OcGregR2D[index]= dbe->book1D(MeName.str(),MeTitle.str(),36,0.5,36.5);
151 
152  MeName.str("");
153  MeName<<"OcGregDistroR3D"<<d;
154  MeTitle.str("");
155  MeTitle<<"Occupancy Distribution for Station "<<d<<" Ring 3";
156  OcGregR3D[index]= dbe->book1D(MeName.str(),MeTitle.str(),36,0.5,36.5);
157 
158  MeName.str("");
159  MeName<<"ExGregDistroR2D"<<d;
160  MeTitle.str("");
161  MeTitle<<"Expected Distribution for Station "<<d<<" Ring 2";
162  ExGregR2D[index]= dbe->book1D(MeName.str(),MeTitle.str(),36,0.5,36.5);
163 
164 
165  MeName.str("");
166  MeName<<"ExGregDistroR3D"<<d;
167  MeTitle.str("");
168  MeTitle<<"Expected Distribution for Station "<<d<<" Ring 3";
169  ExGregR3D[index]= dbe->book1D(MeName.str(),MeTitle.str(),36,0.5,36.5);
170 
171  index++;
172  }
173 }
174 
176 
178 
179  LogDebug("rpcefficiencysecond")<<"Getting the RPC Geometry";
180 
182  iSetup.get<MuonGeometryRecord>().get(rpcGeo);
183 
184 
185  //Setting Labels in Summary Label.
186  std::stringstream binLabel;
187  for (int w = -2; w<=2 ;w++){
188  for(int x = 1;x<=12;x++){
189  for(int y = 1;y<=21;y++){
190  WheelSummary[w+2]->setBinContent(x,y,-1);
191  }
192  }
193  }
194 
195 
196  for (int d = 0 ; d<(numberOfDisks_*2); d++){
197  for(int x = 1;x<=36;x++){
198  for(int y = 1;y<=6;y++){
199  DiskSummary[d]->setBinContent(x,y,-1);
200  }
201  }
202  }
203 
204 
205  binLabel.str("");
206  int indexWheel[5];
207  for(int j=0;j<5;j++){
208  indexWheel[j]=0;
209  }
210 
211  int indexDisk[10];
212  for(int j=0;j<10;j++){
213  indexDisk[j]=0;
214  }
215 
216 
217  for(TrackingGeometry::DetContainer::const_iterator it=rpcGeo->dets().begin(); it!=rpcGeo->dets().end();it++){
218 
219 
220 
221  if(dynamic_cast< const RPCChamber* >( *it ) != 0 ){
222  const RPCChamber* ch = dynamic_cast< const RPCChamber* >( *it );
223  std::vector< const RPCRoll*> roles = (ch->rolls());
224  for(std::vector<const RPCRoll*>::const_iterator r = roles.begin();r != roles.end(); ++r){
225  RPCDetId rpcId = (*r)->id();
226  RPCGeomServ rpcsrv(rpcId);
227  // std::string nameRoll = rpcsrv.name();
228  std::string camera = rpcsrv.name();
229 
230  //Breaking down the geometry
231  int region = rpcId.region();
232  int wheel = rpcId.ring(); int ring = rpcId.ring();
233  int sector = rpcId.sector();
234  int station = rpcId.station();
235  int geolayer = rpcId.layer();
236 
237 
238  if(region==0){//Barrel
239 
240  std::stringstream meIdRPC, meIdDT; //, bxDistroId;
241 
242  meIdRPC<<folderPath<<"MuonSegEff/RPCDataOccupancyFromDT_"<<rpcId.rawId();
243  meIdDT<<folderPath<<"MuonSegEff/ExpectedOccupancyFromDT_"<<rpcId.rawId();
244 
245  histoRPC = dbe->get(meIdRPC.str());
246  histoDT = dbe->get(meIdDT.str());
247 
248  int NumberWithOutPrediction=0;
249  double p = 0.;
250  double o = 0.;
251  double ef = 0.;
252  double er = 0.;
253  double buffef = 0.;
254  double buffer = 0.;
255  double sumbuffef = 0.;
256  double sumbuffer = 0.;
257  double averageeff = 0.;
258  int NumberStripsPointed = 0;
259 
260 
261  if(histoRPC && histoDT){
262 
263  for(int i=1;i<=int((*r)->nstrips());++i){
264 
265  if(histoDT->getBinContent(i)!=0){
266  LogDebug("rpcefficiencysecond")<<"Inside the If";
267 
268  buffef = double(histoRPC->getBinContent(i))/double(histoDT->getBinContent(i));
269  buffer = sqrt(buffef*(1.-buffef)/double(histoDT->getBinContent(i)));
270 
271  sumbuffef=sumbuffef+buffef;
272  sumbuffer = sumbuffer + buffer*buffer;
273  NumberStripsPointed++;
274  }else{
275  NumberWithOutPrediction++;
276  }
277  LogDebug("rpcefficiencysecond")<<"Strip="<<i<<" RPC="<<histoRPC->getBinContent(i)<<" DT="<<histoDT->getBinContent(i)<<" buffef="<<buffef<<" buffer="<<buffer<<" sumbuffef="<<sumbuffef<<" sumbuffer="<<sumbuffer<<" NumberStripsPointed="<<NumberStripsPointed<<" NumberWithOutPrediction"<<NumberWithOutPrediction;
278  }
279 
280  p=histoDT->getTH1F()->Integral();
281  o=histoRPC->getTH1F()->Integral();
282 
283  if(NumberStripsPointed!=0){
284  averageeff = (sumbuffef/double(NumberStripsPointed))*100.;
285  //averageerr = sqrt(sumbuffer/double(NumberStripsPointed))*100.;
286  }
287 
288  }
289 
290  // int Ring = rpcId.ring();
291 
292  if(p!=0){
293  ef = double(o)/double(p);
294  er = sqrt(ef*(1.-ef)/double(p));
295  }
296 
297  ef=ef*100;
298  er=er*100;
299 
300 
301  //Filling azimutal Wheel Histograms
302 
303  int layer = 0;
304 
305  if(station==1&&geolayer==1) layer = 1;
306  else if(station==1&&geolayer==2) layer = 2;
307  else if(station==2&&geolayer==1) layer = 3;
308  else if(station==2&&geolayer==2) layer = 4;
309  else if(station==3) layer = 5;
310  else if(station==4) layer = 6;
311 
312 
313  ExsectorEffW[wheel+2]->Fill(sector,p); OcsectorEffW[wheel+2]->Fill(sector,o);
314  ExpLayerW[wheel+2]->Fill(layer, p); ObsLayerW[wheel+2]->Fill(layer, o);
315 
316 
317 
318  LogDebug("rpcefficiencysecond")<<"Pigi "<<camera<<" "<<rpcsrv.shortname()<<" "<<(*r)->id();
319 
320  if(p > 100){//We need at least 100 predictions to fill the summary plot
321 
322  int xBin,yBin;
323  xBin= (*r)->id().sector();
324  rpcdqm::utils rollNumber;
325  yBin = rollNumber.detId2RollNr((*r)->id());
326 
327 
328  WheelSummary[wheel+2]->setBinContent(xBin,yBin,averageeff);
329  }
330 
331 
332  EffDistroW[wheel+2]->Fill(averageeff);
333  indexWheel[wheel+2]++;
334  EffGlobW[wheel+2]->setBinContent(indexWheel[wheel+2],ef);
335  EffGlobW[wheel+2]->setBinError(indexWheel[wheel+2],er);
336  EffGlobW[wheel+2]->setBinLabel(indexWheel[wheel+2],camera,1);
337 
338  }else{//EndCap
339 
340 
341  std::stringstream meIdRPC,meIdCSC; //, bxDistroId;
342  std::string meIdPRO;
343 
344 
345  meIdRPC<<folderPath<<"MuonSegEff/RPCDataOccupancyFromCSC_"<<rpcId.rawId();
346  meIdCSC<<folderPath<<"MuonSegEff/ExpectedOccupancyFromCSC_"<<rpcId.rawId();
347 
348  meIdPRO = "Profile_"+ rpcsrv.name();
349 
350  histoRPC= dbe->get(meIdRPC.str());
351  histoCSC= dbe->get(meIdCSC.str());
352  //BXDistribution = dbe->get(bxDistroId.str());
353 
354  int NumberWithOutPrediction=0;
355  double p = 0;
356  double o = 0;
357  double ef =0;
358  double er =0;
359  double buffef = 0;
360  double buffer = 0;
361  double sumbuffef = 0;
362  double sumbuffer = 0;
363  double averageeff = 0;
364  int NumberStripsPointed = 0;
365 
366 
367  if(histoRPC && histoCSC) {
368  LogDebug("rpcefficiencysecond")<<rpcsrv.name();
369 
370  for(int i=1;i<=int((*r)->nstrips());++i){
371  if(histoCSC->getBinContent(i)!=0){
372  LogDebug("rpcefficiencysecond")<<"Inside the If";
373  buffef = double(histoRPC->getBinContent(i))/double(histoCSC->getBinContent(i));
374  // meMap[meIdPRO]->setBinContent(i,buffef);
375  buffer = sqrt(buffef*(1.-buffef)/double(histoCSC->getBinContent(i)));
376  // meMap[meIdPRO]->setBinError(i,buffer);
377  sumbuffef=sumbuffef+buffef;
378  sumbuffer = sumbuffer + buffer*buffer;
379  NumberStripsPointed++;
380  }else{
381  NumberWithOutPrediction++;
382  }
383 
384  LogDebug("rpcefficiencysecond")<<"Strip="<<i<<" RPC="<<histoRPC->getBinContent(i)<<" CSC="<<histoCSC->getBinContent(i)<<" buffef="<<buffef<<" buffer="<<buffer<<" sumbuffef="<<sumbuffef<<" sumbuffer="<<sumbuffer<<" NumberStripsPointed="<<NumberStripsPointed<<" NumberWithOutPrediction"<<NumberWithOutPrediction;
385  }
386 
387  p=histoCSC->getTH1F()->Integral();
388  o=histoRPC->getTH1F()->Integral();
389 
390  if(NumberStripsPointed!=0){
391  averageeff = (sumbuffef/double(NumberStripsPointed))*100.;
392  }
393 
394  }
395 
396  int Disk = station*region;
397 
398  if(Disk > numberOfDisks_ || Disk<-numberOfDisks_){continue;} //remove strange disk numbers!!
399 
400  int dIndex = -1;
401  if(region == -1){
402  dIndex = Disk+numberOfDisks_;
403  }else if(region == 1){
404  dIndex = Disk+numberOfDisks_-1;
405  }
406 
407  if( dIndex<0 || dIndex>= numberOfDisks_*2){continue;} //extra check on disk numering
408 
409  if(p!=0){
410  ef = double(o)/double(p);
411  er = sqrt(ef*(1.-ef)/double(p));
412  }
413 
414  ef=ef*100;
415  er=er*100;
416 
417  if(ring==2){
418  if (ExGregR2D[dIndex] && OcGregR2D[dIndex] ){
419  ExGregR2D[dIndex]->Fill(rpcsrv.segment(),p);
420  OcGregR2D[dIndex]->Fill(rpcsrv.segment(),o);
421  }
422  } else if(ring==3){
423  ExGregR3D[dIndex]->Fill(rpcsrv.segment(),p);
424  OcGregR3D[dIndex]->Fill(rpcsrv.segment(),o);
425  }
426 
427  if(p > 100){ //We need at least 100 predictions to fill the summary plot
428  RPCGeomServ RPCServ((*r)->id());
429  int xBin = RPCServ.segment();
430  int yBin= (ring-1)*3-(*r)->id().roll()+1;
431 
432  DiskSummary[dIndex]->setBinContent(xBin, yBin, averageeff);
433  }
434 
435 
436  if ( EffDistroD[dIndex]){EffDistroD[dIndex]->Fill(averageeff);}
437  indexDisk[dIndex]++;
438  if ( EffGlobD[dIndex]){
439  EffGlobD[dIndex]->setBinContent(indexDisk[dIndex],ef);
440  EffGlobD[dIndex]->setBinError(indexDisk[dIndex],er);
441  EffGlobD[dIndex]->setBinLabel(indexDisk[dIndex],camera,1);
442 
443  }
444  }
445  }
446  }
447  }
448 
449  double eff,N,err;
450  int k;
451 
452  for (int d = 0; d<(numberOfDisks_*2); d++){
453  for(k=1;k<=36;k++){
454 
455  err=0; eff=0; N=ExGregR2D[d]->getBinContent(k);
456  if(N!=0.){ eff = OcGregR2D[d]->getBinContent(k)/N; err=sqrt(eff*(1-eff)/N);}
457  GregR2D[d]->setBinContent(k,eff); GregR2D[d]->setBinError(k,err);
458 
459  err=0; eff=0; N=ExGregR3D[d]->getBinContent(k);
460  if(N!=0.){ eff = OcGregR3D[d]->getBinContent(k)/N; err=sqrt(eff*(1-eff)/N);}
461  GregR3D[d]->setBinContent(k,eff); GregR3D[d]->setBinError(k,err);
462  }
463  }
464 
465 
466  for (int w =-2; w<=2; w++){
467  for(k=1;k<=12;k++){
468  err=0; eff=0; N=ExsectorEffW[w+2]->getBinContent(k);
469  if(N!=0.){ eff = OcsectorEffW[w+2]->getBinContent(k)/N;err=sqrt(eff*(1-eff)/N);}
470  sectorEffW[w+2]->setBinContent(k,eff); sectorEffW[w+2]->setBinError(k,err);
471  }
472  }
473 
474 
475  //Ranges for Both
476  //Barrel
477 
478  for (int w=-2; w<=2; w++){
479  EffGlobW[w+2]->setAxisRange(-4.,100.,2);
480  EffGlobW[w+2]->setAxisTitle("%",2);
481  }
482 
483  for (int d=0; d<(numberOfDisks_*2); d++){
484  EffGlobD[d]->setAxisRange(-4.,100.,2);
485  EffGlobD[d]->setAxisTitle("%",2);
486  }
487 
488 
489  if(SaveFile){
490  dbe->save(NameFile);
491  }
492 
493 }
494 
496 
#define LogDebug(id)
MonitorElement * EffGlobW[5]
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
virtual void analyze(const edm::Event &, const edm::EventSetup &)
void setBinContent(int binx, double content)
set content of bin (1-D)
RPCEfficiencySecond(const edm::ParameterSet &)
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:942
MonitorElement * histoCSC
void setAxisRange(double xmin, double xmax, int axis=1)
set x-, y- or z-axis range (axis=1, 2, 3 respectively)
tuple yBin
Definition: cuy.py:891
MonitorElement * OcGregR3D[10]
MonitorElement * histoDT
void setBinLabel(int bin, const std::string &label, int axis=1)
set bin label for x, y or z axis (axis=1, 2, 3 respectively)
MonitorElement * OcsectorEffW[5]
void labelXAxisSegment(MonitorElement *myMe)
Definition: utils.h:250
void labelYAxisRoll(MonitorElement *myMe, int region, int ring, bool useRollInfo)
Definition: utils.h:268
MonitorElement * sectorEffW[5]
void Fill(long long x)
uint32_t rawId() const
get the raw id
Definition: DetId.h:43
virtual void endRun(const edm::Run &, const edm::EventSetup &)
int iEvent
Definition: GenABIO.cc:230
virtual std::string name()
Definition: RPCGeomServ.cc:20
void labelYAxisRing(MonitorElement *myMe, int numberOfRings, bool useRollInfo)
Definition: utils.h:295
MonitorElement * ObsLayerW[5]
int ring() const
Definition: RPCDetId.h:72
T sqrt(T t)
Definition: SSEVec.h:48
MonitorElement * EffGlobD[10]
const std::vector< const RPCRoll * > & rolls() const
Return the Rolls.
Definition: RPCChamber.cc:68
int j
Definition: DBlmapReader.cc:9
MonitorElement * EffDistroD[10]
MonitorElement * ExpLayerW[5]
void setBinError(int binx, double error)
set uncertainty on content of bin (1-D)
MonitorElement * OcGregR2D[10]
MonitorElement * get(const std::string &path) const
get ME from full pathname (e.g. &quot;my/long/dir/my_histo&quot;)
Definition: DQMStore.cc:1696
virtual std::string shortname()
Definition: RPCGeomServ.cc:206
void labelXAxisSector(MonitorElement *myMe)
Definition: utils.h:236
int k[5][pyjets_maxn]
int layer() const
Definition: RPCDetId.h:108
virtual int segment()
Definition: RPCGeomServ.cc:467
int detId2RollNr(const RPCDetId &_id)
Definition: utils.h:18
#define N
Definition: blowfish.cc:9
MonitorElement * WheelSummary[5]
const T & get() const
Definition: EventSetup.h:55
TH1F * getTH1F(void) const
MonitorElement * ExGregR2D[10]
int sector() const
Sector id: the group of chambers at same phi (and increasing r)
Definition: RPCDetId.h:102
MonitorElement * DiskSummary[10]
virtual void beginRun(const edm::Run &, const edm::EventSetup &iSetup)
double getBinContent(int binx) const
get content of bin (1-D)
void save(const std::string &filename, const std::string &path="", const std::string &pattern="", const std::string &rewrite="", const uint32_t run=0, const uint32_t lumi=0, SaveReferenceTag ref=SaveWithReference, int minStatus=dqm::qstatus::STATUS_OK, const std::string &fileupdate="RECREATE", const bool resetMEsAfterWriting=false)
Definition: DQMStore.cc:2490
MonitorElement * ExsectorEffW[5]
T w() const
Definition: DDAxes.h:10
MonitorElement * GregR2D[10]
MonitorElement * ExGregR3D[10]
MonitorElement * book2D(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Book 2D histogram.
Definition: DQMStore.cc:1070
void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
MonitorElement * histoRPC
MonitorElement * GregR3D[10]
MonitorElement * EffDistroW[5]
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:655
Definition: Run.h:41
int region() const
Region id: 0 for Barrel, +/-1 For +/- Endcap.
Definition: RPCDetId.h:63
int station() const
Definition: RPCDetId.h:96