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
SiStripBaselineAnalyzer Class Reference

#include <Validation/SiStripAnalyzer/src/SiStripBaselineAnalyzer.cc>

Inheritance diagram for SiStripBaselineAnalyzer:
edm::EDAnalyzer

Public Member Functions

 SiStripBaselineAnalyzer (const edm::ParameterSet &)
 
 ~SiStripBaselineAnalyzer ()
 
- Public Member Functions inherited from edm::EDAnalyzer
 EDAnalyzer ()
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 

Private Member Functions

virtual void analyze (const edm::Event &, const edm::EventSetup &)
 
virtual void beginJob ()
 
virtual void endJob ()
 

Private Attributes

uint16_t actualModule_
 
TCanvas * Canvas_
 
edm::Service< TFileServicefs_
 
TH1F * h1APVCM_
 
TH1F * h1BadAPVperEvent_
 
TH1F * h1Baseline_
 
TH1F * h1Clusters_
 
TH1F * h1Pedestals_
 
TH1F * h1ProcessedRawDigis_
 
uint16_t nModuletoDisplay_
 
std::vector< int > pedestals
 
edm::ESHandle< SiStripPedestalspedestalsHandle
 
uint32_t peds_cache_id
 
bool plotAPVCM_
 
bool plotBaseline_
 
bool plotBaselinePoints_
 
bool plotClusters_
 
bool plotPedestals_
 
bool plotRawDigi_
 
edm::InputTag srcAPVCM_
 
edm::InputTag srcBaseline_
 
edm::InputTag srcBaselinePoints_
 
edm::InputTag srcProcessedRawDigi_
 
std::auto_ptr
< SiStripPedestalsSubtractor
subtractorPed_
 
std::vector< TH1F > vBaselineHisto_
 
std::vector< TH1F > vBaselinePointsHisto_
 
std::vector< TH1F > vClusterHisto_
 
std::vector< TH1F > vProcessedRawDigiHisto_
 

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)
 
static void prevalidate (ConfigurationDescriptions &)
 
- Protected Member Functions inherited from edm::EDAnalyzer
CurrentProcessingContext const * currentContext () const
 

Detailed Description

Description: <one line="" class="" summary>="">

Implementation: <Notes on="" implementation>="">

Definition at line 75 of file SiStripBaselineAnalyzer.cc.

Constructor & Destructor Documentation

SiStripBaselineAnalyzer::SiStripBaselineAnalyzer ( const edm::ParameterSet conf)
explicit

Definition at line 124 of file SiStripBaselineAnalyzer.cc.

References SiStripRawProcessingFactory::create_SubtractorPed(), fs_, edm::ParameterSet::getParameter(), h1APVCM_, h1BadAPVperEvent_, h1Pedestals_, TFileDirectory::make(), nModuletoDisplay_, plotAPVCM_, plotBaseline_, plotBaselinePoints_, plotClusters_, plotPedestals_, plotRawDigi_, srcAPVCM_, srcBaseline_, srcBaselinePoints_, srcProcessedRawDigi_, and subtractorPed_.

