CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Types | Private Member Functions | Private Attributes | Static Private Attributes
RPCChamberQuality Class Reference

#include <RPCChamberQuality.h>

Inheritance diagram for RPCChamberQuality:
edm::EDAnalyzer

Public Member Functions

void analyze (const edm::Event &iEvent, const edm::EventSetup &c)
 
void beginJob ()
 
void beginLuminosityBlock (edm::LuminosityBlock const &lumiSeg, edm::EventSetup const &context)
 
void beginRun (const edm::Run &r, const edm::EventSetup &c)
 
void endLuminosityBlock (edm::LuminosityBlock const &lumiSeg, edm::EventSetup const &c)
 
void endRun (const edm::Run &r, const edm::EventSetup &c)
 
 RPCChamberQuality (const edm::ParameterSet &ps)
 
virtual ~RPCChamberQuality ()
 
- Public Member Functions inherited from edm::EDAnalyzer
 EDAnalyzer ()
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 

Private Types

enum  chamberQualityState {
  GoodState = 1, OffState =2, NoisyStripState = 3, NoisyRollState = 4,
  PartiallyDeadState =5, DeadState =6, BadShapeState =7
}
 

Private Member Functions

void performeClientOperation (std::string, int, MonitorElement *)
 

Private Attributes

DQMStoredbe_
 
int minEvents
 
int numberOfDisks_
 
int numLumBlock_
 
std::string prefixDir_
 
int prescaleFactor_
 

Static Private Attributes

static const std::string regions_ [3] = {"EndcapNegative","Barrel","EndcapPositive"}
 
static const std::string xLabels_ [7] = {"Good", "OFF", "Nois.St","Nois.Ch","Part.Dead","Dead","Bad.Shape"}
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
typedef WorkerT< EDAnalyzerWorkerType
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDAnalyzer
CurrentProcessingContext const * currentContext () const
 

Detailed Description

Definition at line 19 of file RPCChamberQuality.h.

Member Enumeration Documentation

Constructor & Destructor Documentation

RPCChamberQuality::RPCChamberQuality ( const edm::ParameterSet ps)

Definition at line 24 of file RPCChamberQuality.cc.

References edm::ParameterSet::getUntrackedParameter(), minEvents, numberOfDisks_, prefixDir_, and prescaleFactor_.

24  {
25  edm::LogVerbatim ("rpceventsummary") << "[RPCChamberQuality]: Constructor";
26 
27  prescaleFactor_ = ps.getUntrackedParameter<int>("PrescaleFactor", 9);
28  prefixDir_ = ps.getUntrackedParameter<std::string>("RPCGlobalFolder", "RPC/RecHits/SummaryHistograms");
29  minEvents = ps.getUntrackedParameter<int>("MinimumRPCEvents", 10000);
30  numberOfDisks_ = ps.getUntrackedParameter<int>("NumberOfEndcapDisks", 3);
31 }
T getUntrackedParameter(std::string const &, T const &) const
std::string prefixDir_
RPCChamberQuality::~RPCChamberQuality ( )
virtual

Definition at line 33 of file RPCChamberQuality.cc.

References dbe_.

33  {
34  edm::LogVerbatim ("rpceventsummary") << "[RPCChamberQuality]: Destructor ";
35  dbe_=0;
36 }

Member Function Documentation

void RPCChamberQuality::analyze ( const edm::Event iEvent,
const edm::EventSetup c 
)
virtual

Implements edm::EDAnalyzer.

Definition at line 125 of file RPCChamberQuality.cc.

125 {}
void RPCChamberQuality::beginJob ( void  )
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 38 of file RPCChamberQuality.cc.

References dbe_, and cmsCodeRules.cppFunctionSkipper::operator.

38  {
39  edm::LogVerbatim ("rpceventsummary") << "[RPCChamberQuality]: Begin job ";
41 }
void RPCChamberQuality::beginLuminosityBlock ( edm::LuminosityBlock const &  lumiSeg,
edm::EventSetup const &  context 
)
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 123 of file RPCChamberQuality.cc.

