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 Member Functions | Private Attributes
DTNoiseCalibration Class Reference

#include <DTNoiseCalibration.h>

Inheritance diagram for DTNoiseCalibration:
edm::EDAnalyzer

Public Member Functions

void analyze (const edm::Event &e, const edm::EventSetup &c)
 Analyze. More...
 
virtual void beginRun (const edm::Run &run, const edm::EventSetup &setup)
 
 DTNoiseCalibration (const edm::ParameterSet &ps)
 Constructor. More...
 
void endJob ()
 Endjob. More...
 
virtual ~DTNoiseCalibration ()
 Destructor. More...
 
- Public Member Functions inherited from edm::EDAnalyzer
 EDAnalyzer ()
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 

Private Member Functions

std::string getLayerName (const DTLayerId &lId) const
 Get the name of the layer. More...
 
std::string getSuperLayerName (const DTSuperLayerId &dtSLId) const
 Get the name of the superLayer. More...
 

Private Attributes

bool cosmicRun
 
int counter
 
std::string dbLabel
 
bool debug
 
std::string digiLabel
 
edm::ESHandle< DTGeometrydtGeom
 
bool fastAnalysis
 
TH1F * hTDCTriggerWidth
 
float kFactor
 
int nevents
 
edm::ParameterSet parameters
 
int sect
 
std::map< DTLayerId, int > skippedPlot
 
TFile * theFile
 
std::map< DTLayerId, TH2F * > theHistoEvtPerWireMap
 
std::map< DTLayerId, TH1F * > theHistoOccupancyMap
 
double theOffset
 
int TriggerWidth
 variables to set by configuration file More...
 
float tTrig
 tTrig from the DB More...
 
edm::ESHandle< DTTtrigtTrigMap
 
float tTrigRMS
 
float upperLimit
 
int wh
 

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 33 of file DTNoiseCalibration.h.

Constructor & Destructor Documentation

DTNoiseCalibration::DTNoiseCalibration ( const edm::ParameterSet ps)

Constructor.

Definition at line 44 of file DTNoiseCalibration.cc.

References gather_cfg::cout, debug, ExpressReco_HICollisions_FallBack::digiLabel, edm::ParameterSet::getUntrackedParameter(), ExpressReco_HICollisions_FallBack::parameters, and interactiveExample::theFile.

44  {
45 
46  cout << "[DTNoiseCalibration]: Constructor" <<endl;
47 
48  // Get the debug parameter for verbose output
49  debug = ps.getUntrackedParameter<bool>("debug");
50 
51  // Get the label to retrieve digis from the event
52  digiLabel = ps.getUntrackedParameter<string>("digiLabel");
53 
54  // The analysis type
55  fastAnalysis = ps.getUntrackedParameter<bool>("fastAnalysis", true);
56  // The wheel & sector interested for the time-dependent analysis
57  wh = ps.getUntrackedParameter<int>("wheel", 0);
58  sect = ps.getUntrackedParameter<int>("sector", 6);
59 
60  // The trigger mode
61  cosmicRun = ps.getUntrackedParameter<bool>("cosmicRun", false);
62 
63  // The trigger width (if noise run)
64  TriggerWidth = ps.getUntrackedParameter<int>("TriggerWidth");
65 
66  //get the offset to look for the noise
67  theOffset = ps.getUntrackedParameter<double>("theOffset",500.);
68 
69  // The root file which will contain the histos
70  string rootFileName = ps.getUntrackedParameter<string>("rootFileName");
71  theFile = new TFile(rootFileName.c_str(), "RECREATE");
72  theFile->cd();
73 
74  dbLabel = ps.getUntrackedParameter<string>("dbLabel", "");
75 
76  parameters=ps;
77 
78 }
T getUntrackedParameter(std::string const &, T const &) const
edm::ParameterSet parameters
int TriggerWidth
variables to set by configuration file
tuple cout
Definition: gather_cfg.py:41
DTNoiseCalibration::~DTNoiseCalibration ( )
virtual

Destructor.

Definition at line 403 of file DTNoiseCalibration.cc.

References gather_cfg::cout, nevents, and interactiveExample::theFile.

403  {
404 
405  cout << "DTNoiseCalibration: analyzed " << nevents << " events" <<endl;
406  theFile->Close();
407 
408 }
tuple cout
Definition: gather_cfg.py:41