124  {
125 
126  srcBaseline_ = conf.getParameter<edm::InputTag>( "srcBaseline" );
127  srcBaselinePoints_ = conf.getParameter<edm::InputTag>( "srcBaselinePoints" );
128  srcProcessedRawDigi_ = conf.getParameter<edm::InputTag>( "srcProcessedRawDigi" );
129  srcAPVCM_ = conf.getParameter<edm::InputTag>( "srcAPVCM" );
131  nModuletoDisplay_ = conf.getParameter<uint32_t>( "nModuletoDisplay" );
132  plotClusters_ = conf.getParameter<bool>( "plotClusters" );
133  plotBaseline_ = conf.getParameter<bool>( "plotBaseline" );
134  plotBaselinePoints_ = conf.getParameter<bool>( "plotBaselinePoints" );
135  plotRawDigi_ = conf.getParameter<bool>( "plotRawDigi" );
136  plotAPVCM_ = conf.getParameter<bool>( "plotAPVCM" );
137  plotPedestals_ = conf.getParameter<bool>( "plotPedestals" );
138 
139  h1BadAPVperEvent_ = fs_->make<TH1F>("BadAPV/Event","BadAPV/Event", 2001, -0.5, 2000.5);
140  h1BadAPVperEvent_->SetXTitle("# Modules with Bad APVs");
141  h1BadAPVperEvent_->SetYTitle("Entries");
142  h1BadAPVperEvent_->SetLineWidth(2);
143  h1BadAPVperEvent_->SetLineStyle(2);
144 
145  h1APVCM_ = fs_->make<TH1F>("APV CM","APV CM", 2048, -1023.5, 1023.5);
146  h1APVCM_->SetXTitle("APV CM [adc]");
147  h1APVCM_->SetYTitle("Entries");
148  h1APVCM_->SetLineWidth(2);
149  h1APVCM_->SetLineStyle(2);
150 
151  h1Pedestals_ = fs_->make<TH1F>("Pedestals","Pedestals", 2048, -1023.5, 1023.5);
152  h1Pedestals_->SetXTitle("Pedestals [adc]");
153  h1Pedestals_->SetYTitle("Entries");
154  h1Pedestals_->SetLineWidth(2);
155  h1Pedestals_->SetLineStyle(2);
156 
157 
158 
159 }
T getParameter(std::string const &) const
edm::Service< TFileService > fs_
T * make() const
make new ROOT object
std::auto_ptr< SiStripPedestalsSubtractor > subtractorPed_
static std::auto_ptr< SiStripPedestalsSubtractor > create_SubtractorPed(const edm::ParameterSet &)
SiStripBaselineAnalyzer::~SiStripBaselineAnalyzer ( )

Definition at line 162 of file SiStripBaselineAnalyzer.cc.

163 {
164 
165 
166 
167 }

Member Function Documentation

void SiStripBaselineAnalyzer::analyze ( const edm::Event e,
const edm::EventSetup es 
)
privatevirtual

Implements edm::EDAnalyzer.

Definition at line 170 of file SiStripBaselineAnalyzer.cc.

References actualModule_, ecalMGPA::adc(), sistrip::APV, edm::DetSet< T >::begin(), edm::DetSetVector< T >::begin(), edmNew::DetSetVector< T >::begin(), combineCards::bins, gather_cfg::cout, edm::DetSet< T >::end(), edm::DetSetVector< T >::end(), edmNew::DetSetVector< T >::end(), edm::EventID::event(), event(), fs_, edm::EventSetup::get(), edm::Event::getByLabel(), h1APVCM_, h1BadAPVperEvent_, h1Baseline_, h1Clusters_, h1Pedestals_, h1ProcessedRawDigis_, i, edm::EventBase::id(), edmNew::DetSetVector< T >::id(), TFileDirectory::mkdir(), nModuletoDisplay_, pedestals, pedestalsHandle, peds_cache_id, plotAPVCM_, plotBaseline_, plotBaselinePoints_, plotClusters_, plotPedestals_, plotRawDigi_, edm::EventID::run(), DTTTrigCorrFirst::run, gather_cfg::runs, srcAPVCM_, srcBaseline_, srcProcessedRawDigi_, strip(), and subtractorPed_.