123 {}
void RPCChamberQuality::beginRun ( const edm::Run r,
const edm::EventSetup c 
)
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 43 of file RPCChamberQuality.cc.

References DQMStore::book1D(), DQMStore::book2D(), dbe_, DQMStore::get(), MonitorElement::getName(), rpcdqm::utils::labelXAxisSector(), rpcdqm::utils::labelXAxisSegment(), rpcdqm::utils::labelYAxisRing(), rpcdqm::utils::labelYAxisRoll(), numberOfDisks_, prefixDir_, csvReporter::r, regions_, DQMStore::removeElement(), MonitorElement::setBinLabel(), DQMStore::setCurrentFolder(), ExpressReco_HICollisions_FallBack::x, and xLabels_.

43  {
44  edm::LogVerbatim ("rpceventsummary") << "[RPCChamberQuality]: Begin run";
45 
46  // init_ = false;
47 
48  MonitorElement* me;
50 
51  std::stringstream histoName;
52 
53  rpcdqm::utils rpcUtils;
54 
55  for (int r = 0 ; r < 3; r++){
56 
57  histoName.str("");
58  histoName<<"RPCChamberQuality_"<<regions_[r];
59  me = dbe_->get(prefixDir_+"/"+ histoName.str());
60  if (0!=me) dbe_->removeElement(me->getName());
61  me = dbe_->book1D(histoName.str().c_str(), histoName.str().c_str(), 7, 0.5, 7.5);
62 
63  for (int x = 1; x <8 ; x++) me->setBinLabel(x, xLabels_[x-1]);
64  }
65 
66 
67  histoName.str("");
68  histoName<<"RPC_System_Quality_Overview";
69  me = dbe_->get(prefixDir_+"/"+ histoName.str());
70  if (0!=me) dbe_->removeElement(me->getName());
71  me = dbe_->book2D(histoName.str().c_str(), histoName.str().c_str(), 7, 0.5, 7.5, 3, 0.5, 3.5);
72  me->setBinLabel(1, "E+", 2);
73  me->setBinLabel(2, "B", 2);
74  me->setBinLabel(3, "E-", 2);
75 
76  for (int x = 1; x <8 ; x++) me->setBinLabel(x, xLabels_[x-1]);
77 
78  for(int w=-2; w<3;w++){//Loop on wheels
79 
80  histoName.str("");
81  histoName<<"RPCChamberQuality_Roll_vs_Sector_Wheel"<<w;
82  me = dbe_->get(prefixDir_+"/"+ histoName.str());
83  if (0!=me) dbe_->removeElement(me->getName());
84  me = dbe_->book2D(histoName.str().c_str(), histoName.str().c_str(), 12, 0.5, 12.5, 21, 0.5, 21.5);
85 
86  rpcUtils.labelXAxisSector( me);
87  rpcUtils.labelYAxisRoll(me, 0, w);
88 
89  histoName.str("");
90  histoName<<"RPCChamberQuality_Distribution_Wheel"<<w;
91  me=0;
92  me = dbe_->get(prefixDir_+"/"+ histoName.str());
93  if (0!=me ) dbe_->removeElement(me->getName());
94  me = dbe_->book1D(histoName.str().c_str(), histoName.str().c_str(), 7, 0.5, 7.5);
95 
96  for (int x = 1; x <8; x++) me->setBinLabel(x, xLabels_[x-1]);
97  }//end loop on wheels
98 
99  for(int d= -numberOfDisks_; d<= numberOfDisks_ ; d++) { // Loop on disk
100  if(d==0) continue;
101  histoName.str("");
102  histoName<<"RPCChamberQuality_Ring_vs_Segment_Disk"<<d; // 2D histo for RPC Qtest
103  me = 0;
104  me = dbe_->get(prefixDir_+"/"+ histoName.str());
105  if (0!=me) {
106  dbe_->removeElement(me->getName());
107  }
108  me = dbe_->book2D(histoName.str().c_str(), histoName.str().c_str(), 36, 0.5, 36.5, 6, 0.5, 6.5);
109  rpcUtils.labelXAxisSegment(me);
110  rpcUtils.labelYAxisRing(me, 2);
111 
112  histoName.str("");
113  histoName<<"RPCChamberQuality_Distribution_Disk"<<d;
114  me=0;
115  me = dbe_->get(prefixDir_+"/"+ histoName.str());
116  if (0!=me ) dbe_->removeElement(me->getName());
117  me = dbe_->book1D(histoName.str().c_str(), histoName.str().c_str(), 7, 0.5, 7.5);
118 
119  for (int x = 1; x <8 ; x++) me->setBinLabel(x, xLabels_[x-1]);
120  }
121 }
const std::string & getName(void) const
get name of ME
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:519
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)
void labelXAxisSegment(MonitorElement *myMe)
Definition: utils.h:221
static const std::string xLabels_[7]
void removeElement(const std::string &name)
Definition: DQMStore.cc:2338
MonitorElement * get(const std::string &path) const
get ME from full pathname (e.g. &quot;my/long/dir/my_histo&quot;)
Definition: DQMStore.cc:1270
void labelXAxisSector(MonitorElement *myMe)
Definition: utils.h:206
static const std::string regions_[3]
std::string prefixDir_
void labelYAxisRoll(MonitorElement *myMe, int region, int ring)
Definition: utils.h:239
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:647
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:237
void labelYAxisRing(MonitorElement *myMe, int numberOfRings)
Definition: utils.h:266
void RPCChamberQuality::endLuminosityBlock ( edm::LuminosityBlock const &  lumiSeg,
edm::EventSetup const &  c 
)
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 192 of file RPCChamberQuality.cc.