Member Function Documentation

void DTNoiseCalibration::analyze ( const edm::Event e,
const edm::EventSetup c 
)
virtual

Analyze.

Implements edm::EDAnalyzer.

Definition at line 100 of file DTNoiseCalibration.cc.

References DTSuperLayerId::chamberId(), DTTopology::channels(), cmsDriverOptions::counter, DTTimeUnits::counts, gather_cfg::cout, debug, ExpressReco_HICollisions_FallBack::digiLabel, DTTopology::firstChannel(), edm::Event::getByLabel(), DTTopology::lastChannel(), nevents, ExpressReco_HICollisions_FallBack::parameters, edm::second(), DTChamberId::sector(), DTLayerId::superlayerId(), interactiveExample::theFile, and DTChamberId::wheel().

100  {
101  nevents++;
102  if(debug)
103  cout<<"nevents: "<<nevents<<endl;
104 
105  // Get the digis from the event
107  e.getByLabel(digiLabel, dtdigis);
108 
109  TH1F *hOccupancyHisto;
110  TH2F *hEvtPerWireH;
111  string Histo2Name;
112 
113  // LOOP OVER ALL THE DIGIS OF THE EVENT
115  for (dtLayerId_It=dtdigis->begin(); dtLayerId_It!=dtdigis->end(); ++dtLayerId_It){
116  for (DTDigiCollection::const_iterator digiIt = ((*dtLayerId_It).second).first;
117  digiIt!=((*dtLayerId_It).second).second; ++digiIt){
118 
119  //Check the TDC trigger width
120  int tdcTime = (*digiIt).countsTDC();
121  if(!cosmicRun){
122  if(debug)
123  cout<<"tdcTime (ns): "<<(tdcTime*25)/32<<endl;
124  if(((tdcTime*25)/32)>TriggerWidth){
125  cout<<"***Error*** : your digi has a tdcTime (ns) higher than the TDC trigger width :"<<(tdcTime*25)/32<<endl;
126  abort();
127  }
128  }
129 
130  if((!fastAnalysis &&
131  (*dtLayerId_It).first.superlayerId().chamberId().wheel()==wh &&
132  (*dtLayerId_It).first.superlayerId().chamberId().sector()==sect) ||
133  fastAnalysis)
134  hTDCTriggerWidth->Fill(tdcTime);
135 
136  // Set the window of interest if the run is triggered by cosmics
137  if ( parameters.getUntrackedParameter<bool>("readDB", true) ) {
138  tTrigMap->get( ((*dtLayerId_It).first).superlayerId(), tTrig, tTrigRMS, kFactor, DTTimeUnits::counts );
139 
141  }
142  else {
143  tTrig = parameters.getUntrackedParameter<int>("defaultTtrig", 4000);
145  }
146 
147  if((cosmicRun && (*digiIt).countsTDC()<upperLimit) || (!cosmicRun) ){
148 
149  if(debug && cosmicRun)
150  cout<<"tdcTime (ns): "<<((*digiIt).countsTDC()*25)/32<<" --- TriggerWidth (ns): "<<(upperLimit*25)/32<<endl;
151 
152  // Get the number of wires
153  const DTLayerId dtLId = (*dtLayerId_It).first;
154  const DTTopology& dtTopo = dtGeom->layer(dtLId)->specificTopology();
155  const int nWires = dtTopo.channels();
156  const int firstWire = dtTopo.firstChannel();
157  const int lastWire = dtTopo.lastChannel();
158 
159  // book the occupancy histos
160  theFile->cd();
161  if((!fastAnalysis &&
162  dtLId.superlayerId().chamberId().wheel()==wh &&
163  dtLId.superlayerId().chamberId().sector()==sect) ||
164  fastAnalysis){
165  hOccupancyHisto = theHistoOccupancyMap[dtLId];
166  if(hOccupancyHisto == 0) {
167  string HistoName = "DigiOccupancy_" + getLayerName(dtLId);
168  theFile->cd();
169  hOccupancyHisto = new TH1F(HistoName.c_str(), HistoName.c_str(), nWires, firstWire, lastWire+1);
170  if(debug)
171  cout << " New Occupancy Histo: " << hOccupancyHisto->GetName() << endl;
172  theHistoOccupancyMap[dtLId] = hOccupancyHisto;
173  }
174  hOccupancyHisto->Fill((*digiIt).wire());
175  }
176 
177  // book the digi event plot every 1000 events if the analysis is not "fast" and if is the correct sector
178  if(!fastAnalysis &&
179  dtLId.superlayerId().chamberId().wheel()==wh &&
180  dtLId.superlayerId().chamberId().sector()==sect) {
181  if(theHistoEvtPerWireMap.find(dtLId) == theHistoEvtPerWireMap.end() ||
182  (theHistoEvtPerWireMap.find(dtLId) != theHistoEvtPerWireMap.end() &&
183  skippedPlot[dtLId] != counter)){
184  skippedPlot[dtLId] = counter;
185  stringstream toAppend; toAppend << counter;
186  Histo2Name = "DigiPerWirePerEvent_" + getLayerName(dtLId) + "_" + toAppend.str();
187  theFile->cd();
188  hEvtPerWireH = new TH2F(Histo2Name.c_str(), Histo2Name.c_str(), 1000,0.5,1000.5,nWires, firstWire, lastWire+1);
189  if(hEvtPerWireH){
190  if(debug)
191  cout << " New Histo with the number of digi per evt per wire: " << hEvtPerWireH->GetName() << endl;
192  theHistoEvtPerWireMap[dtLId]=hEvtPerWireH;
193  }
194  }
195  }
196  }
197  }
198  }
199 
200  //Fill the plot of the number of digi per event per wire
201  std::map<int,int > DigiPerWirePerEvent;
202  // LOOP OVER ALL THE CHAMBERS
203  vector<DTChamber*>::const_iterator ch_it = dtGeom->chambers().begin();
204  vector<DTChamber*>::const_iterator ch_end = dtGeom->chambers().end();
205  for (; ch_it != ch_end; ++ch_it) {
206  DTChamberId ch = (*ch_it)->id();
207  vector<const DTSuperLayer*>::const_iterator sl_it = (*ch_it)->superLayers().begin();
208  vector<const DTSuperLayer*>::const_iterator sl_end = (*ch_it)->superLayers().end();
209  // Loop over the SLs
210  for(; sl_it != sl_end; ++sl_it) {
211  DTSuperLayerId sl = (*sl_it)->id();
212  vector<const DTLayer*>::const_iterator l_it = (*sl_it)->layers().begin();
213  vector<const DTLayer*>::const_iterator l_end = (*sl_it)->layers().end();
214  // Loop over the Ls
215  for(; l_it != l_end; ++l_it) {
216  DTLayerId layerId = (*l_it)->id();
217 
218  // Get the number of wires
219  const DTTopology& dtTopo = dtGeom->layer(layerId)->specificTopology();
220  const int firstWire = dtTopo.firstChannel();
221  const int lastWire = dtTopo.lastChannel();
222 
223  if (theHistoEvtPerWireMap.find(layerId) != theHistoEvtPerWireMap.end() &&
224  skippedPlot[layerId] == counter) {
225 
226  for (int wire=firstWire; wire<=lastWire; wire++) {
227  DigiPerWirePerEvent[wire]= 0;
228  }
229  // loop over all the digis of the event
230  DTDigiCollection::Range layerDigi= dtdigis->get(layerId);
231  for (DTDigiCollection::const_iterator digi = layerDigi.first;
232  digi!=layerDigi.second;
233  ++digi){
234  if((cosmicRun && (*digi).countsTDC()<upperLimit) || (!cosmicRun))
235  DigiPerWirePerEvent[(*digi).wire()]+=1;
236  }
237  // fill the digi event histo
238  for (int wire=firstWire; wire<=lastWire; wire++) {
239  theFile->cd();
240  int histoEvents = nevents - (counter*1000);
241  theHistoEvtPerWireMap[layerId]->Fill(histoEvents,wire,DigiPerWirePerEvent[wire]);
242  }
243  }
244  } //Loop Ls
245  } //Loop SLs
246  } //Loop chambers
247 
248 
249  if(nevents % 1000 == 0) {
250  counter++;
251  // save the digis event plot on file
252  for(map<DTLayerId, TH2F* >::const_iterator lHisto = theHistoEvtPerWireMap.begin();
253  lHisto != theHistoEvtPerWireMap.end();
254  lHisto++) {
255  theFile->cd();
256  if((*lHisto).second)
257  (*lHisto).second->Write();
258  }
259  theHistoEvtPerWireMap.clear();
260  }
261 
262 }
const int lastChannel() const
Returns the wire number of the last wire.
Definition: DTTopology.h:82
T getUntrackedParameter(std::string const &, T const &) const
const int channels() const
Returns the number of wires in the layer.
Definition: DTTopology.h:77
edm::ParameterSet parameters
DTChamberId chamberId() const
Return the corresponding ChamberId.
const int firstChannel() const
Returns the wire number of the first wire.
Definition: DTTopology.h:80
DTSuperLayerId superlayerId() const
Return the corresponding SuperLayerId.
Definition: DTLayerId.h:61
U second(std::pair< T, U > const &p)
edm::ESHandle< DTGeometry > dtGeom
std::string getLayerName(const DTLayerId &lId) const
Get the name of the layer.
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:359
float tTrig
tTrig from the DB
int TriggerWidth
variables to set by configuration file
std::vector< DigiType >::const_iterator const_iterator
std::string HistoName
int sector() const
Definition: DTChamberId.h:63
std::map< DTLayerId, TH2F * > theHistoEvtPerWireMap
std::pair< const_iterator, const_iterator > Range
tuple cout
Definition: gather_cfg.py:41
edm::ESHandle< DTTtrig > tTrigMap
std::map< DTLayerId, int > skippedPlot
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:47
std::map< DTLayerId, TH1F * > theHistoOccupancyMap
void DTNoiseCalibration::beginRun ( const edm::Run run,
const edm::EventSetup setup 
)
virtual