171 {
172  using namespace edm;
173  if(plotPedestals_&&actualModule_ ==0){
174  uint32_t p_cache_id = es.get<SiStripPedestalsRcd>().cacheIdentifier();
175  if(p_cache_id != peds_cache_id) {
177  peds_cache_id = p_cache_id;
178  }
179 
180 
181  std::vector<uint32_t> detIdV;
182  pedestalsHandle->getDetIds(detIdV);
183 
184  for(uint32_t i=0; i < detIdV.size(); ++i){
185  pedestals.clear();
186  SiStripPedestals::Range pedestalsRange = pedestalsHandle->getRange(detIdV[i]);
187  pedestals.resize((pedestalsRange.second- pedestalsRange.first)*8/10);
188  pedestalsHandle->allPeds(pedestals, pedestalsRange);
189  for(uint32_t it=0; it < pedestals.size(); ++it) h1Pedestals_->Fill(pedestals[it]);
190  }
191  }
192 
193  if(plotAPVCM_){
195  edm::InputTag CMLabel("siStripZeroSuppression:APVCM");
196  e.getByLabel(srcAPVCM_,moduleCM);
197 
198  edm::DetSetVector<SiStripProcessedRawDigi>::const_iterator itCMDetSetV =moduleCM->begin();
199  for (; itCMDetSetV != moduleCM->end(); ++itCMDetSetV){
201  for(;itCM != itCMDetSetV->end(); ++itCM) h1APVCM_->Fill(itCM->adc());
202  }
203  }
204 
205  if(!plotRawDigi_) return;
206  subtractorPed_->init(es);
207 
208 
209 
211  e.getByLabel(srcProcessedRawDigi_,moduleRawDigi);
212 
214  if(plotBaseline_) e.getByLabel(srcBaseline_, moduleBaseline);
215 
216  edm::Handle<edm::DetSetVector<SiStripDigi> > moduleBaselinePoints;
217  if(plotBaselinePoints_) e.getByLabel(srcBaseline_, moduleBaselinePoints);
218 
220  if(plotClusters_){
221  edm::InputTag clusLabel("siStripClusters");
222  e.getByLabel(clusLabel, clusters);
223  }
224 
225  char detIds[20];
226  char evs[20];
227  char runs[20];
228 
229 
230  TFileDirectory sdProcessedRawDigis_= fs_->mkdir("ProcessedRawDigis");
231  TFileDirectory sdBaseline_= fs_->mkdir("Baseline");
232  TFileDirectory sdBaselinePoints_= fs_->mkdir("BaselinePoints");
233  TFileDirectory sdClusters_= fs_->mkdir("Clusters");
234 
235 
237  if(plotBaseline_) itDSBaseline = moduleBaseline->begin();
238  edm::DetSetVector<SiStripRawDigi>::const_iterator itRawDigis = moduleRawDigi->begin();
239 
240  uint32_t NBabAPVs = moduleRawDigi->size();
241  std::cout<< "Number of module with HIP in this event: " << NBabAPVs << std::endl;
242  h1BadAPVperEvent_->Fill(NBabAPVs);
243 
244  for (; itRawDigis != moduleRawDigi->end(); ++itRawDigis) {
245  if(actualModule_ > nModuletoDisplay_) return;
246  uint32_t detId = itRawDigis->id;
247 
248  if(plotBaseline_){
249  //std::cout << "bas id: " << itDSBaseline->id << " raw id: " << detId << std::endl;
250  if(itDSBaseline->id != detId){
251  std::cout << "Collections out of Synch. Something of fishy is going on ;-)" << std::endl;
252  return;
253  }
254  }
255 
256 
257  actualModule_++;
258  uint32_t event = e.id().event();
259  uint32_t run = e.id().run();
260  //std::cout << "processing module N: " << actualModule_<< " detId: " << detId << " event: "<< event << std::endl;
261 
262 
263 
264  edm::DetSet<SiStripRawDigi>::const_iterator itRaw = itRawDigis->begin();
265  bool restAPV[6] = {0,0,0,0,0,0};
266  int strip =0, totADC=0;
267  int minAPVRes = 7, maxAPVRes = -1;
268  for(;itRaw != itRawDigis->end(); ++itRaw, ++strip){
269  float adc = itRaw->adc();
270  totADC+= adc;
271  if(strip%127 ==0){
272  //std::cout << "totADC " << totADC << std::endl;
273  int APV = strip/128;
274  if(totADC!= 0){
275  restAPV[APV] = true;
276  totADC =0;
277  if(APV>maxAPVRes) maxAPVRes = APV;
278  if(APV<minAPVRes) minAPVRes = APV;
279  }
280  }
281  }
282 
283  uint16_t bins =768;
284  float minx = -0.5, maxx=767.5;
285  if(minAPVRes !=7){
286  minx = minAPVRes * 128 -0.5;
287  maxx = maxAPVRes * 128 + 127.5;
288  bins = maxx-minx;
289  }
290 
291  sprintf(detIds,"%ul", detId);
292  sprintf(evs,"%ul", event);
293  sprintf(runs,"%ul", run);
294  char* dHistoName = Form("Id:%s_run:%s_ev:%s",detIds, runs, evs);
295  h1ProcessedRawDigis_ = sdProcessedRawDigis_.make<TH1F>(dHistoName,dHistoName, bins, minx, maxx);
296 
297  if(plotBaseline_){
298  h1Baseline_ = sdBaseline_.make<TH1F>(dHistoName,dHistoName, bins, minx, maxx);
299  h1Baseline_->SetXTitle("strip#");
300  h1Baseline_->SetYTitle("ADC");
301  h1Baseline_->SetMaximum(1024.);
302  h1Baseline_->SetMinimum(-300.);
303  h1Baseline_->SetLineWidth(2);
304  h1Baseline_->SetLineStyle(2);
305  h1Baseline_->SetLineColor(2);
306  }
307 
308  if(plotClusters_){
309  h1Clusters_ = sdClusters_.make<TH1F>(dHistoName,dHistoName, bins, minx, maxx);
310 
311  h1Clusters_->SetXTitle("strip#");
312  h1Clusters_->SetYTitle("ADC");
313  h1Clusters_->SetMaximum(1024.);
314  h1Clusters_->SetMinimum(-300.);
315  h1Clusters_->SetLineWidth(2);
316  h1Clusters_->SetLineStyle(2);
317  h1Clusters_->SetLineColor(3);
318  }
319 
320  h1ProcessedRawDigis_->SetXTitle("strip#");
321  h1ProcessedRawDigis_->SetYTitle("ADC");
322  h1ProcessedRawDigis_->SetMaximum(1024.);
323  h1ProcessedRawDigis_->SetMinimum(-300.);
324  h1ProcessedRawDigis_->SetLineWidth(2);
325 
326  std::vector<int16_t> ProcessedRawDigis(itRawDigis->size());
327  subtractorPed_->subtract( *itRawDigis, ProcessedRawDigis);
328 
330  if(plotBaseline_) itBaseline = itDSBaseline->begin();
331  std::vector<int16_t>::const_iterator itProcessedRawDigis;
332 
333  strip =0;
334  for(itProcessedRawDigis = ProcessedRawDigis.begin();itProcessedRawDigis != ProcessedRawDigis.end(); ++itProcessedRawDigis){
335  if(restAPV[strip/128]){
336  float adc = *itProcessedRawDigis;
337  h1ProcessedRawDigis_->Fill(strip, adc);
338  if(plotBaseline_){
339  h1Baseline_->Fill(strip, itBaseline->adc());
340  ++itBaseline;
341  }
342  }
343  ++strip;
344  }
345 
346  if(plotBaseline_) ++itDSBaseline;
347  if(plotClusters_){
348  edmNew::DetSetVector<SiStripCluster>::const_iterator itClusters = clusters->begin();
349  for ( ; itClusters != clusters->end(); ++itClusters ){
350  for ( edmNew::DetSet<SiStripCluster>::const_iterator clus = itClusters->begin(); clus != itClusters->end(); ++clus){
351  if(itClusters->id() == detId){
352  int firststrip = clus->firstStrip();
353  //std::cout << "Found cluster in detId " << detId << " " << firststrip << " " << clus->amplitudes().size() << " -----------------------------------------------" << std::endl;
354  strip=0;
355  for( std::vector<uint8_t>::const_iterator itAmpl = clus->amplitudes().begin(); itAmpl != clus->amplitudes().end(); ++itAmpl){
356  h1Clusters_->Fill(firststrip+strip, *itAmpl);
357  ++strip;
358  }
359  }
360  }
361  }
362  }
363  }
364 
365 }
int adc(sample_type sample)
get the ADC sample (12 bits)
RunNumber_t run() const
Definition: EventID.h:42
iterator end()
Definition: DetSet.h:61
EventNumber_t event() const
Definition: EventID.h:44
int i
Definition: DBlmapReader.cc:9
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
const_iterator begin() const
void strip(std::string &input, const std::string &blanks=" \n\t")
Definition: stringTools.cc:16
data_type const * const_iterator
Definition: DetSetNew.h:25
std::pair< ContainerIterator, ContainerIterator > Range
id_type id(size_t cell) const
tuple runs
Definition: gather_cfg.py:87
const_iterator end() const
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
edm::Service< TFileService > fs_
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
iterator end()
Return the off-the-end iterator.
Definition: DetSetVector.h:356
iterator begin()
Definition: DetSet.h:60
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
const T & get() const
Definition: EventSetup.h:55
edm::EventID id() const
Definition: EventBase.h:56
tuple cout
Definition: gather_cfg.py:121
iterator begin()
Return an iterator to the first DetSet.
Definition: DetSetVector.h:341
edm::ESHandle< SiStripPedestals > pedestalsHandle
collection_type::const_iterator const_iterator
Definition: DetSet.h:34
collection_type::const_iterator const_iterator
Definition: DetSetVector.h:106
std::auto_ptr< SiStripPedestalsSubtractor > subtractorPed_
void SiStripBaselineAnalyzer::beginJob ( void  )
privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 369 of file SiStripBaselineAnalyzer.cc.