192 { }
void RPCChamberQuality::endRun ( const edm::Run r,
const edm::EventSetup c 
)
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 127 of file RPCChamberQuality.cc.

References dbe_, python.tagInventory::entries, DQMStore::get(), MonitorElement::getEntries(), i, minEvents, NULL, performeClientOperation(), prefixDir_, csvReporter::r, regions_, MonitorElement::Reset(), MonitorElement::setBinContent(), edmLumisInFiles::summary, and ExpressReco_HICollisions_FallBack::x.

127  {
128  edm::LogVerbatim ("rpceventsummary") <<"[RPCChamberQuality]: End Job, performing DQM client operation";
129 
130  MonitorElement * RpcEvents = NULL;
131  std::stringstream meName;
132 
133 
134  meName.str("");
135  meName<<prefixDir_<<"/RPCEvents";
136  int rpcEvents=0;
137  RpcEvents = dbe_->get(meName.str());
138 
139  if(RpcEvents) rpcEvents= (int)RpcEvents->getEntries();
140 
141  if(rpcEvents >= minEvents){
142 
143  MonitorElement * summary[3];
144 
145  for(int r = 0 ; r < 3 ; r++) {
146  meName.str("");
147  meName<<prefixDir_<<"/RPCChamberQuality_"<<RPCChamberQuality::regions_[r];
148  summary[r] = dbe_ -> get(meName.str());
149 
150  if( summary[r] != 0 ) summary[r]->Reset();
151  }
152 
153  //Barrel
154  for (int wheel=-2; wheel<3; wheel++) { // loop by Wheels
155  meName.str("");
156  meName<<"Roll_vs_Sector_Wheel"<<wheel;
157 
158  this->performeClientOperation(meName.str(), 0 , summary[1]);
159  } // loop by Wheels
160 
161 
162  // Endcap
163  for(int i=-3; i<4; i++) {//loop on Disks
164  if(i==0) continue;
165 
166  meName.str("");
167  meName<<"Ring_vs_Segment_Disk"<<i;
168 
169  if(i<0) this->performeClientOperation(meName.str(), -1 , summary[0]);
170  else this->performeClientOperation(meName.str(), 1 , summary[2]);
171  }//loop on Disks
172 
173  MonitorElement * RpcOverview = NULL;
174  meName.str("");
175  meName<<prefixDir_<<"/RPC_System_Quality_Overview";
176  RpcOverview = dbe_ -> get(meName.str());
177  RpcOverview->Reset();
178 
179  if(RpcOverview) {//Fill Overview ME
180  for(int r = 0 ; r< 3; r++) {
181  if (summary[r] == 0 ) continue;
182  double entries = summary[r]->getEntries();
183  if(entries == 0) continue;
184  for (int x = 1; x <= 7; x++) {
185  RpcOverview->setBinContent(x,r+1,(summary[r]->getBinContent(x)/entries));
186  }
187  }
188  } //loop by LimiBloks
189  }
190 }
int i
Definition: DBlmapReader.cc:9
void setBinContent(int binx, double content)
set content of bin (1-D)
#define NULL
Definition: scimark2.h:8
double getEntries(void) const
get # of entries
void performeClientOperation(std::string, int, MonitorElement *)
MonitorElement * get(const std::string &path) const
get ME from full pathname (e.g. &quot;my/long/dir/my_histo&quot;)
Definition: DQMStore.cc:1270
static const std::string regions_[3]
std::string prefixDir_
void Reset(void)
reset ME (ie. contents, errors, etc)
void RPCChamberQuality::performeClientOperation ( std::string  MESufix,
int  region,
MonitorElement quality 
)
private

