CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
EcalMatacqAnalyzer.cc
Go to the documentation of this file.
1 /*
2  * \class EcalMatacqAnalyzer
3  *
4  * primary author: Gautier Hamel De Monchenault - CEA/Saclay
5  * author: Julie Malcles - CEA/Saclay
6  */
7 
8 
9 #include <TFile.h>
10 #include <TTree.h>
11 #include <TProfile.h>
12 #include <TChain.h>
13 #include <vector>
14 
15 #include <CalibCalorimetry/EcalLaserAnalyzer/plugins/EcalMatacqAnalyzer.h>
16 
17 #include <sstream>
18 #include <iostream>
19 #include <iomanip>
20 
21 #include <FWCore/MessageLogger/interface/MessageLogger.h>
22 #include <FWCore/Utilities/interface/Exception.h>
23 
24 #include <DataFormats/EcalDetId/interface/EcalDetIdCollections.h>
25 #include <DataFormats/EcalRawData/interface/EcalRawDataCollections.h>
26 
27 #include <FWCore/Framework/interface/Event.h>
28 #include <FWCore/Framework/interface/MakerMacros.h>
29 #include <FWCore/ParameterSet/interface/ParameterSet.h>
30 #include <DataFormats/EcalDigi/interface/EcalDigiCollections.h>
31 
32 #include <CalibCalorimetry/EcalLaserAnalyzer/interface/TMatacq.h>
33 #include <CalibCalorimetry/EcalLaserAnalyzer/interface/TMom.h>
34 #include <CalibCalorimetry/EcalLaserAnalyzer/interface/TMTQ.h>
35 
36 using namespace std;
37 
38 //========================================================================
39 EcalMatacqAnalyzer::EcalMatacqAnalyzer(const edm::ParameterSet& iConfig)
40  //========================================================================
41  :
42 
43 iEvent(0),
44 
45 // framework parameters with default values
46 
47 _presample( iConfig.getUntrackedParameter< double >( "nPresamples", 6.7 ) ),
48 _nsamplesaftmax(iConfig.getUntrackedParameter< unsigned int >( "nSamplesAftMax", 80 ) ),
49 _noiseCut( iConfig.getUntrackedParameter< unsigned int >( "noiseCut", 7 ) ),
50 _parabnbefmax( iConfig.getUntrackedParameter< unsigned int >( "paraBeforeMax", 8 ) ),
51 _parabnaftmax( iConfig.getUntrackedParameter< unsigned int >( "paraAfterMax", 7 ) ),
52 _thres( iConfig.getUntrackedParameter< unsigned int >( "threshold", 10 ) ),
53 _lowlev( iConfig.getUntrackedParameter< unsigned int >( "lowLevel", 20 ) ),
54 _highlev( iConfig.getUntrackedParameter< unsigned int >( "highLevel", 80 ) ),
55 _nevlasers( iConfig.getUntrackedParameter< unsigned int >( "nEventLaser", 600 ) ),
56 _timebefmax( iConfig.getUntrackedParameter< unsigned int >( "timeBefMax", 100 ) ),
57 _timeaftmax( iConfig.getUntrackedParameter< unsigned int >( "timeAftMax", 250 ) ),
58 _cutwindow( iConfig.getUntrackedParameter< double >( "cutWindow", 0.1 ) ),
59 _nsamplesshape( iConfig.getUntrackedParameter< unsigned int >( "nSamplesShape", 250 ) ),
60 _presampleshape(iConfig.getUntrackedParameter< unsigned int >( "nPresamplesShape",50 ) ),
61 _slide( iConfig.getUntrackedParameter< unsigned int >( "nSlide", 100 ) ),
62 _fedid( iConfig.getUntrackedParameter< int >( "fedID", -999 ) ),
63 _debug( iConfig.getUntrackedParameter< int >( "debug", 0 ) ),
64 nSides(NSIDES), lightside(0), runType(-1), runNum(0),
65 event(0), color(-1), maxsamp(0), nsamples(0), tt(0)
66 
67  //========================================================================
68 {
69 
70 
71  //now do what ever initialization is needed
72 
73  resdir_ = iConfig.getUntrackedParameter<std::string>("resDir");
74 
75  digiCollection_ = iConfig.getParameter<std::string>("digiCollection");
76  digiProducer_ = iConfig.getParameter<std::string>("digiProducer");
77 
78  eventHeaderCollection_ = iConfig.getParameter<std::string>("eventHeaderCollection");
79  eventHeaderProducer_ = iConfig.getParameter<std::string>("eventHeaderProducer");
80 
81 }
82 
83 //========================================================================
84 EcalMatacqAnalyzer::~EcalMatacqAnalyzer(){
85 //========================================================================
86 
87 // do anything here that needs to be done at desctruction time
88 // (e.g. close files, deallocate resources etc.)
89 
90 }
91 
92 
93 
94 //========================================================================
95 void EcalMatacqAnalyzer::beginJob() {
96 //========================================================================
97 
98 // Define temporary file name
99 
100  sampfile=resdir_;
101  sampfile+="/TmpTreeMatacqAnalyzer.root";
102 
103  sampFile = new TFile(sampfile.c_str(),"RECREATE");
104 
105 
106  // declaration of the tree to fill
107 
108  tree = new TTree("MatacqTree","MatacqTree");
109 
110 
111  //List of branches
112 
113  tree->Branch( "event", &event, "event/I" );
114  tree->Branch( "color", &color , "color/I" );
115  tree->Branch( "matacq", &matacq , "matacq[2560]/D" );
116  tree->Branch( "nsamples", &nsamples , "nsamples/I" );
117  tree->Branch( "maxsamp", &maxsamp , "maxsamp/I" );
118  tree->Branch( "tt", &tt , "tt/D" );
119  tree->Branch( "lightside", &lightside , "lightside/I" );
120 
121  tree->SetBranchAddress( "event", &event );
122  tree->SetBranchAddress( "color", &color );
123  tree->SetBranchAddress( "matacq", matacq );
124  tree->SetBranchAddress( "nsamples", &nsamples );
125  tree->SetBranchAddress( "maxsamp", &maxsamp );
126  tree->SetBranchAddress( "tt", &tt );
127  tree->SetBranchAddress( "lightside", &lightside );
128 
129 
130 
131  // Define output results files' names
132 
133  std::stringstream namefile;
134  namefile << resdir_ <<"/MATACQ.root";
135  outfile=namefile.str();
136 
137 
138  // Laser events counter
139  laserEvents=0;
140  isThereMatacq=false;
141 
142 }
143 
144 
145 //========================================================================
146 void EcalMatacqAnalyzer:: analyze( const edm::Event & e, const edm::EventSetup& c){
147 //========================================================================
148 
149  ++iEvent;
150 
151  if (_debug==1 )std::cout << "-- debug test -- Entering Analyze -- event= "<<iEvent<< std::endl;
152 
153  // retrieving MATACQ :
154  edm::Handle<EcalMatacqDigiCollection> pmatacqDigi;
155  const EcalMatacqDigiCollection* matacqDigi=0;
156  try {
157  e.getByLabel(digiProducer_,digiCollection_, pmatacqDigi);
158  matacqDigi=pmatacqDigi.product();
159  if (_debug==1 )std::cout << "-- debug test -- Matacq Digis Found -- "<< std::endl;
160 
161  }catch ( std::exception& ex ) {
162  std::cerr << "Error! can't get the product EcalMatacqDigi producer:" << digiProducer_.c_str()<<" collection:"<<digiCollection_.c_str() << std::endl;
163  if (_debug==1 )std::cout << "-- debug test -- No Matacq Digis Found -- "<< std::endl;
164  return;
165  }
166 
167  // retrieving DCC header
168 
169  edm::Handle<EcalRawDataCollection> pDCCHeader;
170  const EcalRawDataCollection* DCCHeader=0;
171  try {
172  e.getByLabel(eventHeaderProducer_,eventHeaderCollection_, pDCCHeader);
173  DCCHeader=pDCCHeader.product();
174  }catch ( std::exception& ex ) {
175  std::cerr << "Error! can't get the product EcalRawData producer:" << eventHeaderProducer_.c_str()<<" collection:"<<eventHeaderCollection_.c_str() << std::endl;
176  return;
177  }
178 
179  // ====================================
180  // Decode Basic DCCHeader Information
181  // ====================================
182 
183  if (_debug==1) std::cout <<"-- debug test -- Before header -- "<< std::endl;
184 
185  for ( EcalRawDataCollection::const_iterator headerItr= DCCHeader->begin();headerItr != DCCHeader->end();
186  ++headerItr ) {
187 
188  EcalDCCHeaderBlock::EcalDCCEventSettings settings = headerItr->getEventSettings();
189  color = (int) settings.wavelength;
190  if( color<0 ) return;
191 
192  // Get run type and run number
193 
194  int fed = headerItr->fedId();
195 
196  if(fed!=_fedid && _fedid!=-999) continue;
197 
198  runType=headerItr->getRunType();
199  runNum=headerItr->getRunNumber();
200  event=headerItr->getLV1();
201 
202  if (_debug==1) std::cout <<"-- debug test -- runtype:"<<runType<<" event:"<<event <<" runNum:"<<runNum <<std::endl;
203 
204  dccID=headerItr->getDccInTCCCommand();
205  fedID=headerItr->fedId();
206  lightside=headerItr->getRtHalf();
207 
208  //assert (lightside<2 && lightside>=0);
209 
210  if( lightside!=1 && lightside!=0 ) {
211  std::cout << "Unexpected lightside: "<< lightside<<" for event "<<iEvent << std::endl;
212  return;
213  }
214  if (_debug==1) {
215  std::cout <<"-- debug test -- Inside header before fed cut -- color="<<color<< ", dcc="<<dccID<<", fed="<< fedID<<", lightside="<<lightside<<", runType="<<runType<< std::endl;
216  }
217 
218  // take event only if the fed corresponds to the DCC in TCC
219  if( 600+dccID != fedID ) continue;
220 
221  if (_debug==1) {
222  std::cout <<"-- debug test -- Inside header after fed cut -- color="<<color<< ", dcc="<<dccID<<", fed="<< fedID<<", lightside="<<lightside<<", runType="<<runType<< std::endl;
223  }
224 
225  // Cut on runType
226 
227  if ( runType!=EcalDCCHeaderBlock::LASER_STD && runType!=EcalDCCHeaderBlock::LASER_GAP && runType!=EcalDCCHeaderBlock::LASER_POWER_SCAN && runType!=EcalDCCHeaderBlock::LASER_DELAY_SCAN ) return;
228 
229 
230  std::vector<int>::iterator iter= find( colors.begin(), colors.end(), color );
231  if( iter==colors.end() ){
232  colors.push_back( color );
233  std::cout <<" new color found "<< color<<" "<< colors.size()<< std::endl;
234  }
235 
236  }
237 
238  if (_debug==1) std::cout <<"-- debug test -- Before digis -- Event:"<<iEvent<< std::endl;
239 
240  // Count laser events
241  laserEvents++;
242 
243 
244 
245 // ===========================
246 // Decode Matacq Information
247 // ===========================
248 
249  int iCh=0;
250  double max=0;
251 
252  for(EcalMatacqDigiCollection::const_iterator it = matacqDigi->begin(); it!=matacqDigi->end(); ++it){ // Loop on matacq channel
253 
254  //
255  const EcalMatacqDigi& digis = *it;
256 
257  //if(digis.size()==0 || iCh>=N_channels) continue;
258  if (_debug==1) {
259  std::cout <<"-- debug test -- Inside digis -- digi size="<<digis.size()<< std::endl;
260  }
261 
262  if(digis.size()== 0 ) continue;
263  else isThereMatacq=true;
264 
265  max=0;
266  maxsamp=0;
267  nsamples=digis.size();
268  tt=digis.tTrig();
269 
270  for(int i=0; i<digis.size(); ++i){ // Loop on matacq samples
271  matacq[i]=-digis.adcCount(i);
272  if(matacq[i]>max) {
273  max=matacq[i];
274  maxsamp=i;
275  }
276  }
277  if (_debug==1) {
278  std::cout <<"-- debug test -- Inside digis -- nsamples="<<nsamples<< ", max="<<max<< std::endl;
279  }
280 
281  iCh++;
282  }
283 
284  if (_debug==1)std::cout <<"-- debug test -- After digis -- Event: "<<iEvent<< std::endl;
285  tree->Fill();
286 
287 } // analyze
288 
289 
290 //========================================================================
291 void EcalMatacqAnalyzer::endJob()
292 {
293 
294  // Don't do anything if there is no events
295  if( !isThereMatacq ) {
296 
297  std::cout << "\n\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
298  std::cout << "\t+=+ WARNING! NO MATACQ +=+" << std::endl;
299  std::cout << "\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
300 
301 
302  // Remove temporary file
303  FILE *test;
304  test = fopen(sampfile.c_str(),"r");
305  if (test){
306  std::stringstream del2;
307  del2 << "rm " <<sampfile;
308  system(del2.str().c_str());
309  }
310  return;
311  }
312 
313  assert( colors.size()<= nColor );
314  unsigned int nCol=colors.size();
315 
316  for(unsigned int iCol=0;iCol<nCol;iCol++){
317  for(unsigned int iSide=0;iSide<nSide;iSide++){
318  MTQ[iCol][iSide]=new TMTQ();
319  }
320  }
321 
322  outFile = new TFile(outfile.c_str(),"RECREATE");
323 
324  TProfile *shapeMat = new TProfile("shapeLaser","shapeLaser",_nsamplesshape,-0.5,double(_nsamplesshape)-0.5);
325  TProfile *shapeMatTmp = new TProfile("shapeLaserTmp","shapeLaserTmp",_timeaftmax+_timebefmax,-0.5,double(_timeaftmax+_timebefmax)-0.5);
326 
327  std::cout << "\n\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
328  std::cout << "\t+=+ Analyzing MATACQ data +=+" << std::endl;
329 
330 
331  //
332  // create output ntuple
333  //
334 
335  mtqShape = new TTree("MatacqShape","MatacqShape");
336 
337  // list of branches
338  // keep Patrice's notations
339 
340  mtqShape->Branch( "event", &event, "event/I" );
341  mtqShape->Branch( "color", &color, "color/I" );
342  mtqShape->Branch( "status", &status, "status/I" );
343  mtqShape->Branch( "peak", &peak , "peak/D" );
344  mtqShape->Branch( "sigma", &sigma , "sigma/D" );
345  mtqShape->Branch( "fit", &fit , "fit/D" );
346  mtqShape->Branch( "ampl", &ampl , "ampl/D" );
347  mtqShape->Branch( "trise", &trise , "trise/D" );
348  mtqShape->Branch( "fwhm", &fwhm , "fwhm/D" );
349  mtqShape->Branch( "fw20", &fw20 , "fw20/D" );
350  mtqShape->Branch( "fw80", &fw80 , "fw80/D" );
351  mtqShape->Branch( "ped", &ped , "ped/D" );
352  mtqShape->Branch( "pedsig", &pedsig , "pedsig/D" );
353  mtqShape->Branch( "ttrig", &ttrig , "ttrig/D" );
354  mtqShape->Branch( "sliding", &sliding , "sliding/D" );
355 
356  mtqShape->SetBranchAddress( "event", &event );
357  mtqShape->SetBranchAddress( "color", &color );
358  mtqShape->SetBranchAddress( "status", &status );
359  mtqShape->SetBranchAddress( "peak", &peak );
360  mtqShape->SetBranchAddress( "sigma", &sigma );
361  mtqShape->SetBranchAddress( "fit", &fit );
362  mtqShape->SetBranchAddress( "ampl", &ampl );
363  mtqShape->SetBranchAddress( "fwhm", &fwhm );
364  mtqShape->SetBranchAddress( "fw20", &fw20 );
365  mtqShape->SetBranchAddress( "fw80", &fw80 );
366  mtqShape->SetBranchAddress( "trise", &trise );
367  mtqShape->SetBranchAddress( "ped", &ped );
368  mtqShape->SetBranchAddress( "pedsig", &pedsig );
369  mtqShape->SetBranchAddress( "ttrig", &ttrig );
370  mtqShape->SetBranchAddress( "sliding", &sliding );
371 
372 
373  unsigned int endsample;
374  unsigned int presample;
375 
376 
377 
378  // loop over the entries of the tree
379  TChain* fChain = (TChain*)tree;
380  Long64_t nentries = fChain->GetEntriesFast();
381  Long64_t nbytes = 0, nb = 0;
382 
383  for( Long64_t jentry=0; jentry<nentries; jentry++ )
384  {
385  // load the event
386  Long64_t ientry = fChain->LoadTree(jentry);
387  if (ientry < 0) break;
388  nb = fChain->GetEntry( jentry ); nbytes += nb;
389 
390  status = 0;
391  peak = -1;
392  sigma = 0;
393  fit = -1;
394  ampl = -1;
395  trise = -1;
396  ttrig = tt;
397  fwhm = 0;
398  fw20 = 0;
399  fw80 = 0;
400  ped = 0;
401  pedsig = 0;
402  sliding = 0;
403 
404 
405  if (_debug==1)std::cout <<"-- debug test -- inside loop 1 -- jentry:"<<jentry<<" over nentries="<<nentries<< std::endl;
406 
407  // create the object for Matacq data analysis
408 
409  endsample = maxsamp+_nsamplesaftmax;
410  presample=int(_presample*nsamples/100.);
411  TMatacq* mtq = new TMatacq( nsamples, presample, endsample,
412  _noiseCut, _parabnbefmax, _parabnaftmax,
413  _thres, _lowlev, _highlev,
414  _nevlasers , _slide);
415 
416  if (_debug==1)std::cout <<"-- debug test -- inside loop 2 -- "<< std::endl;
417 
418  // analyse the Matacq data
419  if( mtq->rawPulseAnalysis( nsamples, &matacq[0] )==0 )
420  {
421  status = 1;
422  ped = mtq->getBaseLine();
423  pedsig = mtq->getsigBaseLine();
424 
425  if (_debug==1)std::cout <<"-- debug test -- inside loop 3 -- ped:"<<ped << std::endl;
426  if( mtq->findPeak()==0 )
427  {
428  peak = mtq->getTimpeak();
429  sigma = mtq->getsigTimpeak();
430  }
431  if (_debug==1)std::cout <<"-- debug test -- inside loop 4 -- peak:"<<peak<< std::endl;
432  if( mtq->doFit()==0 )
433  {
434  fit = mtq->getTimax();
435  ampl = mtq->getAmpl();
436  fwhm = mtq->getFwhm();
437  fw20 = mtq->getWidth20();
438  fw80 = mtq->getWidth80();
439  sliding = mtq->getSlide();
440  }
441  if (_debug==1)std::cout <<"-- debug test -- inside loop 4 -- ampl:"<<ampl<< std::endl;
442  if( mtq->compute_trise()==0 )
443  {
444  trise = mtq->getTrise();
445  }
446  if (_debug==1)std::cout <<"-- debug test -- inside loop 5 -- trise:"<<trise<< std::endl;
447  }
448 
449  if (_debug==1)std::cout <<"-- debug test -- inside loop 6 -- status:"<<status<< std::endl;
450 
451  if( status == 1 && mtq->findPeak()==0 ){
452 
453  int firstS=int(peak-double(_timebefmax));
454  int lastS=int(peak+double(_timeaftmax));
455 
456 
457  // Fill histo if there are enough samples
458  if (_debug==1)std::cout <<"-- debug test -- inside loop 7 -- firstS:"<<firstS<<", nsamples:"<< nsamples<< std::endl;
459 
460  if(firstS>=0 && lastS<=nsamples){
461 
462  for (int i=firstS;i<lastS;i++){
463  shapeMatTmp->Fill(double(i)-firstS,matacq[i]);
464  }
465 
466 
467  }else{ // else extrapolate
468 
469 
470  int firstSBis;
471 
472  if(firstS<0){ // fill first bins with 0
473 
474  double thisped;
475  thisped=(matacq[0]+matacq[1]+matacq[2]+matacq[4]+matacq[5])/5.0;
476 
477  for(int i=firstS;i<0;i++){
478  shapeMatTmp->Fill(double(i)-firstS,thisped);
479 
480  }
481  firstSBis=0;
482 
483  }else{
484  firstSBis=firstS;
485  }
486 
487  if(lastS>nsamples){
488 
489  for(int i=firstSBis;i<int(nsamples);i++){
490  shapeMatTmp->Fill(double(i)-firstS,matacq[i]);
491  }
492 
493  //extrapolate with expo tail
494 
495  double expb=0.998;
496  double matacqval=expb*matacq[nsamples-1];
497 
498  for(int i=nsamples;i<lastS;i++){
499  shapeMatTmp->Fill(double(i)-firstS,matacqval);
500  matacqval*=expb;
501  }
502 
503  }else{
504  for (int i=firstSBis;i<lastS;i++){
505  shapeMatTmp->Fill(double(i)-firstS,matacq[i]);
506  }
507  }
508  }
509 
510  }
511  if (_debug==1)std::cout <<"-- debug test -- inside loop 8"<< std::endl;
512 
513  // get back color
514 
515  int iCol=nCol;
516  for(unsigned int i=0;i<nCol;i++){
517  if(color==colors[i]) {
518  iCol=i;
519  i=nCol;
520  }
521  }
522  if (_debug==1)std::cout <<"-- debug test -- inside loop 8bis color:"<<color<<" iCol:"<<iCol<<" nCol:"<< nCol<< std::endl;
523 
524  // fill TMTQ
525 
526  if(status == 1 && mtq->findPeak()==0 && mtq->doFit()==0 && mtq->compute_trise()==0 ) MTQ[iCol][lightside]->addEntry(peak, sigma, fit, ampl, trise, fwhm, fw20, fw80, ped, pedsig, sliding);
527 
528  // fill the output tree
529 
530  if (_debug==1)std::cout <<"-- debug test -- inside loop 9"<< std::endl;
531  mtqShape->Fill();
532 
533  // clean up
534  delete mtq;
535  }
536 
537  if (_debug==1)std::cout <<"-- debug test -- after loop "<< std::endl;
538  sampFile->Close();
539 
540  double Peak[6], Sigma[6], Fit[6], Ampl[6], Trise[6], Fwhm[6], Fw20[6], Fw80[6], Ped[6], Pedsig[6], Sliding[6];
541  int Side;
542 
543  for (unsigned int iColor=0;iColor<nCol;iColor++){
544 
545  std::stringstream nametree;
546  nametree <<"MatacqCol"<<colors[iColor];
547  meanTree[iColor]= new TTree(nametree.str().c_str(),nametree.str().c_str());
548  meanTree[iColor]->Branch( "side", &Side , "Side/I" );
549  meanTree[iColor]->Branch( "peak", &Peak , "Peak[6]/D" );
550  meanTree[iColor]->Branch( "sigma", &Sigma , "Sigma[6]/D" );
551  meanTree[iColor]->Branch( "fit", &Fit , "Fit[6]/D" );
552  meanTree[iColor]->Branch( "ampl", &Ampl , "Ampl[6]/D" );
553  meanTree[iColor]->Branch( "trise", &Trise , "Trise[6]/D" );
554  meanTree[iColor]->Branch( "fwhm", &Fwhm , "Fwhm[6]/D" );
555  meanTree[iColor]->Branch( "fw20", &Fw20 , "Fw20[6]/D" );
556  meanTree[iColor]->Branch( "fw80", &Fw80 , "Fw80[6]/D" );
557  meanTree[iColor]->Branch( "ped", &Ped , "Ped[6]/D" );
558  meanTree[iColor]->Branch( "pedsig", &Pedsig , "Pedsig[6]/D" );
559  meanTree[iColor]->Branch( "sliding", &Sliding , "Sliding[6]/D" );
560 
561  meanTree[iColor]->SetBranchAddress( "side", &Side );
562  meanTree[iColor]->SetBranchAddress( "peak", Peak );
563  meanTree[iColor]->SetBranchAddress( "sigma", Sigma );
564  meanTree[iColor]->SetBranchAddress( "fit", Fit );
565  meanTree[iColor]->SetBranchAddress( "ampl", Ampl );
566  meanTree[iColor]->SetBranchAddress( "fwhm", Fwhm );
567  meanTree[iColor]->SetBranchAddress( "fw20", Fw20 );
568  meanTree[iColor]->SetBranchAddress( "fw80", Fw80 );
569  meanTree[iColor]->SetBranchAddress( "trise", Trise );
570  meanTree[iColor]->SetBranchAddress( "ped", Ped );
571  meanTree[iColor]->SetBranchAddress( "pedsig", Pedsig );
572  meanTree[iColor]->SetBranchAddress( "sliding", Sliding );
573 
574  }
575 
576  for(unsigned int iCol=0;iCol<nCol;iCol++){
577  for(unsigned int iSide=0;iSide<nSides;iSide++){
578 
579  Side=iSide;
580  std::vector<double> val[TMTQ::nOutVar];
581 
582  for(int iVar=0;iVar<TMTQ::nOutVar;iVar++){
583  val[iVar] = MTQ[iCol][iSide]->get(iVar);
584 
585  for(unsigned int i=0;i<val[iVar].size();i++){
586 
587  switch (iVar){
588 
589  case TMTQ::iPeak: Peak[i]=val[iVar].at(i);
590  case TMTQ::iSigma: Sigma[i]=val[iVar].at(i);
591  case TMTQ::iFit: Fit[i]=val[iVar].at(i);
592  case TMTQ::iAmpl: Ampl[i]=val[iVar].at(i);
593  case TMTQ::iFwhm: Fwhm[i]=val[iVar].at(i);
594  case TMTQ::iFw20: Fw20[i]=val[iVar].at(i);
595  case TMTQ::iFw80: Fw80[i]=val[iVar].at(i);
596  case TMTQ::iTrise: Trise[i]=val[iVar].at(i);
597  case TMTQ::iPed: Ped[i]=val[iVar].at(i);
598  case TMTQ::iPedsig: Pedsig[i]=val[iVar].at(i);
599  case TMTQ::iSlide: Sliding[i]=val[iVar].at(i);
600  }
601  }
602  }
603  meanTree[iCol]->Fill();
604  if (_debug==1)std::cout <<"-- debug test -- inside final loop "<< std::endl;
605  }
606  }
607 
608  // Calculate maximum with pol 2
609 
610  int im = shapeMatTmp->GetMaximumBin();
611  double q1=shapeMatTmp->GetBinContent(im-1);
612  double q2=shapeMatTmp->GetBinContent(im);
613  double q3=shapeMatTmp->GetBinContent(im+1);
614 
615  double a2=(q3+q1)/2.0-q2;
616  double a1=q2-q1+a2*(1-2*im);
617  double a0=q2-a1*im-a2*im*im;
618 
619  double am=a0-a1*a1/(4*a2);
620 
621 
622  // Compute pedestal
623 
624  double bl=0;
625  for (unsigned int i=1; i<_presampleshape+1;i++){
626  bl+=shapeMatTmp->GetBinContent(i);
627  }
628  bl/=_presampleshape;
629 
630  // Compute and save laser shape
631 
632  if (_debug==1)std::cout <<"-- debug test -- computing shape "<< std::endl;
633 
634  int firstBin=0;
635  double height=0.0;
636 
637  for (unsigned int i=_timebefmax; i>_presampleshape;i--){
638  height=shapeMatTmp->GetBinContent(i)-bl;
639 
640  if(height<(am-bl)*_cutwindow){
641  firstBin=i;
642  i=_presampleshape;
643  }
644  }
645 
646  unsigned int lastBin=firstBin+_nsamplesshape;
647 
648  for(unsigned int i=firstBin;i<lastBin;i++){
649  shapeMat->Fill(i-firstBin,shapeMatTmp->GetBinContent(i)-bl);
650  }
651 
652  mtqShape->Write();
653  for (unsigned int iColor=0;iColor<nCol;iColor++){
654  meanTree[iColor]->Write();
655  }
656  if (_debug==1)std::cout <<"-- debug test -- writing "<< std::endl;
657  shapeMat->Write();
658 
659  // close the output file
660  outFile->Close();
661 
662  // Remove temporary file
663  FILE *test;
664  test = fopen(sampfile.c_str(),"r");
665  if (test){
666  std::stringstream del2;
667  del2 << "rm " <<sampfile;
668  system(del2.str().c_str());
669  }
670 
671  std::cout << "\t+=+ .................... done +=+" << std::endl;
672  std::cout << "\t+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+=+" << std::endl;
673 
674 }
675 
676 DEFINE_FWK_MODULE(EcalMatacqAnalyzer);
677 
TMatacq::getWidth80
double getWidth80()
Definition: TMatacq.h:76
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
i
int i
Definition: DBlmapReader.cc:9
TMTQ
Definition: TMTQ.h:8
NSIDES
#define NSIDES
Definition: EcalLaserAnalyzer.h:44
EcalMatacqAnalyzer::_presampleshape
unsigned int _presampleshape
Definition: EcalMatacqAnalyzer.h:50
EcalMatacqAnalyzer::analyze
virtual void analyze(const edm::Event &e, const edm::EventSetup &c)
Definition: EcalMatacqAnalyzer.cc:146
EcalMatacqAnalyzer::EcalMatacqAnalyzer
EcalMatacqAnalyzer(const edm::ParameterSet &iConfig)
Definition: EcalMatacqAnalyzer.cc:39
EcalMatacqAnalyzer::colors
std::vector< int > colors
Definition: EcalMatacqAnalyzer.h:93
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
EcalMatacqAnalyzer::_nsamplesshape
unsigned int _nsamplesshape
Definition: EcalMatacqAnalyzer.h:49
edm::SortedCollection< EcalDCCHeaderBlock >::const_iterator
std::vector< EcalDCCHeaderBlock >::const_iterator const_iterator
Definition: SortedCollection.h:82
EcalMatacqDigi::tTrig
float tTrig() const
Definition: EcalMatacqDigi.h:106
TMatacq::compute_trise
int compute_trise()
Definition: TMatacq.cc:330
TMatacq::getTimax
double getTimax()
Definition: TMatacq.h:71
EcalMatacqAnalyzer::eventHeaderProducer_
std::string eventHeaderProducer_
Definition: EcalMatacqAnalyzer.h:59
EcalMatacqAnalyzer::MTQ
TMTQ * MTQ[nColor][nSide]
Definition: EcalMatacqAnalyzer.h:90
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
TMatacq::getBaseLine
double getBaseLine()
Definition: TMatacq.h:64
TMatacq::getsigBaseLine
double getsigBaseLine()
Definition: TMatacq.h:65
q2
double q2[4]
Definition: TauolaWrapper.h:88
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
TMatacq::getFwhm
double getFwhm()
Definition: TMatacq.h:74
EcalMatacqAnalyzer::meanTree
TTree * meanTree[nColor]
Definition: EcalMatacqAnalyzer.h:91
TMatacq::getsigTimpeak
double getsigTimpeak()
Definition: TMatacq.h:68
EcalMatacqDigi::adcCount
const float adcCount(const int &i) const
Definition: EcalMatacqDigi.h:54
EcalMatacqAnalyzer::_timebefmax
unsigned int _timebefmax
Definition: EcalMatacqAnalyzer.h:46
iEvent
int iEvent
Definition: GenABIO.cc:230
EcalMatacqAnalyzer::eventHeaderCollection_
std::string eventHeaderCollection_
Definition: EcalMatacqAnalyzer.h:58
TMatacq::getSlide
double getSlide()
Definition: TMatacq.h:77
Fit
Definition: Fit.h:34
TMTQ::get
std::vector< double > get(int)
Definition: TMTQ.cc:71
TMatacq::doFit
int doFit()
Definition: TMatacq.cc:201
TMatacq::getAmpl
double getAmpl()
Definition: TMatacq.h:70
TMatacq::rawPulseAnalysis
int rawPulseAnalysis(Int_t, Double_t *)
Definition: TMatacq.cc:59
event
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::Event::getByLabel
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:390
EcalMatacqAnalyzer::beginJob
virtual void beginJob()
Definition: EcalMatacqAnalyzer.cc:95
edm::SortedCollection::end
const_iterator end() const
Definition: SortedCollection.h:310
EcalMatacqAnalyzer::_parabnaftmax
unsigned int _parabnaftmax
Definition: EcalMatacqAnalyzer.h:41
environment_file_cff.runType
tuple runType
Definition: environment_file_cff.py:30
TMatacq::getWidth20
double getWidth20()
Definition: TMatacq.h:75
TMatacq::findPeak
int findPeak()
Definition: TMatacq.cc:127
q1
double q1[4]
Definition: TauolaWrapper.h:87
run_regression.test
tuple test
Definition: run_regression.py:386
edm::Handle::product
T const * product() const
Definition: Handle.h:81
dtNoiseDBValidation_cfg.cerr
tuple cerr
Definition: dtNoiseDBValidation_cfg.py:22
TMatacq::getTrise
double getTrise()
Definition: TMatacq.h:73
create_public_lumi_plots.color
list color
Definition: create_public_lumi_plots.py:1800
EcalMatacqAnalyzer::_timeaftmax
unsigned int _timeaftmax
Definition: EcalMatacqAnalyzer.h:47
EcalMatacqDigi::size
int size() const
Definition: EcalMatacqDigi.h:76
TMatacq
Definition: TMatacq.h:9
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
EcalMatacqAnalyzer::digiCollection_
std::string digiCollection_
Definition: EcalMatacqAnalyzer.h:56
EcalMatacqAnalyzer::_nsamplesaftmax
unsigned int _nsamplesaftmax
Definition: EcalMatacqAnalyzer.h:38
edm::SortedCollection::begin
const_iterator begin() const
Definition: SortedCollection.h:303
TMatacq::getTimpeak
double getTimpeak()
Definition: TMatacq.h:67
EcalMatacqAnalyzer::_parabnbefmax
unsigned int _parabnbefmax
Definition: EcalMatacqAnalyzer.h:40
Generated for CMSSW Reference Manual by doxygen 1.8.5