References actualModule_.

370 {
371 
372 
373 actualModule_ =0;
374 
375 
376 }
void SiStripBaselineAnalyzer::endJob ( void  )
privatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 380 of file SiStripBaselineAnalyzer.cc.

380  {
381 
382 }

Member Data Documentation

uint16_t SiStripBaselineAnalyzer::actualModule_
private

Definition at line 120 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and beginJob().

TCanvas* SiStripBaselineAnalyzer::Canvas_
private

Definition at line 113 of file SiStripBaselineAnalyzer.cc.

edm::Service<TFileService> SiStripBaselineAnalyzer::fs_
private

Definition at line 103 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

TH1F* SiStripBaselineAnalyzer::h1APVCM_
private

Definition at line 110 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

TH1F* SiStripBaselineAnalyzer::h1BadAPVperEvent_
private

Definition at line 105 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

TH1F* SiStripBaselineAnalyzer::h1Baseline_
private

Definition at line 108 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze().

TH1F* SiStripBaselineAnalyzer::h1Clusters_
private

Definition at line 109 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze().

TH1F* SiStripBaselineAnalyzer::h1Pedestals_
private

Definition at line 111 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

TH1F* SiStripBaselineAnalyzer::h1ProcessedRawDigis_
private

Definition at line 107 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze().