Reimplemented from edm::EDAnalyzer.

Definition at line 80 of file DTNoiseCalibration.cc.

References cmsDriverOptions::counter, gather_cfg::cout, edm::EventSetup::get(), nevents, and ExpressReco_HICollisions_FallBack::parameters.

80  {
81 
82  cout <<"[DTNoiseCalibration]: BeginJob"<<endl;
83  nevents = 0;
84  counter = 0;
85 
86  // Get the DT Geometry
87  setup.get<MuonGeometryRecord>().get(dtGeom);
88 
89  // tTrig
90  if (parameters.getUntrackedParameter<bool>("readDB", true))
91  setup.get<DTTtrigRcd>().get(dbLabel,tTrigMap);
92 
93  // TDC time distribution
94  int numBin = (TriggerWidth*(32/25))/50;
95  hTDCTriggerWidth = new TH1F("TDC_Time_Distribution", "TDC_Time_Distribution", numBin, 0, TriggerWidth*(32/25));
96 
97 }
T getUntrackedParameter(std::string const &, T const &) const
edm::ParameterSet parameters
edm::ESHandle< DTGeometry > dtGeom
int TriggerWidth
variables to set by configuration file
const T & get() const
Definition: EventSetup.h:55
tuple cout
Definition: gather_cfg.py:41
edm::ESHandle< DTTtrig > tTrigMap
void DTNoiseCalibration::endJob ( void  )
virtual