Definition at line 195 of file RPCChamberQuality.cc.

References BadShapeState, dbe_, DeadState, MonitorElement::Fill(), DQMStore::get(), MonitorElement::getBinContent(), GoodState, VarParsing::mult, NoisyRollState, NoisyStripState, NULL, OffState, PartiallyDeadState, pos, prefixDir_, MonitorElement::Reset(), ExpressReco_HICollisions_FallBack::x, and ExpressReco_HICollisions_FallBack::y.

Referenced by endRun().

195  {
196 
197 
198  MonitorElement * RCQ=NULL;
199  MonitorElement * RCQD=NULL;
200 
201  MonitorElement * DEAD=NULL;
202  MonitorElement * CLS=NULL;
203  MonitorElement * MULT=NULL;
204  MonitorElement * NoisySt=NULL;
205  MonitorElement * Chip=NULL;
206  MonitorElement * HV=NULL;
208  std::stringstream meName;
209 
210  meName.str("");
211  meName<<prefixDir_<<"/RPCChamberQuality_"<<MESufix;
212  RCQ = dbe_ -> get(meName.str());
213  if (RCQ) RCQ->Reset();
214 
215 
216  int pos = MESufix.find_last_of("_");
217  meName.str("");
218  meName<<prefixDir_<<"/RPCChamberQuality_Distribution"<<MESufix.substr(pos);
219  RCQD = dbe_ -> get(meName.str());
220  if (RCQD) RCQD->Reset();
221 
222  //get HV Histo
223  meName.str("");
224  meName<<prefixDir_<<"/HVStatus_"<<MESufix;
225  HV = dbe_ -> get(meName.str());
226  //get LV Histo
227  meName.str("");
228  meName<<prefixDir_<<"/LVStatus_"<<MESufix;
229  LV = dbe_ -> get(meName.str());
230  //Dead
231  meName.str("");
232  meName << prefixDir_<<"/DeadChannelFraction_"<<MESufix;
233  DEAD = dbe_->get(meName.str());
234  //ClusterSize
235  meName.str("");
236  meName<<prefixDir_<<"/ClusterSizeIn1Bin_"<<MESufix;
237  CLS = dbe_ -> get(meName.str());
238  //NoisyStrips
239  meName.str("");
240  meName<<prefixDir_<<"/RPCNoisyStrips_"<<MESufix;
241  NoisySt = dbe_ -> get(meName.str());
242  //Multiplicity
243  meName.str("");
244  meName<<prefixDir_<<"/NumberOfDigi_Mean_"<<MESufix;
245  MULT = dbe_ -> get(meName.str());
246  //Asymetry
247  meName.str("");
248  meName<<prefixDir_<<"/AsymmetryLeftRight_"<<MESufix;
249  Chip = dbe_ -> get(meName.str());
250 
251  int xBinMax, yBinMax;
252 
253  if (region != 0) xBinMax = 37;
254  else xBinMax = 13;
255 
256  for(int x=1; x<xBinMax; x++) {
257  if (region != 0 ) {
258  yBinMax = 7;
259  }else {
260  if(x==4) yBinMax=22;
261  else if(x==9 || x==11) yBinMax=16;
262  else yBinMax=18;
263  }
264  for(int y=1; y<yBinMax; y++) {
265  int hv=1;
266  int lv=1;
267  float dead =0;
268  float firstbin= 0;
269  float noisystrips = 0;
270  float mult = 0;
271  float asy = 0;
272  chamberQualityState chamberState = GoodState;
273 
274  if(HV) hv = (int)HV ->getBinContent(x,y);
275  if(LV) lv = (int)LV ->getBinContent(x,y);
276 
277  if( hv!=1 || lv!=1) {
278  chamberState = OffState;
279  }else {
280  if(DEAD) dead= DEAD -> getBinContent(x,y);
281  if(dead>=0.80 ) {
282  chamberState = DeadState;
283  }else if (0.33<=dead && dead<0.80 ){
284  chamberState = PartiallyDeadState;
285  }else {
286  if(CLS ) firstbin = CLS -> getBinContent(x,y);
287  if(firstbin >= 0.88) {
288  chamberState = NoisyStripState;
289  } else {
290  if(NoisySt) noisystrips = NoisySt -> getBinContent(x,y);
291  if (noisystrips > 0){
292  chamberState = NoisyStripState;
293  }else {
294  if(MULT) mult = MULT -> getBinContent(x,y);
295  if(mult>=6) {
296  chamberState = NoisyRollState;
297  }else {
298  if (Chip) asy = Chip->getBinContent(x,y);
299  if(asy>0.35) {
300  chamberState = BadShapeState;
301  }else {
302  chamberState = GoodState;
303  }
304  }
305  }
306  }
307  }
308  }
309  if (RCQ) RCQ -> setBinContent(x,y, chamberState);
310  if (RCQD) RCQD -> Fill(chamberState);
311  if (quality) quality->Fill(chamberState);
312  }
313  }
314  return;
315 }
#define NULL
Definition: scimark2.h:8
math::XYZTLorentzVectorD LV
void Fill(long long x)
MonitorElement * get(const std::string &path) const
get ME from full pathname (e.g. &quot;my/long/dir/my_histo&quot;)
Definition: DQMStore.cc:1270
std::string prefixDir_
double getBinContent(int binx) const
get content of bin (1-D)
void Reset(void)
reset ME (ie. contents, errors, etc)