uint16_t SiStripBaselineAnalyzer::nModuletoDisplay_
private

Definition at line 119 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

std::vector<int> SiStripBaselineAnalyzer::pedestals
private

Definition at line 88 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze().

edm::ESHandle<SiStripPedestals> SiStripBaselineAnalyzer::pedestalsHandle
private

Definition at line 87 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze().

uint32_t SiStripBaselineAnalyzer::peds_cache_id
private

Definition at line 89 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze().

bool SiStripBaselineAnalyzer::plotAPVCM_
private

Definition at line 95 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

bool SiStripBaselineAnalyzer::plotBaseline_
private

Definition at line 92 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

bool SiStripBaselineAnalyzer::plotBaselinePoints_
private

Definition at line 93 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

bool SiStripBaselineAnalyzer::plotClusters_
private

Definition at line 91 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

bool SiStripBaselineAnalyzer::plotPedestals_
private

Definition at line 96 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

bool SiStripBaselineAnalyzer::plotRawDigi_
private

Definition at line 94 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

edm::InputTag SiStripBaselineAnalyzer::srcAPVCM_
private

Definition at line 100 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

edm::InputTag SiStripBaselineAnalyzer::srcBaseline_
private

Definition at line 98 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

edm::InputTag SiStripBaselineAnalyzer::srcBaselinePoints_
private

Definition at line 99 of file SiStripBaselineAnalyzer.cc.

Referenced by SiStripBaselineAnalyzer().

edm::InputTag SiStripBaselineAnalyzer::srcProcessedRawDigi_
private

Definition at line 101 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

std::auto_ptr<SiStripPedestalsSubtractor> SiStripBaselineAnalyzer::subtractorPed_
private

Definition at line 86 of file SiStripBaselineAnalyzer.cc.

Referenced by analyze(), and SiStripBaselineAnalyzer().

std::vector<TH1F> SiStripBaselineAnalyzer::vBaselineHisto_
private

Definition at line 115 of file SiStripBaselineAnalyzer.cc.

std::vector<TH1F> SiStripBaselineAnalyzer::vBaselinePointsHisto_
private

Definition at line 116 of file SiStripBaselineAnalyzer.cc.

std::vector<TH1F> SiStripBaselineAnalyzer::vClusterHisto_
private

Definition at line 117 of file SiStripBaselineAnalyzer.cc.

std::vector<TH1F> SiStripBaselineAnalyzer::vProcessedRawDigiHisto_
private

Definition at line 114 of file SiStripBaselineAnalyzer.cc.