Endjob.

Reimplemented from edm::EDAnalyzer.

Definition at line 265 of file DTNoiseCalibration.cc.

References newFWLiteAna::bin, DTTimeUnits::counts, gather_cfg::cout, debug, DTTopology::firstChannel(), DTTopology::lastChannel(), nevents, ExpressReco_HICollisions_FallBack::parameters, record, DTStatusFlag::setCellNoise(), crabStatusFromReport::statusMap, and interactiveExample::theFile.

265  {
266 
267  cout << "[DTNoiseCalibration] endjob called!" <<endl;
268 
269  // save the TDC digi plot
270  theFile->cd();
271  hTDCTriggerWidth->Write();
272 
273  // save on file the occupancy histo and write the list of noisy cells
274  double TriggerWidth_s=0;
276  for(map<DTLayerId, TH1F*>::const_iterator lHisto = theHistoOccupancyMap.begin();
277  lHisto != theHistoOccupancyMap.end();
278  lHisto++) {
279  if(cosmicRun){
280  if ( parameters.getUntrackedParameter<bool>("readDB", true) )
281  tTrigMap->get( ((*lHisto).first).superlayerId(), tTrig, tTrigRMS, kFactor,
283  else tTrig = parameters.getUntrackedParameter<int>("defaultTtrig", 4000);
284  double TriggerWidth_ns = ((tTrig-theOffset)*25)/32;
285  TriggerWidth_s = TriggerWidth_ns/1e9;
286  }
287  if(!cosmicRun)
288  TriggerWidth_s = double(TriggerWidth/1e9);
289  if(debug)
290  cout<<"TriggerWidth (s): "<<TriggerWidth_s<<" TotEvents: "<<nevents<<endl;
291  double normalization = 1/double(nevents*TriggerWidth_s);
292  if((*lHisto).second){
293  (*lHisto).second->Scale(normalization);
294  theFile->cd();
295  (*lHisto).second->Write();
296  const DTTopology& dtTopo = dtGeom->layer((*lHisto).first)->specificTopology();
297  const int firstWire = dtTopo.firstChannel();
298  const int lastWire = dtTopo.lastChannel();
299  for(int bin=firstWire; bin<=lastWire; bin++){
300  //from definition of "noisy cell"
301  if((*lHisto).second->GetBinContent(bin)>500){
302  DTWireId wireID((*lHisto).first, bin);
303  statusMap->setCellNoise(wireID,1);
304  }
305  }
306  }
307  }
308  cout << "Writing Noise Map object to DB!" << endl;
309  string record = "DTStatusFlagRcd";
310  DTCalibDBUtils::writeToDB<DTStatusFlag>(record, statusMap);
311 
312 
313  /*
314  //save the digi event plot per SuperLayer
315  bool histo=false;
316  map<DTSuperLayerId, vector<int> > maxPerSuperLayer;
317  int numPlot = (nevents/1000);
318  int num=0;
319  // loop over the numPlot
320  for(int i=0; i<numPlot; i++){
321  vector<DTChamber*>::const_iterator ch_it = dtGeom->chambers().begin();
322  vector<DTChamber*>::const_iterator ch_end = dtGeom->chambers().end();
323  // Loop over the chambers
324  for (; ch_it != ch_end; ++ch_it) {
325  DTChamberId ch = (*ch_it)->id();
326  vector<const DTSuperLayer*>::const_iterator sl_it = (*ch_it)->superLayers().begin();
327  vector<const DTSuperLayer*>::const_iterator sl_end = (*ch_it)->superLayers().end();
328  // Loop over the SLs
329  for(; sl_it != sl_end; ++sl_it) {
330  DTSuperLayerId sl = (*sl_it)->id();
331  vector<const DTLayer*>::const_iterator l_it = (*sl_it)->layers().begin();
332  vector<const DTLayer*>::const_iterator l_end = (*sl_it)->layers().end();
333  double dummy = pow(10.,10.);
334  maxPerSuperLayer[sl].push_back(0);
335  // Loop over the Ls
336  for(; l_it != l_end; ++l_it) {
337  DTLayerId layerId = (*l_it)->id();
338 
339  if (theHistoEvtPerWireMap.find(layerId) != theHistoEvtPerWireMap.end() &&
340  theHistoEvtPerWireMap[layerId].size() > i){
341  if (theHistoEvtPerWireMap[layerId][i]->GetMaximum(dummy)>maxPerSuperLayer[sl][i])
342  maxPerSuperLayer[sl][i] = theHistoEvtPerWireMap[layerId][i]->GetMaximum(dummy);
343  }
344  }
345  } // loop over SLs
346  } // loop over chambers
347  } // loop over numPlot
348 
349  // loop over the numPlot
350  for(int i=0; i<numPlot; i++){
351  vector<DTChamber*>::const_iterator chamber_it = dtGeom->chambers().begin();
352  vector<DTChamber*>::const_iterator chamber_end = dtGeom->chambers().end();
353  // Loop over the chambers
354  for (; chamber_it != chamber_end; ++chamber_it) {
355  DTChamberId ch = (*chamber_it)->id();
356  vector<const DTSuperLayer*>::const_iterator sl_it = (*chamber_it)->superLayers().begin();
357  vector<const DTSuperLayer*>::const_iterator sl_end = (*chamber_it)->superLayers().end();
358  // Loop over the SLs
359  for(; sl_it != sl_end; ++sl_it) {
360  DTSuperLayerId sl = (*sl_it)->id();
361  vector<const DTLayer*>::const_iterator l_it = (*sl_it)->layers().begin();
362  vector<const DTLayer*>::const_iterator l_end = (*sl_it)->layers().end();
363 
364  stringstream num; num << i;
365  string canvasName = "c" + getSuperLayerName(sl) + "_" + num.str();
366  TCanvas c1(canvasName.c_str(),canvasName.c_str(),600,780);
367  TLegend *leg=new TLegend(0.5,0.6,0.7,0.8);
368  for(; l_it != l_end; ++l_it) {
369  DTLayerId layerId = (*l_it)->id();
370 
371  if (theHistoEvtPerWireMap.find(layerId) != theHistoEvtPerWireMap.end() &&
372  theHistoEvtPerWireMap[layerId].size() > i){
373 
374  string TitleName = "DigiPerWirePerEvent_" + getSuperLayerName(sl) + "_" + num.str();
375  theHistoEvtPerWireMap[layerId][i]->SetTitle(TitleName.c_str());
376  stringstream layer; layer << layerId.layer();
377  string legendHisto = "layer " + layer.str();
378  leg->AddEntry(theHistoEvtPerWireMap[layerId][i],legendHisto.c_str(),"L");
379  theHistoEvtPerWireMap[layerId][i]->SetMaximum(maxPerSuperLayer[sl][i]);
380  if(histo==false)
381  theHistoEvtPerWireMap[layerId][i]->Draw("lego");
382  else
383  theHistoEvtPerWireMap[layerId][i]->Draw("same , lego");
384  theHistoEvtPerWireMap[layerId][i]->SetLineColor(layerId.layer());
385  histo=true;
386  }
387  } // loop over Ls
388  if(histo){
389  leg->Draw("same");
390  theFile->cd();
391  c1.Write();
392  }
393  histo=false;
394  } // loop over SLs
395  } // loop over chambers
396  } // loop over numPlot
397  */
398 
399 }
const int lastChannel() const
Returns the wire number of the last wire.
Definition: DTTopology.h:82
T getUntrackedParameter(std::string const &, T const &) const
JetCorrectorParameters::Record record
Definition: classes.h:11
edm::ParameterSet parameters
const int firstChannel() const
Returns the wire number of the first wire.
Definition: DTTopology.h:80
int setCellNoise(int wheelId, int stationId, int sectorId, int slId, int layerId, int cellId, bool flag)
edm::ESHandle< DTGeometry > dtGeom
float tTrig
tTrig from the DB
int TriggerWidth
variables to set by configuration file
tuple cout
Definition: gather_cfg.py:41
edm::ESHandle< DTTtrig > tTrigMap
std::map< DTLayerId, TH1F * > theHistoOccupancyMap
string DTNoiseCalibration::getLayerName ( const DTLayerId lId) const
private

Get the name of the layer.

Definition at line 412 of file DTNoiseCalibration.cc.

References DTSuperLayerId::chamberId(), DTLayerId::layer(), DTChamberId::sector(), DTChamberId::station(), relativeConstraints::station, DTSuperLayerId::superlayer(), DTLayerId::superlayerId(), and DTChamberId::wheel().

412  {
413 
414  const DTSuperLayerId dtSLId = lId.superlayerId();
415  const DTChamberId dtChId = dtSLId.chamberId();
416  stringstream Layer; Layer << lId.layer();
417  stringstream superLayer; superLayer << dtSLId.superlayer();
418  stringstream wheel; wheel << dtChId.wheel();
419  stringstream station; station << dtChId.station();
420  stringstream sector; sector << dtChId.sector();
421 
422  string LayerName =
423  "W" + wheel.str()
424  + "_St" + station.str()
425  + "_Sec" + sector.str()
426  + "_SL" + superLayer.str()
427  + "_L" + Layer.str();
428 
429  return LayerName;
430 
431 }
DTChamberId chamberId() const
Return the corresponding ChamberId.
int layer() const
Return the layer number.
Definition: DTLayerId.h:55
DTSuperLayerId superlayerId() const
Return the corresponding SuperLayerId.
Definition: DTLayerId.h:61
int superlayer() const
Return the superlayer number (deprecated method name)
int sector() const
Definition: DTChamberId.h:63
int station() const
Return the station number.
Definition: DTChamberId.h:53
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:47
string DTNoiseCalibration::getSuperLayerName ( const DTSuperLayerId dtSLId) const
private

Get the name of the superLayer.

Definition at line 434 of file DTNoiseCalibration.cc.

References DTSuperLayerId::chamberId(), DTChamberId::sector(), DTChamberId::station(), relativeConstraints::station, DTSuperLayerId::superlayer(), and DTChamberId::wheel().

434  {
435 
436  const DTChamberId dtChId = dtSLId.chamberId();
437  stringstream superLayer; superLayer << dtSLId.superlayer();
438  stringstream wheel; wheel << dtChId.wheel();
439  stringstream station; station << dtChId.station();
440  stringstream sector; sector << dtChId.sector();
441 
442  string SuperLayerName =
443  "W" + wheel.str()
444  + "_St" + station.str()
445  + "_Sec" + sector.str()
446  + "_SL" + superLayer.str();
447 
448  return SuperLayerName;
449 
450 }
DTChamberId chamberId() const
Return the corresponding ChamberId.
int superlayer() const
Return the superlayer number (deprecated method name)
int sector() const
Definition: DTChamberId.h:63
int station() const
Return the station number.
Definition: DTChamberId.h:53
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:47

Member Data Documentation

bool DTNoiseCalibration::cosmicRun
private

Definition at line 63 of file DTNoiseCalibration.h.

int DTNoiseCalibration::counter
private

Definition at line 58 of file DTNoiseCalibration.h.

std::string DTNoiseCalibration::dbLabel
private

Definition at line 75 of file DTNoiseCalibration.h.

bool DTNoiseCalibration::debug
private

Definition at line 56 of file DTNoiseCalibration.h.

std::string DTNoiseCalibration::digiLabel
private

Definition at line 78 of file DTNoiseCalibration.h.

edm::ESHandle<DTGeometry> DTNoiseCalibration::dtGeom
private

Definition at line 86 of file DTNoiseCalibration.h.

bool DTNoiseCalibration::fastAnalysis
private

Definition at line 64 of file DTNoiseCalibration.h.

TH1F* DTNoiseCalibration::hTDCTriggerWidth
private

Definition at line 83 of file DTNoiseCalibration.h.

float DTNoiseCalibration::kFactor
private

Definition at line 71 of file DTNoiseCalibration.h.

int DTNoiseCalibration::nevents
private

Definition at line 57 of file DTNoiseCalibration.h.

edm::ParameterSet DTNoiseCalibration::parameters
private
int DTNoiseCalibration::sect
private

Definition at line 66 of file DTNoiseCalibration.h.

std::map<DTLayerId, int> DTNoiseCalibration::skippedPlot
private

Definition at line 101 of file DTNoiseCalibration.h.

TFile* DTNoiseCalibration::theFile
private

Definition at line 92 of file DTNoiseCalibration.h.

std::map<DTLayerId, TH2F*> DTNoiseCalibration::theHistoEvtPerWireMap
private

Definition at line 95 of file DTNoiseCalibration.h.

std::map<DTLayerId, TH1F*> DTNoiseCalibration::theHistoOccupancyMap
private

Definition at line 98 of file DTNoiseCalibration.h.

double DTNoiseCalibration::theOffset
private

Definition at line 73 of file DTNoiseCalibration.h.

int DTNoiseCalibration::TriggerWidth
private

variables to set by configuration file

Definition at line 61 of file DTNoiseCalibration.h.

float DTNoiseCalibration::tTrig
private

tTrig from the DB

Definition at line 69 of file DTNoiseCalibration.h.

edm::ESHandle<DTTtrig> DTNoiseCalibration::tTrigMap
private

Definition at line 89 of file DTNoiseCalibration.h.

float DTNoiseCalibration::tTrigRMS
private

Definition at line 70 of file DTNoiseCalibration.h.

float DTNoiseCalibration::upperLimit
private

Definition at line 62 of file DTNoiseCalibration.h.

int DTNoiseCalibration::wh
private

Definition at line 65 of file DTNoiseCalibration.h.