Member Data Documentation

DQMStore* RPCChamberQuality::dbe_
private
int RPCChamberQuality::minEvents
private

Definition at line 51 of file RPCChamberQuality.h.

Referenced by endRun(), and RPCChamberQuality().

int RPCChamberQuality::numberOfDisks_
private

Definition at line 46 of file RPCChamberQuality.h.

Referenced by beginRun(), and RPCChamberQuality().

int RPCChamberQuality::numLumBlock_
private

Definition at line 52 of file RPCChamberQuality.h.

std::string RPCChamberQuality::prefixDir_
private

Definition at line 41 of file RPCChamberQuality.h.

Referenced by beginRun(), endRun(), performeClientOperation(), and RPCChamberQuality().

int RPCChamberQuality::prescaleFactor_
private

Definition at line 45 of file RPCChamberQuality.h.

Referenced by RPCChamberQuality().

const std::string RPCChamberQuality::regions_ = {"EndcapNegative","Barrel","EndcapPositive"}
staticprivate

Definition at line 43 of file RPCChamberQuality.h.

Referenced by beginRun(), and endRun().

const std::string RPCChamberQuality::xLabels_ = {"Good", "OFF", "Nois.St","Nois.Ch","Part.Dead","Dead","Bad.Shape"}
staticprivate

Definition at line 42 of file RPCChamberQuality.h.

Referenced by beginRun().