CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HOCalibAnalyzer.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: HOCalibAnalyzer
4 // Class: HOCalibAnalyzer
5 //
13 //
14 // Original Author: Gobinda Majumder
15 // Created: Sat Jul 7 09:51:31 CEST 2007
16 // $Id: HOCalibAnalyzer.cc,v 1.12 2012/09/27 20:32:04 wdd Exp $
17 // $Id: HOCalibAnalyzer.cc,v 1.12 2012/09/27 20:32:04 wdd Exp $
18 //
19 //
20 
21 
22 // system include files
23 #include <memory>
24 
25 // user include files
28 
31 
33 //#include "DataFormats/HOCalibHit/interface/HOCalibVariables.h"
35 
39 
40 #include "TMath.h"
41 #include "TFile.h"
42 #include "TH1F.h"
43 #include "TH2F.h"
44 #include "TTree.h"
45 #include "TProfile.h"
46 #include "TPostScript.h"
47 #include "TCanvas.h"
48 #include "TF1.h"
49 #include "TStyle.h"
50 #include "TMinuit.h"
51 #include "TMath.h"
52 
53 #include <string>
54 
55 #include <iostream>
56 #include <fstream>
57 #include <iomanip>
58 //#include <sstream>
59 
60 
61 using namespace std;
62 using namespace edm;
63 
64 
65  //
66  // Look for nearby pixel through eta, phi informations for pixel cross-talk
67  // 1. Look PIXEL code from (eta,phi)
68  // 2. Go to nearby pixel code
69  // 3. Come back to (eta,phi) from pixel code
70  // Though it works, it is a very ugly/crude way to get cross talk, need better algorithms
71  //
72 
73 static const int mapx1[6][3]={{1,4,8}, {12,7,3}, {5,9,13}, {11,6,2}, {16,15,14}, {19,18,17}};
74 
75 static const int mapx2[6][3]={{1,4,8}, {12,7,3}, {5,9,13}, {11,6,2}, {16,15,14}, {-1,-1,-1}};
76 
77 static const int mapx0p[9][2]={{3,1}, {7,4}, {6,5}, {12,8}, {0,0}, {11,9}, {16,13}, {15,14}, {19,17}};
78 static const int mapx0m[9][2]={{17,19}, {14,15}, {13,16}, {9,11}, {0,0}, {8,12}, {5,6}, {4,7}, {1,3}};
79 
80 static const int etamap[4][21]={{-1, 0,3,1, 0,2,3, 1,0,2, -1, 3,1,2, 4,4,4, -1,-1,-1, -1}, //etamap2
81  {-1, 0,3,1, 0,2,3, 1,0,2, -1, 3,1,2, 4,4,4, 5,5,5, -1}, //etamap1
82  {-1, 0,-1,0, 1,2,2, 1,3,5, -1, 5,3,6, 7,7,6, 8,-1,8, -1}, //etamap0p
83  {-1, 8,-1,8, 7,6,6, 7,5,3, -1, 3,5,2, 1,1,2, 0,-1,0, -1}}; //etamap0m
84 
85 static const int phimap[4][21] ={{-1, 0,2,2, 1,0,1, 1,2,1, -1, 0,0,2, 2,1,0, 2,1,0, -1}, //phimap2
86  {-1, 0,2,2, 1,0,1, 1,2,1, -1, 0,0,2, 2,1,0, 2,1,0, -1}, //phimap1
87  {-1, 1,-1,0, 1,1,0, 0,1,1, -1, 0,0,1, 1,0,0, 1,-1,0, -1}, //phimap0p
88  {-1, 0,-1,1, 0,0,1, 1,0,0, -1, 1,1,0, 0,1,1, 0,-1,1, -1}}; //phimap0m
89 //swapped phi map for R0+/R0- (15/03/07)
90 
91 static const int npixleft[21]={0, 0, 1, 2, 0, 4, 5, 6, 0, 8, 0, 0,11, 0,13,14,15, 0,17,18,0};
92 static const int npixrigh[21]={0, 2, 3, 0, 5, 6, 7, 0, 9, 0, 0,12, 0,14,15,16, 0,18,19, 0,0};
93 static const int npixlebt[21]={0, 0, 0, 0, 0, 1, 2, 3, 0, 4, 0, 6, 7, 8, 9, 0,11,13,14,15,0};
94 static const int npixribt[21]={0, 0, 0, 0, 1, 2, 3, 0, 4, 5, 0, 7, 0, 9, 0,11,12,14,15,16,0};
95 static const int npixleup[21]={0, 4, 5, 6, 8, 9, 0,11, 0,13, 0,15,16, 0,17,18,19, 0, 0, 0,0};
96 static const int npixriup[21]={0, 5, 6, 7, 9, 0,11,12,13,14, 0,16, 0,17,18,19, 0, 0, 0, 0,0};
97 
98 static const int netamx=30;
99 static const int nphimx =72;
100  const int nbgpr = 3;
101  const int nsgpr = 7;
102 
105 vector<float>sig_reg[netamx][nphimx+1];
106 vector<float>cro_ssg[netamx][nphimx+1];
107 
108 
109 //#define CORREL
110 //
111 // class decleration
112 //
113 
114 Double_t gausX(Double_t* x, Double_t* par){
115  return par[0]*(TMath::Gaus(x[0], par[1], par[2], kTRUE));
116 }
117 
118 Double_t langaufun(Double_t *x, Double_t *par) {
119 
120  //Fit parameters:
121  //par[0]*par[1]=Width (scale) parameter of Landau density
122  //par[1]=Most Probable (MP, location) parameter of Landau density
123  //par[2]=Total area (integral -inf to inf, normalization constant)
124  //par[3]=Width (sigma) of convoluted Gaussian function
125  //
126  //In the Landau distribution (represented by the CERNLIB approximation),
127  //the maximum is located at x=-0.22278298 with the location parameter=0.
128  //This shift is corrected within this function, so that the actual
129  //maximum is identical to the MP parameter.
130  // /*
131  // Numeric constants
132  Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
133  Double_t mpshift = -0.22278298; // Landau maximum location
134 
135  // Control constants
136  Double_t np = 100.0; // number of convolution steps
137  Double_t sc = 5.0; // convolution extends to +-sc Gaussian sigmas
138 
139  // Variables
140  Double_t xx;
141  Double_t mpc;
142  Double_t fland;
143  Double_t sum = 0.0;
144  Double_t xlow,xupp;
145  Double_t step;
146  Double_t i;
147 
148 
149  // MP shift correction
150  mpc = par[1] - mpshift * par[0]*par[1];
151 
152  // Range of convolution integral
153  xlow = x[0] - sc * par[3];
154  xupp = x[0] + sc * par[3];
155 
156  step = (xupp-xlow) / np;
157 
158  // Convolution integral of Landau and Gaussian by sum
159  for(i=1.0; i<=np/2; i++) {
160  xx = xlow + (i-.5) * step;
161  fland = TMath::Landau(xx,mpc,par[0]*par[1], kTRUE); // / par[0];
162  sum += fland * TMath::Gaus(x[0],xx,par[3]);
163  xx = xupp - (i-.5) * step;
164  fland = TMath::Landau(xx,mpc,par[0]*par[1], kTRUE); // / par[0];
165  sum += fland * TMath::Gaus(x[0],xx,par[3]);
166  }
167 
168  return (par[2] * step * sum * invsq2pi / par[3]);
169 }
170 
171 Double_t totalfunc(Double_t* x, Double_t* par){
172  return gausX(x, par) + langaufun(x, &par[3]);
173 }
174 
175 void fcnbg(Int_t &npar, Double_t* gin, Double_t &f, Double_t* par, Int_t flag) {
176 
177  double fval = -par[0];
178  for (unsigned i=0; i<cro_ssg[ietafit][iphifit].size(); i++) {
179  double xval = (double)cro_ssg[ietafit][iphifit][i];
180  fval +=log(max(1.e-30,par[0]*TMath::Gaus(xval, par[1], par[2], 1)));
181  // fval +=log(par[0]*TMath::Gaus(xval, par[1], par[2], 1));
182  }
183  f = -fval;
184 }
185 
186 void fcnsg(Int_t &npar, Double_t* gin, Double_t &f, Double_t* par, Int_t flag) {
187 
188  double xval[2];
189  double fval = -(par[0]+par[5]);
190  for (unsigned i=0; i<sig_reg[ietafit][iphifit].size(); i++) {
191  xval[0] = (double)sig_reg[ietafit][iphifit][i];
192  fval +=log(totalfunc(xval, par));
193  }
194  f = -fval;
195 }
196 
197 void set_mean(double& x, bool mdigi) {
198  if(mdigi) {
199  x = min(x, 0.5);
200  x = max(x, -0.5);
201  } else {
202  x = min(x, 0.1);
203  x = max(x, -0.1);
204  }
205 }
206 
207 void set_sigma(double& x, bool mdigi) {
208  if(mdigi) {
209  x = min(x, 1.2);
210  x = max(x, -1.2);
211  } else {
212  x = min(x, 0.24);
213  x = max(x, -0.24);
214  }
215 }
216 
217 
218 
220  public:
221  explicit HOCalibAnalyzer(const edm::ParameterSet&);
222  ~HOCalibAnalyzer();
223 
224 
225  private:
226 
227  virtual void beginJob() ;
228  virtual void analyze(const edm::Event&, const edm::EventSetup&);
229  virtual void endJob() ;
230 
231  TFile* theFile;
235 
236  bool m_hotime;
237  bool m_hbtime;
238  bool m_correl;
240  bool m_hbinfo;
243  bool m_figure;
244  bool m_digiInput;
245  bool m_cosmic;
246  bool m_histfit;
248  double m_sigma;
249 
250  static const int ncut = 13;
251  static const int mypow_2_ncut = 8192; // 2^13, should be changed to match ncut
252 
253  int ipass;
254 
255  TTree* T1;
256 
257  TH1F* muonnm;
258  TH1F* muonmm;
259  TH1F* muonth;
260  TH1F* muonph;
261  TH1F* muonch;
262 
263  TH1F* sel_muonnm;
264  TH1F* sel_muonmm;
265  TH1F* sel_muonth;
266  TH1F* sel_muonph;
267  TH1F* sel_muonch;
268 
269  // if (m_hotime) { // #ifdef HOTIME
270  TProfile* hotime[netamx][nphimx];
271  TProfile* hopedtime[netamx][nphimx];
272  // } // if (m_hotime) #endif
273  // if (m_hbtime) { // #ifdef HBTIME
274  TProfile* hbtime[netamx][nphimx];
275  // } // m_hbtime #endif
276 
277  // double pedval[netamx][nphimx]; // ={0};
278 
279  // if (m_correl) { // #ifdef CORREL
280  TH1F* corrsglb[netamx][nphimx];
281  TH1F* corrsgrb[netamx][nphimx];
282  TH1F* corrsglu[netamx][nphimx];
283  TH1F* corrsgru[netamx][nphimx];
284  TH1F* corrsgall[netamx][nphimx];
285 
286  TH1F* corrsgl[netamx][nphimx];
287  TH1F* corrsgr[netamx][nphimx];
288 
289  TH1F* mncorrsglb;
290  TH1F* mncorrsgrb;
291  TH1F* mncorrsglu;
292  TH1F* mncorrsgru;
293  TH1F* mncorrsgall;
294 
295  TH1F* mncorrsgl;
296  TH1F* mncorrsgr;
297 
298  TH1F* rmscorrsglb;
299  TH1F* rmscorrsgrb;
300  TH1F* rmscorrsglu;
301  TH1F* rmscorrsgru;
303 
304  TH1F* rmscorrsgl;
305  TH1F* rmscorrsgr;
306 
307  TH1F* nevcorrsglb;
308  TH1F* nevcorrsgrb;
309  TH1F* nevcorrsglu;
310  TH1F* nevcorrsgru;
312 
313  TH1F* nevcorrsgl;
314  TH1F* nevcorrsgr;
315 
316  // } //m_correl #endif
317  // if (m_checkmap) { // #ifdef CHECKMAP
318  TH1F* corrsgc[netamx][nphimx];
319  TH1F* mncorrsgc;
320  TH1F* rmscorrsgc;
321  TH1F* nevcorrsgc;
322 
323  // } //m_checkmap #endif
324 
325  TH1F* sigrsg[netamx][nphimx+1];
326  TH1F* crossg[netamx][nphimx+1];
327  float invang[netamx][nphimx+1];
328 
329  TH1F* mnsigrsg;
330  TH1F* mncrossg;
331 
332  TH1F* rmssigrsg;
333  TH1F* rmscrossg;
334 
335  TH1F* nevsigrsg;
336  TH1F* nevcrossg;
337 
338  TH1F* ho_sig2p[9];
339  TH1F* ho_sig1p[9];
340  TH1F* ho_sig00[9];
341  TH1F* ho_sig1m[9];
342  TH1F* ho_sig2m[9];
343 
344  // if (m_hbinfo) { // #ifdef HBINFO
345  TH1F* hbhe_sig[9];
346  // } // m_hbinfo #endif
347 
348 
349  // if (m_combined) { //#ifdef COMBINED
350  static const int ringmx=5;
351  static const int sectmx=12;
352  static const int routmx=36;
353  static const int rout12mx=24;
354  static const int neffip=6;
355 
356  // if (m_hotime) { // #ifdef HOTIME
357  TProfile* com_hotime[ringmx][sectmx];
358  TProfile* com_hopedtime[ringmx][sectmx];
359  // } //m_hotime #endif
360  // if (m_hbtime) { #ifdef HBTIME
361  TProfile* com_hbtime[ringmx][sectmx];
362  // } // m_hbtime #endif
363  // if (m_correl) { //#ifdef CORREL
364  TH1F* com_corrsglb[ringmx][sectmx];
365  TH1F* com_corrsgrb[ringmx][sectmx];
366  TH1F* com_corrsglu[ringmx][sectmx];
367  TH1F* com_corrsgru[ringmx][sectmx];
368  TH1F* com_corrsgall[ringmx][sectmx];
369 
370  TH1F* com_corrsgl[ringmx][sectmx];
371  TH1F* com_corrsgr[ringmx][sectmx];
372  // } //m_correl #endif
373  // if (m_checkmap) { #ifdef CHECKMAP
374  TH1F* com_corrsgc[ringmx][sectmx];
375  // } // m_checkmap #endif
376 
377  TH1F* com_sigrsg[ringmx][routmx+1];
378  TH1F* com_crossg[ringmx][routmx+1];
379  float com_invang[ringmx][routmx+1];
380 
381 
382  // } //m_combined #endif
383 
384  TH1F* ped_evt;
385  TH1F* ped_mean;
386  TH1F* ped_width;
387  TH1F* fit_chi;
388  TH1F* sig_evt;
389  TH1F* fit_sigevt;
390  TH1F* fit_bkgevt;
391  TH1F* sig_mean;
392  TH1F* sig_diff;
393  TH1F* sig_width;
394  TH1F* sig_sigma;
395  TH1F* sig_meanerr;
397  TH1F* sig_signf;
398 
401  TH1F* ped_rms;
402  TH1F* sig_rms;
403 
405 
406  TH1F* const_eta[netamx];
407  TH1F* stat_eta[netamx];
408  TH1F* statmn_eta[netamx];
409  TH1F* peak_eta[netamx];
410 
411  TH1F* const_hpdrm[ringmx];
412  // TH1F* stat_hpdrm[ringmx];
413  // TH1F* statmn_hpdrm[ringmx];
414  TH1F* peak_hpdrm[ringmx];
415 
419 
420  TH2F* sig_effi[neffip];
421  TH2F* mean_energy;
422 
423  double fitprm[nsgpr][netamx];
424 
425 
426  TProfile* sigvsevt[15][ncut];
427 
428 
429  // int irun, ievt, itrg1, itrg2, isect, nrecht, nfound, nlost, ndof, nmuon;
430  int irun, ievt, itrg1, itrg2, isect, ndof, nmuon;
431  float trkdr, trkdz, trkvx, trkvy, trkvz, trkmm, trkth, trkph, chisq, therr, pherr, hodx, hody,
432  hoang, htime, hosig[9], hocorsig[18], hocro, hbhesig[9], caloen[3];
433 
434  int Nevents;
435  int nbn;
436  float alow;
437  float ahigh;
438  float binwid;
439  int irunold;
440 
442 
443  // ----------member data ---------------------------
444 
445 };
446 
447 const int HOCalibAnalyzer::ringmx;
448 const int HOCalibAnalyzer::sectmx;
449 const int HOCalibAnalyzer::routmx;
450 const int HOCalibAnalyzer::rout12mx;
451 const int HOCalibAnalyzer::neffip;
452 
453 //
454 // constants, enums and typedefs
455 //
456 
457 //
458 // static data member definitions
459 //
460 
461 //
462 // constructors and destructor
463 //
465  hoCalibVariableCollectionTag(iConfig.getParameter<edm::InputTag>("hoCalibVariableCollectionTag"))
466  // It is very likely you want the following in your configuration
467  // hoCalibVariableCollectionTag = cms.InputTag('hoCalibProducer', 'HOCalibVariableCollection')
468 {
469  //now do what ever initialization is needed
470  ipass = 0;
471  Nevents = 0;
472 
473  theRootFileName = iConfig.getUntrackedParameter<string>("RootFileName", "test.root");
474  theoutputtxtFile = iConfig.getUntrackedParameter<string>("txtFileName", "test.txt");
475  theoutputpsFile = iConfig.getUntrackedParameter<string>("psFileName", "test.ps");
476 
477  m_hbinfo = iConfig.getUntrackedParameter<bool>("hbinfo", false);
478  m_hbtime = iConfig.getUntrackedParameter<bool>("hbtime", false);
479  m_hotime = iConfig.getUntrackedParameter<bool>("hotime", false);
480  m_correl = iConfig.getUntrackedParameter<bool>("correl", false);
481  m_checkmap = iConfig.getUntrackedParameter<bool>("checkmap", false);
482  m_combined = iConfig.getUntrackedParameter<bool>("combined", false);
483  m_constant = iConfig.getUntrackedParameter<bool>("get_constant", false);
484  m_figure = iConfig.getUntrackedParameter<bool>("get_figure", true);
485  m_digiInput = iConfig.getUntrackedParameter<bool>("digiInput", true);
486  m_histfit = iConfig.getUntrackedParameter<bool>("histFit", true);
487  m_pedsuppr = iConfig.getUntrackedParameter<bool>("pedSuppr", true);
488  m_cosmic = iConfig.getUntrackedParameter<bool>("cosmic", true);
489  m_sigma = iConfig.getUntrackedParameter<double>("sigma", 1.0);
490 
492 
493  theFile = new TFile(theRootFileName.c_str(), "RECREATE");
494  theFile->cd();
495 
496  T1 = new TTree("T1", "DT+CSC+HO");
497 
498  T1->Branch("irun",&irun,"irun/I");
499  T1->Branch("ievt",&ievt,"ievt/I");
500  T1->Branch("itrg1",&itrg1,"itrg1/I");
501  T1->Branch("itrg2",&itrg2,"itrg2/I");
502 
503  T1->Branch("isect",&isect,"isect/I");
504  // T1->Branch("nrecht",&nrecht,"nrecht/I");
505  T1->Branch("ndof",&ndof,"ndof/I");
506  T1->Branch("nmuon",&nmuon,"nmuon/I");
507 
508  T1->Branch("trkdr",&trkdr,"trkdr/F");
509  T1->Branch("trkdz",&trkdz,"trkdz/F");
510 
511  T1->Branch("trkvx",&trkvx,"trkvx/F");
512  T1->Branch("trkvy",&trkvy,"trkvy/F");
513  T1->Branch("trkvz",&trkvz,"trkvz/F");
514  T1->Branch("trkmm",&trkmm,"trkmm/F");
515  T1->Branch("trkth",&trkth,"trkth/F");
516  T1->Branch("trkph",&trkph,"trkph/F");
517 
518  T1->Branch("chisq",&chisq,"chisq/F");
519  T1->Branch("therr",&therr,"therr/F");
520  T1->Branch("pherr",&pherr,"pherr/F");
521  T1->Branch("hodx",&hodx,"hodx/F");
522  T1->Branch("hody",&hody,"hody/F");
523  T1->Branch("hoang",&hoang,"hoang/F");
524 
525  T1->Branch("htime",&htime,"htime/F");
526  T1->Branch("hosig",hosig,"hosig[9]/F");
527  T1->Branch("hocro",&hocro,"hocro/F");
528  T1->Branch("hocorsig",hocorsig,"hocorsig[18]/F");
529  T1->Branch("caloen",caloen,"caloen[3]/F");
530 
531  if (m_hbinfo) { // #ifdef HBINFO
532  T1->Branch("hbhesig",hbhesig,"hbhesig[9]/F");
533  } //m_hbinfo #endif
534 
535  muonnm = fs->make<TH1F>("muonnm", "No of muon", 10, -0.5, 9.5);
536  muonmm = fs->make<TH1F>("muonmm", "P_{mu}", 200, -100., 100.);
537  muonth = fs->make<TH1F>("muonth", "{Theta}_{mu}", 180, 0., 180.);
538  muonph = fs->make<TH1F>("muonph", "{Phi}_{mu}", 180, -180., 180.);
539  muonch = fs->make<TH1F>("muonch", "{chi^2}/ndf", 100, 0., 1000.);
540 
541  sel_muonnm = fs->make<TH1F>("sel_muonnm", "No of muon(sel)", 10, -0.5, 9.5);
542  sel_muonmm = fs->make<TH1F>("sel_muonmm", "P_{mu}(sel)", 200, -100., 100.);
543  sel_muonth = fs->make<TH1F>("sel_muonth", "{Theta}_{mu}(sel)", 180, 0., 180.);
544  sel_muonph = fs->make<TH1F>("sel_muonph", "{Phi}_{mu}(sel)", 180, -180., 180.);
545  sel_muonch = fs->make<TH1F>("sel_muonch", "{chi^2}/ndf(sel)", 100, 0., 1000.);
546 
547 
548  int nbin = 50; //40;// 45; //50; //55; //60; //55; //45; //40; //50;
549  alow = -2.0;// -1.85; //-1.90; // -1.95; // -2.0;
550  ahigh = 8.0;// 8.15; // 8.10; // 8.05; // 8.0;
551 
552  if (m_digiInput) {alow = -10.0; ahigh = 40.;}
553 
554  float tmpwid = (ahigh-alow)/nbin;
555  nbn = int(-alow/tmpwid)+1;
556  if (nbn <0) nbn = 0;
557  if (nbn>nbin) nbn = nbin;
558 
559  char name[200];
560  char title[200];
561 
562  cout <<"nbin "<< nbin<<" "<<alow<<" "<<ahigh<<" "<<tmpwid<<" "<<nbn<<endl;
563 
564  for (int i=0; i<15; i++) {
565 
566  sprintf(title, "sigvsndof_ring%i", i+1);
567  sigvsevt[i][0] = fs->make<TProfile>(title, title, 50, 0., 50.,-9., 20.);
568 
569  sprintf(title, "sigvschisq_ring%i", i+1);
570  sigvsevt[i][1] = fs->make<TProfile>(title, title, 50, 0., 30.,-9., 20.);
571 
572  sprintf(title, "sigvsth_ring%i", i+1);
573  sigvsevt[i][2] = fs->make<TProfile>(title, title, 50, .7, 2.4,-9., 20.);
574 
575  sprintf(title, "sigvsph_ring%i", i+1);
576  sigvsevt[i][3] = fs->make<TProfile>(title, title, 50, -2.4, -0.7,-9., 20.);
577 
578  sprintf(title, "sigvstherr_ring%i", i+1);
579  sigvsevt[i][4] = fs->make<TProfile>(title, title, 50, 0., 0.2,-9., 20.);
580 
581  sprintf(title, "sigvspherr_ring%i", i+1);
582  sigvsevt[i][5] = fs->make<TProfile>(title, title, 50, 0., 0.2,-9., 20.);
583 
584  sprintf(title, "sigvsdircos_ring%i", i+1);
585  sigvsevt[i][6] = fs->make<TProfile>(title, title, 50, 0.5, 1.,-9., 20.);
586 
587  sprintf(title, "sigvstrkmm_ring%i", i+1);
588  sigvsevt[i][7] = fs->make<TProfile>(title, title, 50, 0., 50.,-9., 20.);
589 
590  sprintf(title, "sigvsnmuon_ring%i", i+1);
591  sigvsevt[i][8] = fs->make<TProfile>(title, title, 5, 0.5, 5.5,-9., 20.);
592 
593  sprintf(title, "sigvserr_ring%i", i+1);
594  sigvsevt[i][9] = fs->make<TProfile>(title, title, 50, 0., .3, -9., 20.);
595 
596  sprintf(title, "sigvsaccx_ring%i", i+1);
597  sigvsevt[i][10] = fs->make<TProfile>(title, title, 100, -25., 25., -9., 20.);
598 
599  sprintf(title, "sigvsaccy_ring%i", i+1);
600  sigvsevt[i][11] = fs->make<TProfile>(title, title, 100, -25., 25., -9., 20.);
601 
602  sprintf(title, "sigvscalo_ring%i", i+1);
603  sigvsevt[i][12] = fs->make<TProfile>(title, title, 100, 0., 15., -9., 20.);
604  }
605 
606  for (int j=0; j<netamx; j++) {
607  int ieta = (j<15) ? j+1 : 14-j;
608  for (int i=0;i<nphimx+1;i++) {
609  if (i==nphimx) {
610  sprintf(title, "sig_eta%i_allphi", ieta);
611  } else {
612  sprintf(title, "sig_eta%i_phi%i", ieta,i+1);
613  }
614  sigrsg[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
615  if (i==nphimx) {
616  sprintf(title, "ped_eta%i_allphi", ieta);
617  } else {
618  sprintf(title, "ped_eta%i_phi%i", ieta,i+1);
619  }
620  crossg[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
621  }
622 
623  for (int i=0;i<nphimx;i++) {
624  if (m_hotime) { //#ifdef HOTIME
625  sprintf(title, "hotime_eta%i_phi%i", (j<=14) ? j+1 : 14-j, i+1);
626  hotime[j][i] = fs->make<TProfile>(title, title, 10, -0.5, 9.5, -1.0, 30.0);
627 
628  sprintf(title, "hopedtime_eta%i_phi%i", (j<=14) ? j+1 : 14-j, i+1);
629  hopedtime[j][i] = fs->make<TProfile>(title, title, 10, -0.5, 9.5, -1.0, 30.0);
630 
631  } //m_hotime #endif
632  if (m_hbtime) { //#ifdef HBTIME
633  sprintf(title, "hbtime_eta%i_phi%i", (j<=15) ? j+1 : 15-j, i+1);
634  hbtime[j][i] = fs->make<TProfile>(title, title, 10, -0.5, 9.5, -1.0, 30.0);
635  } //m_hbtime #endif
636 
637  if (m_correl) { //#ifdef CORREL
638  sprintf(title, "corrsg_eta%i_phi%i_leftbottom", ieta,i+1);
639  corrsglb[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
640 
641  sprintf(title, "corrsg_eta%i_phi%i_rightbottom", ieta,i+1);
642  corrsgrb[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
643 
644  sprintf(title, "corrsg_eta%i_phi%i_leftup", ieta,i+1);
645  corrsglu[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
646 
647  sprintf(title, "corrsg_eta%i_phi%i_rightup", ieta,i+1);
648  corrsgru[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
649 
650  sprintf(title, "corrsg_eta%i_phi%i_all", ieta,i+1);
651  corrsgall[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
652 
653  sprintf(title, "corrsg_eta%i_phi%i_left", ieta,i+1);
654  corrsgl[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
655 
656  sprintf(title, "corrsg_eta%i_phi%i_right", ieta,i+1);
657  corrsgr[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
658  } //m_correl #endif
659  if (m_checkmap) {// #ifdef CHECKMAP
660  sprintf(title, "corrsg_eta%i_phi%i_centrl", ieta,i+1);
661  corrsgc[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
662  } //m_checkmap #endif
663  }
664  }
665 
666  mnsigrsg = fs->make<TH1F>("mnsigrsg","mnsigrsg", netamx*nphimx+ringmx*routmx, -0.5, netamx*nphimx+ringmx*routmx -0.5);
667  rmssigrsg = fs->make<TH1F>("rmssigrsg","rmssigrsg", netamx*nphimx+ringmx*routmx, -0.5, netamx*nphimx+ringmx*routmx -0.5);
668  nevsigrsg = fs->make<TH1F>("nevsigrsg","nevsigrsg", netamx*nphimx+ringmx*routmx, -0.5, netamx*nphimx+ringmx*routmx -0.5);
669 
670  mncrossg = fs->make<TH1F>("mncrossg","mncrossg", netamx*nphimx+ringmx*routmx, -0.5, netamx*nphimx+ringmx*routmx -0.5);
671  rmscrossg = fs->make<TH1F>("rmscrossg","rmscrossg", netamx*nphimx+ringmx*routmx, -0.5, netamx*nphimx+ringmx*routmx -0.5);
672  nevcrossg = fs->make<TH1F>("nevcrossg","nevcrossg", netamx*nphimx+ringmx*routmx, -0.5, netamx*nphimx+ringmx*routmx -0.5);
673 
674 
675  for (int i=0; i<neffip; i++) {
676  if (i==0) {
677  sprintf(title, "Total projected muon in tower");
678  sprintf(name, "total_evt");
679  } else {
680  sprintf(title, "Efficiency with sig >%i #sigma", i);
681  sprintf(name, "Effi_with_gt%i_sig", i);
682  }
683  sig_effi[i] = fs->make<TH2F>(name, title, netamx+1, -netamx/2-0.5, netamx/2+0.5, nphimx, 0.5, nphimx+0.5);
684  }
685 
686  sprintf(title, "Mean Energy of all towers");
687  sprintf(name, "mean_energy");
688  mean_energy = fs->make<TH2F>(name, title, netamx+1, -netamx/2-0.5, netamx/2+0.5, nphimx, 0.5, nphimx+0.5);
689 
690  if (m_correl) { //#ifdef CORREL
691  mncorrsglb = fs->make<TH1F>("mncorrsglb","mncorrsglb", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
692  rmscorrsglb = fs->make<TH1F>("rmscorrsglb","rmscorrsglb", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
693  nevcorrsglb = fs->make<TH1F>("nevcorrsglb","nevcorrsglb", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
694 
695  mncorrsgrb = fs->make<TH1F>("mncorrsgrb","mncorrsgrb", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
696  rmscorrsgrb = fs->make<TH1F>("rmscorrsgrb","rmscorrsgrb", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
697  nevcorrsgrb = fs->make<TH1F>("nevcorrsgrb","nevcorrsgrb", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
698 
699  mncorrsglu = fs->make<TH1F>("mncorrsglu","mncorrsglu", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
700  rmscorrsglu = fs->make<TH1F>("rmscorrsglu","rmscorrsglu", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
701  nevcorrsglu = fs->make<TH1F>("nevcorrsglu","nevcorrsglu", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
702 
703  mncorrsgru = fs->make<TH1F>("mncorrsgru","mncorrsgru", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
704  rmscorrsgru = fs->make<TH1F>("rmscorrsgru","rmscorrsgru", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
705  nevcorrsgru = fs->make<TH1F>("nevcorrsgru","nevcorrsgru", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
706 
707  mncorrsgall = fs->make<TH1F>("mncorrsgall","mncorrsgall", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
708  rmscorrsgall = fs->make<TH1F>("rmscorrsgall","rmscorrsgall", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
709  nevcorrsgall = fs->make<TH1F>("nevcorrsgall","nevcorrsgall", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
710 
711  mncorrsgl = fs->make<TH1F>("mncorrsgl","mncorrsgl", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
712  rmscorrsgl = fs->make<TH1F>("rmscorrsgl","rmscorrsgl", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
713  nevcorrsgl = fs->make<TH1F>("nevcorrsgl","nevcorrsgl", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
714 
715  mncorrsgr = fs->make<TH1F>("mncorrsgr","mncorrsgr", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
716  rmscorrsgr = fs->make<TH1F>("rmscorrsgr","rmscorrsgr", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
717  nevcorrsgr = fs->make<TH1F>("nevcorrsgr","nevcorrsgr", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
718  } //m_correl #endif
719 
720  if (m_checkmap) { //#ifdef CHECKMAP
721  mncorrsgc = fs->make<TH1F>("mncorrsgc","mncorrsgc", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
722  rmscorrsgc = fs->make<TH1F>("rmscorrsgc","rmscorrsgc", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
723  nevcorrsgc = fs->make<TH1F>("nevcorrsgc","nevcorrsgc", netamx*nphimx+60, -0.5, netamx*nphimx+59.5);
724  } //m_checkmap #endif
725 
726  if (m_combined) { //#ifdef COMBINED
727  for (int j=0; j<ringmx; j++) {
728 
729  for (int i=0;i<routmx+1;i++) {
730  if (j!=2 && i>rout12mx) continue;
731  int phmn = 3*i-1;
732  int phmx = 3*i+1;
733  if (j==2) {phmn = 2*i-1; phmx=2*i;}
734  if (phmn <=0) phmn = nphimx+phmn;
735  if (phmx <=0) phmx = nphimx+phmx;
736 
737  if ((j==2 && i==routmx) || (j!=2 && i==rout12mx)) {
738  sprintf(title, "sig_ring%i_allrm", j-2);
739  sprintf(name, "sig_ring%i_allrm", j-2);
740  } else {
741  sprintf(title, "sig_ring%i_phi%i-%i", j-2,phmn,phmx);
742  sprintf(name, "sig_ring%i_rout%i", j-2,i+1);
743  }
744  com_sigrsg[j][i] = fs->make<TH1F>(name, title, nbin, alow, ahigh);
745  if ((j==2 && i==routmx) || (j!=2 && i==rout12mx)) {
746  sprintf(title, "ped_ring%i_allrm", j-2);
747  sprintf(name, "ped_ring%i_allrm", j-2);
748  } else {
749  sprintf(title, "ped_ring%i_phi%i-%i", j-2,phmn, phmx);
750  sprintf(name, "ped_ring%i_rout%i", j-2,i+1);
751  }
752  com_crossg[j][i] = fs->make<TH1F>(name, title, nbin, alow, ahigh);
753  }
754 
755  for (int i=0;i<sectmx;i++) {
756  if (m_hotime) { //#ifdef HOTIME
757  sprintf(title, "com_hotime_ring%i_sect%i", j-2, i+1);
758  com_hotime[j][i] = fs->make<TProfile>(title, title, 10, -0.5, 9.5, -1.0, 30.0);
759 
760  sprintf(title, "com_hopedtime_ring%i_sect%i", j-2, i+1);
761  com_hopedtime[j][i] = fs->make<TProfile>(title, title, 10, -0.5, 9.5, -1.0, 30.0);
762  } //m_hotime #endif
763  if (m_hbtime){ //#ifdef HBTIME
764  sprintf(title, "_com_hbtime_ring%i_serrct%i", j-2, i+1);
765  com_hbtime[j][i] = fs->make<TProfile>(title, title, 10, -0.5, 9.5, -1.0, 30.0);
766  } //m_hbtime #endif
767 
768  if (m_correl) { //#ifdef CORREL
769  sprintf(title, "com_corrsg_ring%i_sect%i_leftbottom", j-2,i+1);
770  com_corrsglb[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
771 
772  sprintf(title, "com_corrsg_ring%i_sect%i_rightbottom", j-2,i+1);
773  com_corrsgrb[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
774 
775  sprintf(title, "com_corrsg_ring%i_sect%i_leftup", j-2,i+1);
776  com_corrsglu[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
777 
778  sprintf(title, "com_corrsg_ring%i_sect%i_rightup", j-2,i+1);
779  com_corrsgru[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
780 
781  sprintf(title, "com_corrsg_ring%i_sect%i_all", j-2,i+1);
782  com_corrsgall[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
783 
784  sprintf(title, "com_corrsg_ring%i_sect%i_left", j-2,i+1);
785  com_corrsgl[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
786 
787  sprintf(title, "com_corrsg_ring%i_sect%i_right", j-2,i+1);
788  com_corrsgr[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
789  } //m_correl #endif
790 
791  if (m_checkmap) { // #ifdef CHECKMAP
792  sprintf(title, "com_corrsg_ring%i_sect%i_centrl", j-2,i+1);
793  com_corrsgc[j][i] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
794  } //m_checkmap #endif
795  }
796  }
797  } //m_combined #endif
798 
799  for (int i=-1; i<=1; i++) {
800  for (int j=-1; j<=1; j++) {
801  int k = 3*(i+1)+j+1;
802 
803  sprintf(title, "hosct2p_eta%i_phi%i", i, j);
804  ho_sig2p[k] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
805 
806  sprintf(title, "hosct1p_eta%i_phi%i", i, j);
807  ho_sig1p[k] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
808 
809  sprintf(title, "hosct00_eta%i_phi%i", i, j);
810  ho_sig00[k] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
811 
812  sprintf(title, "hosct1m_eta%i_phi%i", i, j);
813  ho_sig1m[k] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
814 
815  sprintf(title, "hosct2m_eta%i_phi%i", i, j);
816  ho_sig2m[k] = fs->make<TH1F>(title, title, nbin, alow, ahigh);
817 
818  if (m_hbinfo) { // #ifdef HBINFO
819  sprintf(title, "hbhesig_eta%i_phi%i", i, j);
820  hbhe_sig[k] = fs->make<TH1F>(title, title, 51, -10.5, 40.5);
821  } //m_hbinfo #endif
822  }
823  }
824 
825  if (m_constant) {
826  ped_evt = fs->make<TH1F>("ped_evt", "ped_evt", netamx*nphimx, -0.5, netamx*nphimx-0.5);
827  ped_mean = fs->make<TH1F>("ped_mean", "ped_mean", netamx*nphimx, -0.5, netamx*nphimx-0.5);
828  ped_width = fs->make<TH1F>("ped_width", "ped_width", netamx*nphimx, -0.5, netamx*nphimx-0.5);
829 
830  fit_chi = fs->make<TH1F>("fit_chi", "fit_chi", netamx*nphimx, -0.5, netamx*nphimx-0.5);
831  sig_evt = fs->make<TH1F>("sig_evt", "sig_evt", netamx*nphimx, -0.5, netamx*nphimx-0.5);
832  fit_sigevt = fs->make<TH1F>("fit_sigevt", "fit_sigevt", netamx*nphimx, -0.5, netamx*nphimx-0.5);
833  fit_bkgevt = fs->make<TH1F>("fit_bkgevt", "fit_bkgevt", netamx*nphimx, -0.5, netamx*nphimx-0.5);
834  sig_mean = fs->make<TH1F>("sig_mean", "sig_mean", netamx*nphimx, -0.5, netamx*nphimx-0.5);
835  sig_diff = fs->make<TH1F>("sig_diff", "sig_diff", netamx*nphimx, -0.5, netamx*nphimx-0.5);
836  sig_width = fs->make<TH1F>("sig_width", "sig_width", netamx*nphimx, -0.5, netamx*nphimx-0.5);
837  sig_sigma = fs->make<TH1F>("sig_sigma", "sig_sigma", netamx*nphimx, -0.5, netamx*nphimx-0.5);
838  sig_meanerr = fs->make<TH1F>("sig_meanerr", "sig_meanerr", netamx*nphimx, -0.5, netamx*nphimx-0.5);
839  sig_meanerrp = fs->make<TH1F>("sig_meanerrp", "sig_meanerrp", netamx*nphimx, -0.5, netamx*nphimx-0.5);
840  sig_signf = fs->make<TH1F>("sig_signf", "sig_signf", netamx*nphimx, -0.5, netamx*nphimx-0.5);
841 
842  ped_statmean = fs->make<TH1F>("ped_statmean", "ped_statmean", netamx*nphimx, -0.5, netamx*nphimx-0.5);
843  sig_statmean = fs->make<TH1F>("sig_statmean", "sig_statmean", netamx*nphimx, -0.5, netamx*nphimx-0.5);
844  ped_rms = fs->make<TH1F>("ped_rms", "ped_rms", netamx*nphimx, -0.5, netamx*nphimx-0.5);
845  sig_rms = fs->make<TH1F>("sig_rms", "sig_rms", netamx*nphimx, -0.5, netamx*nphimx-0.5);
846 
847  const_eta_phi = fs->make<TH2F>("const_eta_phi", "const_eta_phi", netamx+1, -(netamx+1)/2., (netamx+1)/2., nphimx, 0.5, nphimx+0.5);
848 
849  for (int ij=0; ij<netamx; ij++) {
850  int ieta = (ij<15) ? ij+1 : 14-ij;
851  sprintf(title, "Cont_Eta_%i", ieta);
852  const_eta[ij] = fs->make<TH1F>(title, title, nphimx, 0.5, nphimx+0.5);
853 
854  sprintf(title, "Peak_Eta_%i", ieta);
855  peak_eta[ij] = fs->make<TH1F>(title, title, nphimx, 0.5, nphimx+0.5);
856  }
857 
858  for (int ij=0; ij<ringmx; ij++) {
859  int iring = ij-2;
860  int iread = (ij==2) ? routmx : rout12mx;
861  sprintf(title, "Cont_hpdrm_%i", iring);
862  const_hpdrm[ij] = fs->make<TH1F>(title, title, iread, 0.5, iread+0.5);
863 
864  sprintf(title, "Peak_hpdrm_%i", iring);
865  peak_hpdrm[ij] = fs->make<TH1F>(title, title, iread, 0.5, iread+0.5);
866  }
867 
868  mean_phi_hst = fs->make<TH1F>("mean_phi_hst", "mean_phi_hst", netamx+1, -(netamx+1)/2., (netamx+1)/2.);
869  mean_phi_ave = fs->make<TH1F>("mean_phi_ave", "mean_phi_ave", netamx+1, -(netamx+1)/2., (netamx+1)/2.);
870 
871  mean_eta_ave = fs->make<TH1F>("mean_eta_ave", "mean_eta_ave", nphimx, 0.5, nphimx+0.5);
872 
873  } //m_constant
874 
875  for (int ij=0; ij<netamx; ij++) {
876  int ieta = (ij<15) ? ij+1 : 14-ij;
877 
878  sprintf(title, "Stat_Eta_%i", ieta);
879  stat_eta[ij] = fs->make<TH1F>(title, title, nphimx, 0.5, nphimx+0.5);
880 
881  sprintf(title, "#mu(stat)_Eta_%i", ieta);
882  statmn_eta[ij] = fs->make<TH1F>(title, title, nphimx, 0.5, nphimx+0.5);
883  }
884 
885  for (int j=0; j<netamx; j++) {
886  for (int i=0; i<nphimx; i++) {
887  invang[j][i] = 0.0;
888  }
889  }
890  for (int j=0; j<ringmx; j++) {
891  for (int i=0;i<routmx+1;i++) {
892  com_invang[j][i] = 0.0;
893  }
894  }
895 
896 }
897 
898 
900 {
901 
902  // do anything here that needs to be done at desctruction time
903  // (e.g. close files, deallocate resources etc.)
904 
905  theFile->cd();
906  theFile->Write();
907  theFile->Close();
908  cout <<" Selected events # is "<<ipass<<endl;
909 }
910 
911 
912 //
913 // member functions
914 //
915 
916 // ------------ method called to for each event ------------
917 void
919 {
920 
921  // calcualte these once (and avoid the pow(int,int) ambiguities for c++)
922  int mypow_2_0 = 0; // 2^0
923  int mypow_2_1 = 2; // 2^1
924  int mypow_2_2 = 4; // 2^2
925 
926  int mypow_2_3 = 8; // 2^3
927  int mypow_2_4 = 16; // 2^4
928  int mypow_2_5 = 32; // 2^5
929  int mypow_2_6 = 64; // 2^6
930  int mypow_2_7 = 128; // 2^7
931  int mypow_2_8 = 256; // 2^8
932  int mypow_2_9 = 512; // 2^9
933  int mypow_2_10 = 1024; // 2^10
934  int mypow_2_11 = 2048; // 2^11
935  int mypow_2_12 = 4096; // 2^12
936 
937  /*
938  //FIXGM Put this is initialiser
939  int mapx1[6][3]={{1,4,8}, {12,7,3}, {5,9,13}, {11,6,2}, {16,15,14}, {19,18,17}};
940  // int etamap1[21]={-1, 0,3,1, 0,2,3, 1,0,2, -1, 3,1,2, 4,4,4, 5,5,5, -1};
941  // int phimap1[21]={-1, 0,2,2, 1,0,1, 1,2,1, -1, 0,0,2, 2,1,0, 2,1,0,-1};
942 
943  int mapx2[6][3]={{1,4,8}, {12,7,3}, {5,9,13}, {11,6,2}, {16,15,14}, {-1,-1,-1}};
944  // int etamap2[21]={-1, 0,3,1, 0,2,3, 1,0,2, -1, 3,1,2, 4,4,4, -1,-1,-1, -1};
945  // int phimap2[21]={-1, 0,2,2, 1,0,1, 1,2,1, -1, 0,0,2, 2,1,0, 2, 1, 0, -1};
946 
947  int mapx0p[9][2]={{3,1}, {7,4}, {6,5}, {12,8}, {0,0}, {11,9}, {16,13}, {15,14}, {19,17}};
948  int mapx0m[9][2]={{17,19}, {14,15}, {13,16}, {9,11}, {0,0}, {8,12}, {5,6}, {4,7}, {1,3}};
949 
950  // int etamap0p[21]={-1, 0,-1,0, 1,2,2, 1,3,5, -1, 5,3,6, 7,7,6, 8,-1,8, -1};
951  // int phimap0p[21]={-1, 1,-1,0, 1,1,0, 0,1,1, -1, 0,0,1, 1,0,0, 1,-1,0, -1};
952 
953  // int etamap0m[21]={-1, 8,-1,8, 7,6,6, 7,5,3, -1, 3,5,2, 1,1,2, 0,-1,0, -1};
954  // int phimap0m[21]={-1, 0,-1,1, 0,0,1, 1,0,0, -1, 1,1,0, 0,1,1, 0,-1,1, -1};
955 
956  int etamap[4][21]={{-1, 0,3,1, 0,2,3, 1,0,2, -1, 3,1,2, 4,4,4, -1,-1,-1, -1}, //etamap2
957  {-1, 0,3,1, 0,2,3, 1,0,2, -1, 3,1,2, 4,4,4, 5,5,5, -1}, //etamap1
958  {-1, 0,-1,0, 1,2,2, 1,3,5, -1, 5,3,6, 7,7,6, 8,-1,8, -1}, //etamap0p
959  {-1, 8,-1,8, 7,6,6, 7,5,3, -1, 3,5,2, 1,1,2, 0,-1,0, -1}}; //etamap0m
960 
961  int phimap[4][21] ={{-1, 0,2,2, 1,0,1, 1,2,1, -1, 0,0,2, 2,1,0, 2,1,0, -1}, //phimap2
962  {-1, 0,2,2, 1,0,1, 1,2,1, -1, 0,0,2, 2,1,0, 2,1,0, -1}, //phimap1
963  {-1, 1,-1,0, 1,1,0, 0,1,1, -1, 0,0,1, 1,0,0, 1,-1,0, -1}, //phimap0p
964  {-1, 0,-1,1, 0,0,1, 1,0,0, -1, 1,1,0, 0,1,1, 0,-1,1, -1}}; //phimap0m
965  //swapped phi map for R0+/R0- (15/03/07)
966  for (int i=0; i<4; i++) {
967  for (int j=0; j<21; j++) {
968  cout <<"ieta "<<i<<" "<<j<<" "<<etamap[i][j]<<endl;
969  }
970  }
971 
972  // Character convention for R+/-1/2
973  // int npixleft[21] = {-1, F, Q,-1, M, D, J,-1, T,-1, C,-1, R, P, H,-1, N, G};
974  // int npixrigh[21] = { Q, S,-1, D, J, L,-1, K,-1, E,-1, P, H, B,-1, G, A,-1};
975 
976  // int npixlb1[21]={-1,-1,-1,-1, F, Q, S,-1, M, J, L, T, K,-1, C, R, P, H};
977  // int npixrb1[21]={-1,-1,-1, F, Q, S,-1, M, D, L,-1, K,-1, C, E, P, H, B};
978  // int npixlu1[21]={ M, D, J, T, K,-1, C,-1, R, H, B,-1, N, G, A,-1,-1,-1};
979  // int npixru1[21]={ D, J, L, K,-1, C, E, R, P, B,-1, N, G, A,-1,-1,-1,-1};
980 
981  int npixleft[21]={0, 0, 1, 2, 0, 4, 5, 6, 0, 8, 0, 0,11, 0,13,14,15, 0,17,18,0};
982  int npixrigh[21]={0, 2, 3, 0, 5, 6, 7, 0, 9, 0, 0,12, 0,14,15,16, 0,18,19, 0,0};
983  int npixlebt[21]={0, 0, 0, 0, 0, 1, 2, 3, 0, 4, 0, 6, 7, 8, 9, 0,11,13,14,15,0};
984  int npixribt[21]={0, 0, 0, 0, 1, 2, 3, 0, 4, 5, 0, 7, 0, 9, 0,11,12,14,15,16,0};
985  int npixleup[21]={0, 4, 5, 6, 8, 9, 0,11, 0,13, 0,15,16, 0,17,18,19, 0, 0, 0,0};
986  int npixriup[21]={0, 5, 6, 7, 9, 0,11,12,13,14, 0,16, 0,17,18,19, 0, 0, 0, 0,0};
987  */
988 
989  int iaxxx = 0;
990  int ibxxx = 0;
991 
992  Nevents++;
993 
994  using namespace edm;
995 
996  float pival = acos(-1.);
997  irunold = irun = iEvent.id().run();
998  ievt = iEvent.id().event();
999 
1000  // cout <<"Nevents "<<Nevents<<endl;
1001 
1003  bool isCosMu = true;
1004  try {
1005  iEvent.getByLabel(hoCalibVariableCollectionTag, HOCalib);
1006  // iEvent.getByLabel("hoCalibProducer","HOCalibVariableCollection",HOCalib);
1007 
1008  } catch ( cms::Exception &iEvent ) { isCosMu = false; }
1009 
1010  if (isCosMu && (*HOCalib).size() >0 ) {
1011  nmuon = (*HOCalib).size();
1012  if (Nevents%5000==1) cout <<"nmuon event # "<<Nevents<<" Run # "<<iEvent.id().run()<<" Evt # "<<iEvent.id().event()<<" "<<nmuon<<endl;
1013 
1014  for (HOCalibVariableCollection::const_iterator hoC=(*HOCalib).begin(); hoC!=(*HOCalib).end(); hoC++){
1015  itrg1 = (*hoC).trig1;
1016  itrg2 = (*hoC).trig2;
1017 
1018  trkdr = (*hoC).trkdr;
1019  trkdz = (*hoC).trkdz;
1020 
1021  trkvx = (*hoC).trkvx;
1022  trkvy = (*hoC).trkvy;
1023  trkvz = (*hoC).trkvz;
1024 
1025  trkmm = (*hoC).trkmm;
1026  trkth = (*hoC).trkth;
1027  trkph = (*hoC).trkph;
1028 
1029  ndof = (int)(*hoC).ndof;
1030  // nrecht = (int)(*hoC).nrecht;
1031  chisq = (*hoC).chisq;
1032 
1033  therr = (*hoC).therr;
1034  pherr = (*hoC).pherr;
1035  trkph = (*hoC).trkph;
1036 
1037  isect = (*hoC).isect;
1038  hodx = (*hoC).hodx;
1039  hody = (*hoC).hody;
1040  hoang = (*hoC).hoang;
1041  htime = (*hoC).htime;
1042  for (int i=0; i<9; i++) { hosig[i] = (*hoC).hosig[i];} //cout<<"hosig "<<i<<" "<<hosig[i]<<endl;}
1043  for (int i=0; i<18; i++) { hocorsig[i] = (*hoC).hocorsig[i];} // cout<<"hocorsig "<<i<<" "<<hocorsig[i]<<endl;}
1044  hocro = (*hoC).hocro;
1045  for (int i=0; i<3; i++) { caloen[i] = (*hoC).caloen[i];}
1046 
1047  if (m_hbinfo) { for (int i=0; i<9; i++) { hbhesig[i] = (*hoC).hbhesig[i];}} // cout<<"hbhesig "<<i<<" "<<hbhesig[i]<<endl;}}
1048 
1049 
1050  int ipsall=0;
1051  int ips0=0;
1052  int ips1=0;
1053  int ips2=0;
1054  int ips3=0;
1055  int ips4=0;
1056  int ips5=0;
1057  int ips6=0;
1058  int ips7=0;
1059  int ips8=0;
1060  int ips9=0;
1061  int ips10 =0;
1062  int ips11 =0;
1063  int ips12 = 0;
1064 
1065  int iselect3 = 0;
1066  if (ndof >=15 && chisq <30) iselect3 = 1;
1067 
1068  if (isect <0) continue; //FIXGM Is it proper place ?
1069  if (fabs(trkth-pival/2)<0.000001) continue; //22OCT07
1070 
1071  int ieta = int((abs(isect)%10000)/100.)-30; //an offset to acodate -ve eta values
1072  if (abs(ieta)>=16) continue;
1073  int iphi = abs(isect)%100;
1074 
1075  int tmpsect = int((iphi + 1)/6.) + 1;
1076  if (tmpsect >12) tmpsect = 1;
1077 
1078  int iring = 0;
1079  int tmpeta = ieta + 4; //For pixel mapping
1080  if (ieta >=-15 && ieta <=-11) {iring = -2; tmpeta =-11-ieta; } //abs(ieta)-11;}
1081  if (ieta >=-10 && ieta <=-5) {iring = -1; tmpeta =-5-ieta; } // abs(ieta)-5;}
1082  if (ieta >= 5 && ieta <= 10) {iring = 1; tmpeta =ieta-5; }
1083  if (ieta >= 11 && ieta <= 15) {iring = 2; tmpeta =ieta-11;}
1084 
1085  int iring2 = iring + 2;
1086 
1087  int tmprout = int((iphi + 1)/3.) + 1;
1088  int tmproutmx =routmx;
1089  if (iring==0) {
1090  tmprout = int((iphi + 1)/2.) + 1;
1091  if (tmprout >routmx) tmprout = 1;
1092  } else {
1093  if (tmprout >rout12mx) tmprout = 1;
1094  tmproutmx =rout12mx;
1095  }
1096 
1097  // CRUZET1
1098  if (m_cosmic) {
1099  /* GMA temoparily change to increase event size at 3 & 9 O'clock position */
1100  if (abs(ndof) >=20 && abs(ndof) <40) {ips0 = (int)mypow_2_0; ipsall += ips0;}
1101  if (chisq >0 && chisq<15) {ips1 = (int)mypow_2_1; ipsall +=ips1;} //18Jan2008
1102  if (fabs(trkth-pival/2) <21.5) {ips2 = (int)mypow_2_2; ipsall +=ips2;} //No nead for pp evt
1103  if (fabs(trkph+pival/2) <21.5) {ips3 = (int)mypow_2_3; ipsall +=ips3;} //No nead for pp evt
1104 
1105  if (therr <0.02) {ips4 = (int)mypow_2_4; ipsall +=ips4;}
1106  if (pherr <0.0002) {ips5 = (int)mypow_2_5; ipsall +=ips5;}
1107  if (fabs(hoang) >0.30) {ips6 = (int)mypow_2_6; ipsall +=ips6;}
1108  if (fabs(trkmm) >0.100) {ips7 = (int)mypow_2_7; ipsall +=ips7;}
1109  // if (nmuon ==1) {ips8 = (int)mypow_2_8; ipsall +=ips8;}
1110  if (nmuon >=1 && nmuon <=4) {ips8 = (int)mypow_2_8; ipsall +=ips8;}
1111 
1112  if (iring2==2) {
1113  if (fabs(hodx)<100 && fabs(hodx)>2 && fabs(hocorsig[8]) <40 && fabs(hocorsig[8]) >2 )
1114  {ips10=(int)mypow_2_10;ipsall +=ips10;}
1115 
1116  if (fabs(hody)<100 && fabs(hody)>2 && fabs(hocorsig[9]) <40 && fabs(hocorsig[9]) >2 )
1117  {ips11=(int)mypow_2_11;ipsall +=ips11;}
1118 
1119  } else {
1120  if (fabs(hodx)<100 && fabs(hodx)>2 )
1121  {ips10=(int)mypow_2_10;ipsall +=ips10;}
1122 
1123  if (fabs(hody)<100 && fabs(hody)>2)
1124  {ips11=(int)mypow_2_11;ipsall +=ips11;}
1125  }
1126  if (caloen[0] ==0) { ips12=(int)mypow_2_12;ipsall +=ips12;}
1127  } else {
1128  //csa08
1129  if (abs(ndof) >=20 && abs(ndof) <40) {ips0 = (int)mypow_2_0; ipsall += ips0;}
1130  if (chisq >0 && chisq<15) {ips1 = (int)mypow_2_1; ipsall +=ips1;} //18Jan2008
1131  if (fabs(trkth-pival/2) <21.5) {ips2 = (int)mypow_2_2; ipsall +=ips2;} //No nead for pp evt
1132  if (fabs(trkph+pival/2) <21.5) {ips3 = (int)mypow_2_3; ipsall +=ips3;} //No nead for pp evt
1133 
1134  if (therr <0.02) {ips4 = (int)mypow_2_4; ipsall +=ips4;}
1135  if (pherr <0.0002) {ips5 = (int)mypow_2_5; ipsall +=ips5;}
1136  if (fabs(hoang) >0.30) {ips6 = (int)mypow_2_6; ipsall +=ips6;}
1137  if (fabs(trkmm) >4.0) {ips7 = (int)mypow_2_7; ipsall +=ips7;}
1138  if (nmuon >=1 && nmuon <=2) {ips8 = (int)mypow_2_8; ipsall +=ips8;}
1139 
1140  if (iring2==2) {
1141  if (fabs(hodx)<100 && fabs(hodx)>2 && fabs(hocorsig[8]) <40 && fabs(hocorsig[8]) >2 )
1142  {ips10=(int)mypow_2_10;ipsall +=ips10;}
1143 
1144  if (fabs(hody)<100 && fabs(hody)>2 && fabs(hocorsig[9]) <40 && fabs(hocorsig[9]) >2 )
1145  {ips11=(int)mypow_2_11;ipsall +=ips11;}
1146 
1147  } else {
1148  if (fabs(hodx)<100 && fabs(hodx)>2 )
1149  {ips10=(int)mypow_2_10;ipsall +=ips10;}
1150 
1151  if (fabs(hody)<100 && fabs(hody)>2)
1152  {ips11=(int)mypow_2_11;ipsall +=ips11;}
1153  }
1154  // if (m_cosmic || (caloen[0] >0.5 && caloen[0]<5.0)) {ips12=(int)pow_2_12;ipsall +=ips12;}
1155  if (ndof >0 && caloen[0]<5.0) {ips12=(int)mypow_2_12;ipsall +=ips12;}
1156  /* */
1157  }
1158 
1159  if (m_digiInput) {
1160  if (htime >2.5 && htime <6.5) {ips9=(int)mypow_2_9;ipsall +=ips9;}
1161  } else {
1162  if (htime >-40 && htime <60) {ips9=(int)mypow_2_9;ipsall +=ips9;}
1163  }
1164 
1165  if (ipsall-ips0==mypow_2_ncut-mypow_2_0-1) sigvsevt[iring2][0]->Fill(abs(ndof), hosig[4]);
1166  if (ipsall-ips1==mypow_2_ncut-mypow_2_1-1) sigvsevt[iring2][1]->Fill(chisq, hosig[4]);
1167  if (ipsall-ips2==mypow_2_ncut-mypow_2_2-1) sigvsevt[iring2][2]->Fill(trkth, hosig[4]);
1168  if (ipsall-ips3==mypow_2_ncut-mypow_2_3-1) sigvsevt[iring2][3]->Fill(trkph, hosig[4]);
1169  if (ipsall-ips4==mypow_2_ncut-mypow_2_4-1) sigvsevt[iring2][4]->Fill(therr, hosig[4]);
1170  if (ipsall-ips5==mypow_2_ncut-mypow_2_5-1) sigvsevt[iring2][5]->Fill(pherr, hosig[4]);
1171  if (ipsall-ips6==mypow_2_ncut-mypow_2_6-1) sigvsevt[iring2][6]->Fill(hoang, hosig[4]);
1172  if (ipsall-ips7==mypow_2_ncut-mypow_2_7-1) sigvsevt[iring2][7]->Fill(fabs(trkmm), hosig[4]);
1173  if (ipsall-ips8==mypow_2_ncut-mypow_2_8-1) sigvsevt[iring2][8]->Fill(nmuon, hosig[4]);
1174  if (ipsall-ips9==mypow_2_ncut-mypow_2_9-1) sigvsevt[iring2][9]->Fill(htime, hosig[4]);
1175  if (ipsall-ips10==mypow_2_ncut-mypow_2_10-1) sigvsevt[iring2][10]->Fill(hodx, hosig[4]);
1176  if (ipsall-ips11==mypow_2_ncut-mypow_2_11-1) sigvsevt[iring2][11]->Fill(hody, hosig[4]);
1177  if (!m_cosmic) {
1178  if (ipsall-ips12==mypow_2_ncut-mypow_2_12-1) sigvsevt[iring2][12]->Fill(caloen[0], hosig[4]);
1179  }
1180 
1181  sigvsevt[iring2+5][0]->Fill(abs(ndof), hosig[4]);
1182  if (ips0 >0) {
1183  sigvsevt[iring2+5][1]->Fill(chisq, hosig[4]);
1184  if (ips1 >0) {
1185  sigvsevt[iring2+5][2]->Fill(trkth, hosig[4]);
1186  if (ips2 >0) {
1187  sigvsevt[iring2+5][3]->Fill(trkph, hosig[4]);
1188  if (ips3 >0) {
1189  sigvsevt[iring2+5][4]->Fill(therr, hosig[4]);
1190  if (ips4 >0) {
1191  sigvsevt[iring2+5][5]->Fill(pherr, hosig[4]);
1192  if (ips5 >0) {
1193  sigvsevt[iring2+5][6]->Fill(hoang, hosig[4]);
1194  if (ips6 >0) {
1195  sigvsevt[iring2+5][7]->Fill(fabs(trkmm), hosig[4]);
1196  if (ips7 >0) {
1197  sigvsevt[iring2+5][8]->Fill(nmuon, hosig[4]);
1198  if (ips8 >0) {
1199  sigvsevt[iring2+5][9]->Fill(htime, hosig[4]);
1200  if (ips9 >0) {
1201  sigvsevt[iring2+5][10]->Fill(hodx, hosig[4]);
1202  if (ips10>0) {
1203  sigvsevt[iring2+5][11]->Fill(hody, hosig[4]);
1204  if (ips11>0) {
1205  if (!m_cosmic) sigvsevt[iring2+5][12]->Fill(caloen[0], hosig[4]);
1206  }
1207  }
1208  }
1209  }
1210  }
1211  }
1212  }
1213  }
1214  }
1215  }
1216  }
1217  }
1218 
1219  sigvsevt[iring2+10][0]->Fill(abs(ndof), hosig[4]);
1220  sigvsevt[iring2+10][1]->Fill(chisq, hosig[4]);
1221  sigvsevt[iring2+10][2]->Fill(trkth, hosig[4]);
1222  sigvsevt[iring2+10][3]->Fill(trkph, hosig[4]);
1223  sigvsevt[iring2+10][4]->Fill(therr, hosig[4]);
1224  sigvsevt[iring2+10][5]->Fill(pherr, hosig[4]);
1225  sigvsevt[iring2+10][6]->Fill(hoang, hosig[4]);
1226  sigvsevt[iring2+10][7]->Fill(fabs(trkmm), hosig[4]);
1227  sigvsevt[iring2+10][8]->Fill(nmuon, hosig[4]);
1228  sigvsevt[iring2+10][9]->Fill(htime, hosig[4]);
1229  sigvsevt[iring2+10][10]->Fill(hodx, hosig[4]);
1230  sigvsevt[iring2+10][11]->Fill(hody, hosig[4]);
1231  if (!m_cosmic) sigvsevt[iring2+10][12]->Fill(caloen[0], hosig[4]);
1232 
1233  int iselect = (ipsall == mypow_2_ncut - 1) ? 1 : 0;
1234 
1235  if (hocro !=-100.0 && hocro <-50.0) hocro +=100.;
1236 
1237  muonnm->Fill(nmuon);
1238  muonmm->Fill(trkmm);
1239  muonth->Fill(trkth*180/pival);
1240  muonph->Fill(trkph*180/pival);
1241  muonch->Fill(chisq);
1242 
1243  if (iselect==1) {
1244  ipass++;
1245  sel_muonnm->Fill(nmuon);
1246  sel_muonmm->Fill(trkmm);
1247  sel_muonth->Fill(trkth*180/pival);
1248  sel_muonph->Fill(trkph*180/pival);
1249  sel_muonch->Fill(chisq);
1250  }
1251 
1252  if (iselect3) T1->Fill();
1253 
1254  int tmpphi = (iphi + 1)%3; //pixel mapping
1255  int npixel = 0;
1256  int itag = -1;
1257  int iflip = 0;
1258  int fact = 2;
1259 
1260  if (iring ==0) {
1261  tmpphi = (iphi+1)%2;
1262  if (tmpsect==2 || tmpsect==3 || tmpsect==6 || tmpsect==7 || tmpsect==10 || tmpsect==11) {
1263  npixel = mapx0p[tmpeta][tmpphi];
1264  itag = 2;
1265  } else {
1266  npixel = mapx0m[tmpeta][tmpphi];
1267  itag = 3;
1268  }
1269  } else {
1270  fact = 3;
1271  if (tmpsect%2==1) iflip =1;
1272  if (abs(iring)==1) {
1273  npixel = mapx1[tmpeta][(iflip==0) ? tmpphi : abs(tmpphi-2)];
1274  itag = 1;
1275  } else {
1276  npixel = mapx2[tmpeta][(iflip==0) ? tmpphi : abs(tmpphi-2)];
1277  itag = 0;
1278  }
1279  }
1280 
1281  int tmpeta1 = (ieta>0) ? ieta -1 : -ieta +14;
1282 
1283  //Histogram filling for noise study: phi shift according to DTChamberAnalysis
1284  int tmpphi1 = iphi -1 ; //(iphi+6 <=nphimx) ? iphi+5 : iphi+5-nphimx;
1285 
1286  int iselect2=0;
1287  if (hosig[4]!=-100) {
1288  if (m_cosmic) {
1289  if (caloen[2]<=0.0) iselect2=1;
1290  } else {
1291  if (caloen[2]<=3.0) iselect2=1;
1292  }
1293  }
1294 
1295  // cout <<"cosmic "<<hosig[4]<<" "<<caloen[3]<<" "<<int(iselect2)<<" "<<int(m_cosmic)<<endl;
1296 
1297 
1298  if (iselect2==1) {
1299  int tmpphi2 = iphi - 1;
1300  if (!m_digiInput) { tmpphi2 = (iphi+6 <=nphimx) ? iphi+5 : iphi+5-nphimx;}
1301 
1302  int tmpsect2 = int((tmpphi2 + 2)/6.) + 1;
1303  if (tmpsect2 >12) tmpsect2 = 1;
1304 
1305  int tmprout2 = int((tmpphi2 + 2)/3.) + 1;
1306  if (iring==0) {
1307  tmprout2 = int((tmpphi2 + 2)/2.) + 1;
1308  if (tmprout2 >routmx) tmprout2 = 1;
1309  } else {
1310  if (tmprout2 >rout12mx) tmprout2 = 1;
1311  }
1312 
1313  if (cro_ssg[tmpeta1][tmpphi2].size()<4000) {
1314  if (hocro>alow && hocro<ahigh) {
1315  if (!m_histfit) cro_ssg[tmpeta1][tmpphi2].push_back(hocro);
1316  crossg[tmpeta1][tmpphi2]->Fill(hocro);
1317  }
1318  }
1319 
1320  if (tmpphi2 >=0 && tmpphi2 <nphimx) {
1321  crossg[tmpeta1][nphimx]->Fill(hocro);
1322  }
1323  if (m_combined) {
1324  com_crossg[iring2][tmprout2-1]->Fill(hocro);
1325  com_crossg[iring2][tmproutmx]->Fill(hocro);
1326  }
1327  }
1328 
1329  if (iselect==1) {
1330  for (int ij=0; ij<neffip; ij++) {
1331  if (ij==0) {
1332  sig_effi[ij]->Fill(ieta, iphi, 1.);
1333  } else {
1334  if (hosig[4] >ij*m_sigma) {
1335  sig_effi[ij]->Fill(ieta, iphi, 1.);
1336  }
1337  }
1338  }
1339 
1340  tmpphi1 = iphi - 1;
1341 
1342  if (sig_reg[tmpeta1][tmpphi1].size()<4000 ) {
1343  if (hosig[4]>-50&& hosig[4] <15) {
1344  sigrsg[tmpeta1][tmpphi1]->Fill(hosig[4]);
1345  if (!m_histfit && hosig[4]<=ahigh/2.) sig_reg[tmpeta1][tmpphi1].push_back(hosig[4]);
1346  invang[tmpeta1][tmpphi1] += 1./fabs(hoang);
1347  }
1348  }
1349 
1350  if (tmpphi1 >=0 && tmpphi1 <nphimx) { //GREN
1351  sigrsg[tmpeta1][nphimx]->Fill(hosig[4]);
1352  invang[tmpeta1][nphimx] += 1./fabs(hoang);
1353  }
1354 
1355  if (m_combined) { //#ifdef COMBINED
1356  com_sigrsg[iring2][tmprout-1]->Fill(hosig[4]);
1357  com_invang[iring2][tmprout-1] += 1./fabs(hoang);
1358 
1359  com_sigrsg[iring2][tmproutmx]->Fill(hosig[4]);
1360  com_invang[iring2][tmproutmx] += 1./fabs(hoang);
1361  } //m_combined #endif
1362 
1363  if (m_checkmap || m_correl) { //#ifdef CHECKMAP
1364  tmpeta = etamap[itag][npixel];
1365  tmpphi = phimap[itag][npixel];
1366  if (tmpeta>=0 && tmpphi >=0) {
1367  if (iflip !=0) tmpphi = abs(tmpphi-2);
1368  if (int((hocorsig[fact*tmpeta+tmpphi]-hosig[4])*10000)/10000.!=0) {
1369  iaxxx++;
1370  cout<<"iring2xxx "<<irun<<" "<<ievt<<" "<<isect<<" "<<iring<<" "<<tmpsect<<" "<<ieta<<" "<<iphi<<" "<<npixel<<" "<<tmpeta<<" "<<tmpphi<<" "<<tmpeta1<<" "<<tmpphi1<<" itag "<<itag<<" "<<iflip<<" "<<fact<<" "<<hocorsig[fact*tmpeta+tmpphi]<<" "<<fact*tmpeta+tmpphi<<" "<<hosig[4]<<" "<<hodx<<" "<<hody<<endl;
1371 
1372  for (int i=0; i<18; i++) {cout <<" "<<i<<" "<<hocorsig[i];}
1373  cout<<" ix "<<iaxxx<<" "<<ibxxx<<endl;
1374  } else { ibxxx++; }
1375 
1376  corrsgc[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1377  if (m_combined) { //#ifdef COMBINED
1378  com_corrsgc[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1379  } //m_combined #endif
1380  }
1381  } //m_checkmap #endif
1382 
1383  if (m_correl) { //#ifdef CORREL
1384  float allcorsig = 0.0;
1385 
1386  tmpeta = etamap[itag][npixleft[npixel]];
1387  tmpphi = phimap[itag][npixleft[npixel]];
1388 
1389  if (tmpeta>=0 && tmpphi >=0) {
1390  if (iflip !=0) tmpphi = abs(tmpphi-2);
1391  corrsgl[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1392  allcorsig +=hocorsig[fact*tmpeta+tmpphi];
1393  if (m_combined) { //#ifdef COMBINED
1394  com_corrsgl[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1395  } //m_combined #endif
1396  }
1397 
1398  tmpeta = etamap[itag][npixrigh[npixel]];
1399  tmpphi = phimap[itag][npixrigh[npixel]];
1400  if (tmpeta>=0 && tmpphi >=0) {
1401  if (iflip !=0) tmpphi = abs(tmpphi-2);
1402  corrsgr[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1403  allcorsig +=hocorsig[fact*tmpeta+tmpphi];
1404  if (m_combined) { // #ifdef COMBINED
1405  com_corrsgr[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1406  } //m_combined #endif
1407  }
1408 
1409  tmpeta = etamap[itag][npixlebt[npixel]];
1410  tmpphi = phimap[itag][npixlebt[npixel]];
1411  if (tmpeta>=0 && tmpphi >=0) {
1412  if (iflip !=0) tmpphi = abs(tmpphi-2);
1413  corrsglb[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1414  allcorsig +=hocorsig[fact*tmpeta+tmpphi];
1415  if (m_combined){ //#ifdef COMBINED
1416  com_corrsglb[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1417  } //m_combined #endif
1418  }
1419 
1420  tmpeta = etamap[itag][npixribt[npixel]];
1421  tmpphi = phimap[itag][npixribt[npixel]];
1422  if (tmpeta>=0 && tmpphi >=0) {
1423  if (iflip !=0) tmpphi = abs(tmpphi-2);
1424  corrsgrb[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1425  allcorsig +=hocorsig[fact*tmpeta+tmpphi];
1426  if (m_combined) { // #ifdef COMBINED
1427  com_corrsgrb[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1428  } //m_combined #endif
1429  }
1430 
1431  tmpeta = etamap[itag][npixleup[npixel]];
1432  tmpphi = phimap[itag][npixleup[npixel]];
1433  if (tmpeta>=0 && tmpphi >=0) {
1434  if (iflip !=0) tmpphi = abs(tmpphi-2);
1435  corrsglu[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1436  allcorsig +=hocorsig[fact*tmpeta+tmpphi];
1437  if (m_combined) {// #ifdef COMBINED
1438  com_corrsglu[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1439  } //m_combined #endif
1440  }
1441 
1442  tmpeta = etamap[itag][npixriup[npixel]];
1443  tmpphi = phimap[itag][npixriup[npixel]];
1444  if (tmpeta>=0 && tmpphi >=0) {
1445  if (iflip !=0) tmpphi = abs(tmpphi-2);
1446  corrsgru[tmpeta1][tmpphi1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1447  allcorsig +=hocorsig[fact*tmpeta+tmpphi];
1448  if (m_combined) { // #ifdef COMBINED
1449  com_corrsgru[iring2][tmpsect-1]->Fill(hocorsig[fact*tmpeta+tmpphi]);
1450  } //m_combined #endif
1451  }
1452  corrsgall[tmpeta1][tmpphi1]->Fill(allcorsig);
1453  if (m_combined) { // #ifdef COMBINED
1454  com_corrsgall[iring2][tmpsect-1]->Fill(allcorsig);
1455  } //m_combined #endif
1456 
1457 
1458  } //m_correl #endif
1459  for (int k=0; k<9; k++) {
1460  switch (iring)
1461  {
1462  case 2 : ho_sig2p[k]->Fill(hosig[k]); break;
1463  case 1 : ho_sig1p[k]->Fill(hosig[k]); break;
1464  case 0 : ho_sig00[k]->Fill(hosig[k]); break;
1465  case -1 : ho_sig1m[k]->Fill(hosig[k]); break;
1466  case -2 : ho_sig2m[k]->Fill(hosig[k]); break;
1467  }
1468  if (m_hbinfo) { // #ifdef HBINFO
1469  hbhe_sig[k]->Fill(hbhesig[k]);
1470  // cout <<"hbhe "<<k<<" "<<hbhesig[k]<<endl;
1471  } //m_hbinfo #endif
1472  }
1473  } //if (iselect==1)
1474  } //for (HOCalibVariableCollection::const_iterator hoC=(*HOCalib).begin(); hoC!=(*HOCalib).end(); hoC++){
1475 
1476  } //if (isCosMu)
1477 
1478 #ifdef THIS_IS_AN_EVENT_EXAMPLE
1479  Handle<ExampleData> pIn;
1480  iEvent.getByLabel("example",pIn);
1481 #endif
1482 
1483 #ifdef THIS_IS_AN_EVENTSETUP_EXAMPLE
1484  ESHandle<SetupData> pSetup;
1485  iSetup.get<SetupRecord>().get(pSetup);
1486 #endif
1487 }
1488 
1489 
1490 // ------------ method called once each job just before starting event loop ------------
1491 void
1493 {
1494 }
1495 
1496 // ------------ method called once each job just after ending the event loop ------------
1497 void
1499 
1500  theFile->cd();
1501 
1502  for (int i=0; i<nphimx; i++) {
1503  for (int j=0; j<netamx; j++) {
1504 
1505  nevsigrsg->Fill(netamx*i+j, sigrsg[j][i]->GetEntries());
1506  mnsigrsg->Fill(netamx*i+j, sigrsg[j][i]->GetMean());
1507  rmssigrsg->Fill(netamx*i+j, sigrsg[j][i]->GetRMS());
1508 
1509  nevcrossg->Fill(netamx*i+j, crossg[j][i]->GetEntries());
1510  mncrossg->Fill(netamx*i+j, crossg[j][i]->GetMean());
1511  rmscrossg->Fill(netamx*i+j, crossg[j][i]->GetRMS());
1512 
1513  if (m_correl) { //#ifdef CORREL
1514 
1515  nevcorrsglb->Fill(netamx*i+j, corrsglb[j][i]->GetEntries());
1516  mncorrsglb->Fill(netamx*i+j, corrsglb[j][i]->GetMean());
1517  rmscorrsglb->Fill(netamx*i+j, corrsglb[j][i]->GetRMS());
1518 
1519  nevcorrsgrb->Fill(netamx*i+j, corrsgrb[j][i]->GetEntries());
1520  mncorrsgrb->Fill(netamx*i+j, corrsgrb[j][i]->GetMean());
1521  rmscorrsgrb->Fill(netamx*i+j, corrsgrb[j][i]->GetRMS());
1522 
1523  nevcorrsglu->Fill(netamx*i+j, corrsglu[j][i]->GetEntries());
1524  mncorrsglu->Fill(netamx*i+j, corrsglu[j][i]->GetMean());
1525  rmscorrsglu->Fill(netamx*i+j, corrsglu[j][i]->GetRMS());
1526 
1527  nevcorrsgru->Fill(netamx*i+j, corrsgru[j][i]->GetEntries());
1528  mncorrsgru->Fill(netamx*i+j, corrsgru[j][i]->GetMean());
1529  rmscorrsgru->Fill(netamx*i+j, corrsgru[j][i]->GetRMS());
1530 
1531  nevcorrsgall->Fill(netamx*i+j, corrsgall[j][i]->GetEntries());
1532  mncorrsgall->Fill(netamx*i+j, corrsgall[j][i]->GetMean());
1533  rmscorrsgall->Fill(netamx*i+j, corrsgall[j][i]->GetRMS());
1534 
1535  nevcorrsgl->Fill(netamx*i+j, corrsgl[j][i]->GetEntries());
1536  mncorrsgl->Fill(netamx*i+j, corrsgl[j][i]->GetMean());
1537  rmscorrsgl->Fill(netamx*i+j, corrsgl[j][i]->GetRMS());
1538 
1539  nevcorrsgr->Fill(netamx*i+j, corrsgr[j][i]->GetEntries());
1540  mncorrsgr->Fill(netamx*i+j, corrsgr[j][i]->GetMean());
1541  rmscorrsgr->Fill(netamx*i+j, corrsgr[j][i]->GetRMS());
1542  } //m_correl #endif
1543  if (m_checkmap) { //#ifdef CHECKMAP
1544  nevcorrsgc->Fill(netamx*i+j, corrsgc[j][i]->GetEntries());
1545  mncorrsgc->Fill(netamx*i+j, corrsgc[j][i]->GetMean());
1546  rmscorrsgc->Fill(netamx*i+j, corrsgc[j][i]->GetRMS());
1547  } //m_checkmap #endif
1548  }
1549  }
1550 
1551  if (m_combined) { // #ifdef COMBINED
1552  for (int j=0; j<ringmx; j++) {
1553  for (int i=0; i<routmx; i++) {
1554  if (j!=2 && i>=rout12mx) continue;
1555  nevsigrsg->Fill(netamx*nphimx+ringmx*i+j, com_sigrsg[j][i]->GetEntries());
1556  mnsigrsg->Fill(netamx*nphimx+ringmx*i+j, com_sigrsg[j][i]->GetMean());
1557  rmssigrsg->Fill(netamx*nphimx+ringmx*i+j, com_sigrsg[j][i]->GetRMS());
1558 
1559  nevcrossg->Fill(netamx*nphimx+ringmx*i+j, com_crossg[j][i]->GetEntries());
1560  mncrossg->Fill(netamx*nphimx+ringmx*i+j, com_crossg[j][i]->GetMean());
1561  rmscrossg->Fill(netamx*nphimx+ringmx*i+j, com_crossg[j][i]->GetRMS());
1562  }
1563  }
1564 
1565  for (int i=0; i<sectmx; i++) {
1566  for (int j=0; j<ringmx; j++) {
1567  if (m_correl) { // #ifdef CORREL
1568  nevcorrsglb->Fill(netamx*nphimx+ringmx*i+j, com_corrsglb[j][i]->GetEntries());
1569  mncorrsglb->Fill(netamx*nphimx+ringmx*i+j, com_corrsglb[j][i]->GetMean());
1570  rmscorrsglb->Fill(netamx*nphimx+ringmx*i+j, com_corrsglb[j][i]->GetRMS());
1571 
1572  nevcorrsgrb->Fill(netamx*nphimx+ringmx*i+j, com_corrsgrb[j][i]->GetEntries());
1573  mncorrsgrb->Fill(netamx*nphimx+ringmx*i+j, com_corrsgrb[j][i]->GetMean());
1574  rmscorrsgrb->Fill(netamx*nphimx+ringmx*i+j, com_corrsgrb[j][i]->GetRMS());
1575 
1576  nevcorrsglu->Fill(netamx*nphimx+ringmx*i+j, com_corrsglu[j][i]->GetEntries());
1577  mncorrsglu->Fill(netamx*nphimx+ringmx*i+j, com_corrsglu[j][i]->GetMean());
1578  rmscorrsglu->Fill(netamx*nphimx+ringmx*i+j, com_corrsglu[j][i]->GetRMS());
1579 
1580  nevcorrsgru->Fill(netamx*nphimx+ringmx*i+j, com_corrsgru[j][i]->GetEntries());
1581  mncorrsgru->Fill(netamx*nphimx+ringmx*i+j, com_corrsgru[j][i]->GetMean());
1582  rmscorrsgru->Fill(netamx*nphimx+ringmx*i+j, com_corrsgru[j][i]->GetRMS());
1583 
1584  nevcorrsgall->Fill(netamx*nphimx+ringmx*i+j, com_corrsgall[j][i]->GetEntries());
1585  mncorrsgall->Fill(netamx*nphimx+ringmx*i+j, com_corrsgall[j][i]->GetMean());
1586  rmscorrsgall->Fill(netamx*nphimx+ringmx*i+j, com_corrsgall[j][i]->GetRMS());
1587 
1588  nevcorrsgl->Fill(netamx*nphimx+ringmx*i+j, com_corrsgl[j][i]->GetEntries());
1589  mncorrsgl->Fill(netamx*nphimx+ringmx*i+j, com_corrsgl[j][i]->GetMean());
1590  rmscorrsgl->Fill(netamx*nphimx+ringmx*i+j, com_corrsgl[j][i]->GetRMS());
1591 
1592  nevcorrsgr->Fill(netamx*nphimx+ringmx*i+j, com_corrsgr[j][i]->GetEntries());
1593  mncorrsgr->Fill(netamx*nphimx+ringmx*i+j, com_corrsgr[j][i]->GetMean());
1594  rmscorrsgr->Fill(netamx*nphimx+ringmx*i+j, com_corrsgr[j][i]->GetRMS());
1595  } //m_correl #endif
1596  if (m_checkmap) { // #ifdef CHECKMAP
1597  nevcorrsgc->Fill(netamx*nphimx+ringmx*i+j, com_corrsgc[j][i]->GetEntries());
1598  mncorrsgc->Fill(netamx*nphimx+ringmx*i+j, com_corrsgc[j][i]->GetMean());
1599  rmscorrsgc->Fill(netamx*nphimx+ringmx*i+j, com_corrsgc[j][i]->GetRMS());
1600  } //m_checkmap #endif
1601  }
1602  }
1603  } //m_combined #endif
1604 
1605 
1606  for (int i=1; i<neffip; i++) {
1607  sig_effi[i]->Divide(sig_effi[0]);
1608  }
1609  for (int ij=0; ij<netamx; ij++) {
1610  for (int jk = 0; jk <nphimx; jk++) {
1611  int ieta = (ij<15) ? ij+1 : 14-ij;
1612  int iphi = jk+1;
1613  double signal = sigrsg[ij][jk]->GetMean();
1614  mean_energy->Fill(ieta, iphi, signal);
1615  }
1616  }
1617 
1618  int irunold = irun;
1619 
1620 
1621  gStyle->SetOptLogy(0);
1622  gStyle->SetTitleFillColor(10);
1623  gStyle->SetStatColor(10);
1624 
1625  gStyle->SetCanvasColor(10);
1626  gStyle->SetOptStat(0); //1110);
1627  gStyle->SetOptTitle(1);
1628 
1629  gStyle->SetTitleColor(10);
1630  gStyle->SetTitleFontSize(0.09);
1631  gStyle->SetTitleOffset(-0.05);
1632  gStyle->SetTitleBorderSize(1);
1633 
1634  gStyle->SetPadColor(10);
1635  gStyle->SetPadBorderMode(0);
1636  gStyle->SetStatColor(10);
1637  gStyle->SetPadBorderMode(0);
1638  gStyle->SetStatBorderSize(1);
1639  gStyle->SetStatFontSize(.07);
1640 
1641  gStyle->SetStatStyle(1001);
1642  gStyle->SetOptFit(101);
1643  gStyle->SetCanvasColor(10);
1644  gStyle->SetCanvasBorderMode(0);
1645 
1646  gStyle->SetStatX(.99);
1647  gStyle->SetStatY(.99);
1648  gStyle->SetStatW(.45);
1649  gStyle->SetStatH(.16);
1650  gStyle->SetLabelSize(0.075,"XY");
1651  gStyle->SetLabelOffset(0.21,"XYZ");
1652  gStyle->SetTitleSize(0.065,"XY");
1653  gStyle->SetTitleOffset(0.06,"XYZ");
1654  gStyle->SetPadTopMargin(.09);
1655  gStyle->SetPadBottomMargin(0.11);
1656  gStyle->SetPadLeftMargin(0.12);
1657  gStyle->SetPadRightMargin(0.15);
1658  gStyle->SetPadGridX(3);
1659  gStyle->SetPadGridY(3);
1660  gStyle->SetGridStyle(2);
1661  gStyle->SetNdivisions(303,"XY");
1662 
1663  gStyle->SetMarkerSize(0.60);
1664  gStyle->SetMarkerColor(2);
1665  gStyle->SetMarkerStyle(20);
1666  gStyle->SetTitleFontSize(0.07);
1667 
1668  char out_file[200];
1669  int xsiz = 700;
1670  int ysiz = 500;
1671 
1672  TCanvas *c2 = new TCanvas("c2", "Statistics and efficiency", xsiz, ysiz);
1673  c2->Divide(2,1); //(3,2);
1674  for (int ij=0; ij<neffip; ij=ij+3) {
1675  sig_effi[ij]->GetXaxis()->SetTitle("#eta");
1676  sig_effi[ij]->GetXaxis()->SetTitleSize(0.075);
1677  sig_effi[ij]->GetXaxis()->SetTitleOffset(0.65); //0.85
1678  sig_effi[ij]->GetXaxis()->CenterTitle();
1679  sig_effi[ij]->GetXaxis()->SetLabelSize(0.055);
1680  sig_effi[ij]->GetXaxis()->SetLabelOffset(0.001);
1681 
1682  sig_effi[ij]->GetYaxis()->SetTitle("#phi");
1683  sig_effi[ij]->GetYaxis()->SetTitleSize(0.075);
1684  sig_effi[ij]->GetYaxis()->SetTitleOffset(0.9);
1685  sig_effi[ij]->GetYaxis()->CenterTitle();
1686  sig_effi[ij]->GetYaxis()->SetLabelSize(0.055);
1687  sig_effi[ij]->GetYaxis()->SetLabelOffset(0.01);
1688 
1689  c2->cd(int(ij/3.)+1); sig_effi[ij]->Draw("colz");
1690  }
1691  sprintf(out_file, "comb_hosig_evt_%i.jpg",irunold);
1692  c2->SaveAs(out_file);
1693 
1694  gStyle->SetTitleFontSize(0.045);
1695  gStyle->SetPadRightMargin(0.1);
1696  gStyle->SetPadLeftMargin(0.1);
1697  gStyle->SetPadBottomMargin(0.12);
1698 
1699  TCanvas *c1 = new TCanvas("c1", "Mean signal in each tower", xsiz, ysiz);
1700 
1701  mean_energy->GetXaxis()->SetTitle("#eta");
1702  mean_energy->GetXaxis()->SetTitleSize(0.075);
1703  mean_energy->GetXaxis()->SetTitleOffset(0.65); //0.6
1704  mean_energy->GetXaxis()->CenterTitle();
1705  mean_energy->GetXaxis()->SetLabelSize(0.045);
1706  mean_energy->GetXaxis()->SetLabelOffset(0.001);
1707 
1708  mean_energy->GetYaxis()->SetTitle("#phi");
1709  mean_energy->GetYaxis()->SetTitleSize(0.075);
1710  mean_energy->GetYaxis()->SetTitleOffset(0.5);
1711  mean_energy->GetYaxis()->CenterTitle();
1712  mean_energy->GetYaxis()->SetLabelSize(0.045);
1713  mean_energy->GetYaxis()->SetLabelOffset(0.01);
1714 
1715  mean_energy->Draw("colz");
1716  sprintf(out_file, "homean_energy_%i.jpg",irunold);
1717  c1->SaveAs(out_file);
1718 
1719  delete c1;
1720  delete c2;
1721 
1722  gStyle->SetPadBottomMargin(0.14);
1723  gStyle->SetPadLeftMargin(0.17);
1724  gStyle->SetPadRightMargin(0.03);
1725 
1726  gStyle->SetOptStat(1110);
1727 
1728  const int nsample =8;
1729  TF1* gx0[nsample]={0};
1730  TF1* ped0fun[nsample]={0};
1731  TF1* signal[nsample]={0};
1732  TF1* pedfun[nsample]={0};
1733  TF1* sigfun[nsample]={0};
1734  TF1* signalx[nsample]={0};
1735 
1736  TH1F* signall[nsample]={0};
1737  TH1F* pedstll[nsample]={0};
1738 
1739  if (m_constant) {
1740 
1741  gStyle->SetOptFit(101);
1742  gStyle->SetCanvasBorderMode(0);
1743  gStyle->SetPadBorderMode(0);
1744  gStyle->SetStatBorderSize(1);
1745  gStyle->SetStatStyle(1001);
1746  gStyle->SetTitleColor(10);
1747  gStyle->SetTitleFontSize(0.09);
1748  gStyle->SetTitleOffset(-0.05);
1749  gStyle->SetTitleBorderSize(1);
1750 
1751  gStyle->SetCanvasColor(10);
1752  gStyle->SetPadColor(10);
1753  gStyle->SetStatColor(10);
1754  gStyle->SetStatFontSize(.07);
1755  gStyle->SetStatX(0.99);
1756  gStyle->SetStatY(0.99);
1757  gStyle->SetStatW(0.30);
1758  gStyle->SetStatH(0.10);
1759  gStyle->SetTitleSize(0.065,"XYZ");
1760  gStyle->SetLabelSize(0.075,"XYZ");
1761  gStyle->SetLabelOffset(0.012,"XYZ");
1762  gStyle->SetPadGridX(1);
1763  gStyle->SetPadGridY(1);
1764  gStyle->SetGridStyle(3);
1765  gStyle->SetNdivisions(101,"XY");
1766  gStyle->SetOptLogy(0);
1767  int iiter = 0;
1768 
1769 
1770  ofstream file_out(theoutputtxtFile.c_str());
1771  // TPostScript* ps=0;
1772  int ips=111;
1773  TPostScript ps(theoutputpsFile.c_str(),ips);
1774  ps.Range(20,28);
1775 
1776  xsiz = 900; //900;
1777  ysiz = 1200; //600;
1778  TCanvas *c0 = new TCanvas("c0", " Pedestal vs signal", xsiz, ysiz);
1779 
1780  float mean_eta[nphimx];
1781  float mean_phi[netamx];
1782  float rms_eta[nphimx];
1783  float rms_phi[netamx];
1784 
1785  for (int ij=0; ij<nphimx; ij++) {mean_phi[ij] = rms_phi[ij] =0;}
1786  for (int ij=0; ij<netamx; ij++) {mean_eta[ij] = rms_eta[ij] =0;}
1787 
1788  int mxeta = 0;
1789  int mxphi = 0;
1790  int mneta = 0;
1791  int mnphi = 0;
1792 
1793  //iijj = 0 : Merging all ring
1794  // = 1 : Individual HPD
1795  //iijj = 2 : merging all phi
1796  // = 3 : Individual tower
1797 
1798  for (int iijj = 0; iijj <4; iijj++) {
1799  // if ((!mx_combined) && iijj==1) continue; //Use this only for combined data
1800  if (iijj==0){
1801  mxeta = ringmx; mxphi = 1; mneta = 0; mnphi = 0;
1802  } else if (iijj==1) {
1803  mxeta = ringmx; mxphi = routmx;
1804  mneta = 0; mnphi = 0;
1805  } else if (iijj==2) {
1806  mxeta = netamx; mxphi = 1; mneta = 0; mnphi = 0;
1807  } else if (iijj==3) {
1808  mxeta = netamx; mxphi = nphimx;
1809  mneta = 0; mnphi = 0;
1810  }
1811 
1812  for (int jk=mneta; jk<mxeta; jk++) {
1813  for (int ij=mnphi; ij<mxphi; ij++) {
1814  if (iijj==1) continue;
1815  if ((iijj==0 || iijj==1) && jk !=2 && ij >=rout12mx) continue;
1816  int izone = iiter%nsample;
1817 
1818  if (iijj==0) {
1819  int iread = (jk==2) ? routmx : rout12mx;
1820  signall[izone] = (TH1F*)com_sigrsg[jk][iread]->Clone("hnew");
1821  pedstll[izone] = (TH1F*)com_crossg[jk][iread]->Clone("hnew");
1822  } else if (iijj==1) {
1823  signall[izone] = (TH1F*)com_sigrsg[jk][ij]->Clone("hnew");
1824  pedstll[izone] = (TH1F*)com_crossg[jk][ij]->Clone("hnew");
1825  } else if (iijj==2) {
1826  signall[izone] = (TH1F*)sigrsg[jk][nphimx]->Clone("hnew");
1827  pedstll[izone] = (TH1F*)crossg[jk][nphimx]->Clone("hnew");
1828  } else if (iijj==3) {
1829  signall[izone] = (TH1F*)sigrsg[jk][ij]->Clone("hnew");
1830  pedstll[izone] = (TH1F*)crossg[jk][ij]->Clone("hnew");
1831  }
1832 
1833  pedstll[izone]->SetLineWidth(2);
1834  signall[izone]->SetLineWidth(2);
1835  pedstll[izone]->SetLineColor(2);
1836  signall[izone]->SetLineColor(4);
1837  pedstll[izone]->SetNdivisions(506,"XY");
1838  signall[izone]->SetNdivisions(506,"XY");
1839 
1840  signall[izone]->GetXaxis()->SetLabelSize(.065);
1841  signall[izone]->GetYaxis()->SetLabelSize(.06);
1842  if (m_digiInput) {
1843  signall[izone]->GetXaxis()->SetTitle("Signal (fC)");
1844  } else {
1845  signall[izone]->GetXaxis()->SetTitle("Signal (GeV)");
1846  }
1847  signall[izone]->GetXaxis()->SetTitleSize(.065);
1848  signall[izone]->GetXaxis()->CenterTitle();
1849 
1850  if (izone==0) { //iiter%8 ==0) {
1851  ps.NewPage();
1852  c0->Divide(4,4); //c0->Divide(2,4); // c0->Divide(1,2);
1853  }
1854  c0->cd(2*izone+1); // (iiter%8)+1); //c0->cd(iiter%8+1);
1855 
1856  /*
1857  if (iijj==0 && izone==0) {
1858  gStyle->SetOptLogy(1);
1859  gStyle->SetOptStat(0);
1860  gStyle->SetOptFit(0);
1861  c0->Divide(3,2);
1862  }
1863 
1864  if (iijj>0) {
1865  gStyle->SetOptLogy(0);
1866  gStyle->SetOptStat(1110);
1867  gStyle->SetOptFit(101);
1868 
1869  if (iiter==0) {
1870  int ips=111;
1871  ps = new TPostScript(theoutputpsFile.c_str(),ips);
1872  ps.Range(20,28);
1873  xsiz = 900; //900;
1874  ysiz = 1200; //600;
1875  c0 = new TCanvas("c0", " Pedestal vs signal", xsiz, ysiz);
1876  }
1877  if (izone==0) {
1878  ps.NewPage();
1879  c0->Divide(4,4);
1880  }
1881  }
1882  if (iijj==0) {c0->cd(izone+1); } else { c0->cd(2*izone+1);}
1883  */
1884 
1885  float mean = pedstll[izone]->GetMean();
1886  float rms = pedstll[izone]->GetRMS();
1887  if (m_digiInput) {
1888  if (rms <0.6) rms = 0.6;
1889  if (rms >1.2) rms = 1.2;
1890  if (mean >1.2) mean = 1.2;
1891  if (mean <-1.2) mean = -1.2;
1892  } else {
1893  if (rms <0.10) rms = 0.10;
1894  if (rms >0.15) rms=0.15;
1895  if (mean >0.20) mean = 0.20;
1896  if (mean <-0.20) mean = -0.20;
1897  }
1898  float xmn = mean-6.*rms;
1899  float xmx = mean+6.*rms;
1900 
1901  binwid = pedstll[izone]->GetBinWidth(1);
1902  if (xmx > pedstll[izone]->GetXaxis()->GetXmax()) xmx = pedstll[izone]->GetXaxis()->GetXmax()-0.5*binwid;
1903  if (xmn < pedstll[izone]->GetXaxis()->GetXmin()) xmn = pedstll[izone]->GetXaxis()->GetXmin()+0.5*binwid;
1904 
1905  float height = pedstll[izone]->GetEntries();
1906 
1907  double par[nbgpr] ={height, mean, 0.75*rms};
1908 
1909  double gaupr[nbgpr];
1910  double parer[nbgpr];
1911 
1912  ietafit = jk;
1913  iphifit = ij;
1914  pedstll[izone]->GetXaxis()->SetLabelSize(.065);
1915  pedstll[izone]->GetYaxis()->SetLabelSize(.06);
1916 
1917  // if (iijj==0) {
1918  // pedstll[izone]->GetXaxis()->SetRangeUser(alow, ahigh);
1919  // } else {
1920  pedstll[izone]->GetXaxis()->SetRangeUser(xmn, xmx);
1921  // }
1922 
1923  if (m_digiInput) {
1924  if (iijj==0) {
1925  pedstll[izone]->GetXaxis()->SetTitle("Pedestal/Signal (fC)");
1926  } else {
1927  pedstll[izone]->GetXaxis()->SetTitle("Pedestal (fC)");
1928  }
1929  } else {
1930  if (iijj==0) {
1931  pedstll[izone]->GetXaxis()->SetTitle("Pedestal/Signal (GeV)");
1932  } else {
1933  pedstll[izone]->GetXaxis()->SetTitle("Pedestal (GeV)");
1934  }
1935  }
1936  pedstll[izone]->GetXaxis()->SetTitleSize(.065);
1937  pedstll[izone]->GetXaxis()->CenterTitle();
1938  // pedstll[izone]->SetLineWidth(2);
1939 
1940  pedstll[izone]->Draw();
1941  if (m_pedsuppr && !m_digiInput) {
1942  gaupr[0] = 0;
1943  gaupr[1] = 0.0; // pedmean[ietafit][iphifit];
1944  if (m_digiInput) {
1945  gaupr[2] = 0.90; //GMA need from database pedwidth[ietafit][iphifit];
1946  } else {
1947  gaupr[2] = 0.15; //GMA need from database
1948  }
1949  parer[0] = parer[1] = parer[2] = 0;
1950  } else {
1951 
1952  if (pedstll[izone]->GetEntries() >5) {
1953 
1954  if ((iijj!=3) || m_histfit) {
1955  char temp[20];
1956  sprintf(temp, "gx0_%i",izone);
1957  gx0[izone] = new TF1(temp, gausX, xmn, xmx, nbgpr);
1958  gx0[izone]->SetParameters(par);
1959  gx0[izone]->SetLineWidth(1);
1960  pedstll[izone]->Fit(gx0[izone], "R+");
1961 
1962  for (int k=0; k<nbgpr; k++) {
1963  parer[k] = gx0[izone]->GetParError(k);
1964  gaupr[k] = gx0[izone]->GetParameter(k);
1965  }
1966  } else {
1967  double strt[nbgpr] = {height, mean, 0.75*rms};
1968  double step[nbgpr] = {1.0, 0.001, 0.001};
1969  double alowmn[nbgpr] = {0.5*height, mean-rms, 0.3*rms};
1970  double ahighmn[nbgpr] ={1.5*height, mean+rms, 1.5*rms};
1971 
1972  TMinuit *gMinuit = new TMinuit(nbgpr);
1973  gMinuit->SetFCN(fcnbg);
1974 
1975  double arglist[10];
1976  int ierflg = 0;
1977  arglist[0] =0.5;
1978  gMinuit->mnexcm("SET ERR", arglist, 1, ierflg);
1979  char name[100];
1980  for (int k=0; k<nbgpr; k++) {
1981  sprintf(name, "pedpar%i",k);
1982  gMinuit->mnparm(k, name, strt[k], step[k], alowmn[k], ahighmn[k],ierflg);
1983  }
1984 
1985  arglist[0] = 0;
1986  gMinuit->mnexcm("SIMPLEX", arglist, 0, ierflg);
1987 
1988  arglist[0] = 0;
1989  gMinuit->mnexcm("IMPROVE", arglist, 0, ierflg);
1990 
1991  TString chnam;
1992  double parv,err,xlo,xup, plerr, mierr, eparab, gcc;
1993  int iuit;
1994 
1995  for (int k=0; k<nbgpr; k++) {
1996  if (step[k] >-10) {
1997  gMinuit->mnpout(k, chnam, parv, err, xlo, xup, iuit);
1998  gMinuit->mnerrs(k, plerr, mierr, eparab, gcc);
1999  // cout <<"k "<< k<<" "<<chnam<<" "<<parv<<" "<<err<<" "<<xlo<<" "<<xup<<" "<<plerr<<" "<<mierr<<" "<<eparab<<endl;
2000  if (k==0) {
2001  gaupr[k] = parv*binwid;
2002  parer[k] = err*binwid;
2003  } else {
2004  gaupr[k] = parv;
2005  parer[k] = err;
2006  }
2007  }
2008  }
2009 
2010  // gx0[izone]->SetParameters(gaupr);
2011 
2012  char temp[20];
2013  sprintf(temp, "ped0fun_%i",izone);
2014  ped0fun[izone] = new TF1(temp, gausX, xmn, xmx, nbgpr);
2015  ped0fun[izone]->SetParameters(gaupr);
2016  ped0fun[izone]->SetLineColor(3);
2017  ped0fun[izone]->SetLineWidth(1);
2018  ped0fun[izone]->Draw("same");
2019 
2020  delete gMinuit;
2021  }
2022  } else {
2023  for (int k=0; k<nbgpr; k++) {gaupr[k] = par[k]; }
2024  if (m_digiInput) { gaupr[2] = 0.90; } else { gaupr[2] = 0.15;}
2025  }
2026  }
2027  // if (iijj!=0)
2028  c0->cd(2*izone+2);
2029  if (signall[izone]->GetEntries() >5) {
2030  Double_t parall[nsgpr];
2031  double parserr[nsgpr];
2032  double fitres[nsgpr];
2033  double pedht = 0;
2034 
2035  char temp[20];
2036  sprintf(temp, "signal_%i",izone);
2037  xmn = signall[izone]->GetXaxis()->GetXmin();
2038  xmx = 0.5*signall[izone]->GetXaxis()->GetXmax();
2039  signal[izone] = new TF1(temp, totalfunc, xmn, xmx, nsgpr);
2040  xmx *=2.0;
2041  if ((iijj!=3) || m_histfit) {
2042  pedht = (signall[izone]->GetBinContent(nbn-1)+
2043  signall[izone]->GetBinContent(nbn)+
2044  signall[izone]->GetBinContent(nbn+1))/3.;
2045 
2046  if (m_pedsuppr && !m_digiInput) {
2047  parall[1] = 0.0; // pedmean[ietafit][iphifit];
2048  if (m_digiInput) { parall[2] = 0.90; } else { parall[2] = 0.15;}
2049  } else {
2050  for (int i=0; i<nbgpr; i++) {parall[i] = gaupr[i];}
2051  }
2052 
2053  set_mean(parall[1], m_digiInput);
2054  set_sigma(parall[2], m_digiInput);
2055 
2056  parall[0] = 0.9*pedht; //GM for Z-mumu, there is almost no pedestal
2057  parall[3] = 0.14;
2058  double area = binwid*signall[izone]->GetEntries();
2059  parall[5]= area;
2060 
2061  if (iijj==3) {
2062  parall[4] = fitprm[4][jk];
2063  parall[6] = fitprm[6][jk];
2064  } else {
2065  parall[4] = signall[izone]->GetMean();
2066  parall[6]=parall[2];
2067  }
2068 
2069  signal[izone]->SetParameters(parall);
2070  signal[izone]->FixParameter(1, parall[1]);
2071  signal[izone]->FixParameter(2, parall[2]);
2072  signal[izone]->SetParLimits(0, 0.00, 2.0*pedht+0.1);
2073  signal[izone]->FixParameter(3, 0.14);
2074 
2075  signal[izone]->SetParLimits(5, 0.40*area, 1.15*area);
2076  // if (m_histfit) { //GMA
2077  if (iijj==3) {
2078  signal[izone]->SetParLimits(4, 0.2*fitprm[4][jk], 2.0*fitprm[4][jk]);
2079  signal[izone]->SetParLimits(6, 0.2*fitprm[6][jk], 2.0*fitprm[6][jk]);
2080  } else {
2081  if (m_digiInput) {
2082  signal[izone]->SetParLimits(4, 0.6, 6.0);
2083  signal[izone]->SetParLimits(6, 0.60, 3.0);
2084  } else {
2085  signal[izone]->SetParLimits(4, 0.1, 1.0);
2086  signal[izone]->SetParLimits(6, 0.035, 0.3);
2087  }
2088  }
2089  signal[izone]->SetParNames("const", "mean", "sigma","Width","MP","Area","GSigma");
2090  signall[izone]->Fit(signal[izone], "0R+");
2091 
2092  signall[izone]->GetXaxis()->SetRangeUser(xmn,xmx);
2093  for (int k=0; k<nsgpr; k++) {
2094  fitres[k] = fitprm[k][jk] = signal[izone]->GetParameter(k);
2095  parserr[k] = signal[izone]->GetParError(k);
2096  }
2097 
2098  } else {
2099  double pedhtx = 0;
2100  for (unsigned i =0; i<sig_reg[ietafit][iphifit].size(); i++) {
2101  if (sig_reg[ietafit][iphifit][i] >gaupr[1]-3*gaupr[2] && sig_reg[ietafit][iphifit][i]<gaupr[1]+gaupr[2]) pedhtx++;
2102  }
2103 
2104  set_mean(gaupr[1], m_digiInput);
2105  set_sigma(gaupr[2], m_digiInput);
2106 
2107  TString name[nsgpr] = {"const", "mean", "sigma","Width","MP","Area","GSigma"};
2108  double strt[nsgpr] = {0.9*pedhtx, gaupr[1], gaupr[2], fitprm[3][jk], fitprm[4][jk], signall[izone]->GetEntries(), fitprm[6][jk]};
2109  double alowmn[nsgpr] = {0.1*pedhtx-0.1, gaupr[1]-0.1, gaupr[2]-0.1,0.07, 0.2*strt[4], 0.1*strt[5], 0.2*strt[6]};
2110  double ahighmn[nsgpr] ={1.2*pedhtx+0.1, gaupr[1]+0.1, gaupr[2]+0.1,0.20, 2.5*strt[4], 1.5*strt[5], 2.2*strt[6]};
2111  double step[nsgpr] = {1.0, 0.0, 0.0, 0.0, 0.001, 1.0, 0.002};
2112 
2113  TMinuit *gMinuit = new TMinuit(nsgpr);
2114  gMinuit->SetFCN(fcnsg);
2115 
2116  double arglist[10];
2117  int ierflg = 0;
2118  arglist[0] =0.5;
2119  gMinuit->mnexcm("SET ERR", arglist, 1, ierflg);
2120 
2121  for (int k=0; k<nsgpr; k++) {
2122  gMinuit->mnparm(k, name[k], strt[k], step[k], alowmn[k], ahighmn[k],ierflg);
2123  }
2124 
2125  arglist[0] = 0;
2126  gMinuit->mnexcm("SIMPLEX", arglist, 0, ierflg);
2127 
2128  arglist[0] = 0;
2129  gMinuit->mnexcm("IMPROVE", arglist, 0, ierflg);
2130 
2131  TString chnam;
2132  double parv,err,xlo,xup, plerr, mierr, eparab, gcc;
2133  int iuit;
2134 
2135  for (int k=0; k<nsgpr; k++) {
2136  if (step[k] >-10) {
2137  gMinuit->mnpout(k, chnam, parv, err, xlo, xup, iuit);
2138  gMinuit->mnerrs(k, plerr, mierr, eparab, gcc);
2139  if (k==0 || k==5) {
2140  fitres[k] = parv*binwid;
2141  parserr[k]= err*binwid;
2142  } else {
2143  fitres[k] = parv;
2144  parserr[k]= err;
2145  }
2146 
2147  }
2148  }
2149 
2150  delete gMinuit;
2151  }
2152 
2153  // if (iijj==0) {
2154  // signall[izone]->Draw("same");
2155  // } else {
2156  signall[izone]->Draw();
2157  // }
2158 
2159  sprintf(temp, "pedfun_%i",izone);
2160  pedfun[izone] = new TF1(temp, gausX, xmn, xmx, nbgpr);
2161  pedfun[izone]->SetParameters(fitres);
2162  pedfun[izone]->SetLineColor(3);
2163  pedfun[izone]->SetLineWidth(1);
2164  pedfun[izone]->Draw("same");
2165 
2166  sprintf(temp, "signalfun_%i",izone);
2167  sigfun[izone] = new TF1(temp, langaufun, xmn, xmx, nsgpr-nbgpr);
2168  sigfun[izone]->SetParameters(&fitres[3]);
2169  sigfun[izone]->SetLineWidth(1);
2170  sigfun[izone]->SetLineColor(4);
2171  sigfun[izone]->Draw("same");
2172 
2173  sprintf(temp, "total_%i",izone);
2174  signalx[izone] = new TF1(temp, totalfunc, xmn, xmx, nsgpr);
2175  signalx[izone]->SetParameters(fitres);
2176  signalx[izone]->SetLineWidth(1);
2177  signalx[izone]->Draw("same");
2178 
2179  int kl = (jk<15) ? jk+1 : 14-jk;
2180 
2181  cout<<"histinfo"<<iijj<<" fit "
2182  <<std::setw(3)<< kl<<" "
2183  <<std::setw(3)<< ij+1<<" "
2184  <<std::setw(5)<<pedstll[izone]->GetEntries()<<" "
2185  <<std::setw(6)<<pedstll[izone]->GetMean()<<" "
2186  <<std::setw(6)<<pedstll[izone]->GetRMS()<<" "
2187  <<std::setw(5)<<signall[izone]->GetEntries()<<" "
2188  <<std::setw(6)<<signall[izone]->GetMean()<<" "
2189  <<std::setw(6)<<signall[izone]->GetRMS()<<" "
2190  <<std::setw(6)<< signal[izone]->GetChisquare()<<" "
2191  <<std::setw(3)<< signal[izone]->GetNDF()<<endl;
2192 
2193  file_out<<"histinfo"<<iijj<<" fit "
2194  <<std::setw(3)<< kl<<" "
2195  <<std::setw(3)<< ij+1<<" "
2196  <<std::setw(5)<<pedstll[izone]->GetEntries()<<" "
2197  <<std::setw(6)<<pedstll[izone]->GetMean()<<" "
2198  <<std::setw(6)<<pedstll[izone]->GetRMS()<<" "
2199  <<std::setw(5)<<signall[izone]->GetEntries()<<" "
2200  <<std::setw(6)<<signall[izone]->GetMean()<<" "
2201  <<std::setw(6)<<signall[izone]->GetRMS()<<" "
2202  <<std::setw(6)<< signal[izone]->GetChisquare()<<" "
2203  <<std::setw(3)<< signal[izone]->GetNDF()<<endl;
2204 
2205  file_out <<"fitres x"<<iijj<<" "<<kl<<" "<<ij+1<<" "<< fitres[0]<<" "<< fitres[1]<<" "<< fitres[2]<<" "<< fitres[3]<<" "<< fitres[4]<<" "<< fitres[5]<<" "<< fitres[6]<<endl;
2206  file_out <<"parserr"<<iijj<<" "<<kl<<" "<<ij+1<<" "<< parserr[0]<<" "<< parserr[1]<<" "<< parserr[2]<<" "<< parserr[3]<<" "<< parserr[4]<<" "<< parserr[5]<<" "<< parserr[6]<<endl;
2207 
2208  double diff=fitres[4]-fitres[1];
2209  if (diff <=0) diff = 0.000001;
2210  double error=parserr[4]*parserr[4]+parer[2]*parer[2];
2211  error = pow(error,0.5);
2212 
2213  int ieta = (jk<15) ? (15+jk) : (29-jk);
2214  int ifl = nphimx*ieta + ij;
2215 
2216  if (iijj==3) {
2217  ped_evt->Fill(ifl,pedstll[izone]->GetEntries());
2218  ped_mean->Fill(ifl,gaupr[1]);
2219  ped_width->Fill(ifl,gaupr[2]);
2220  fit_chi->Fill(ifl,signal[izone]->GetChisquare());
2221  sig_evt->Fill(ifl, signall[izone]->GetEntries());
2222  fit_sigevt->Fill(ifl, fitres[5]);
2223  fit_bkgevt->Fill(ifl, fitres[0]*sqrt(2*acos(-1.))*gaupr[2]);
2224  sig_mean->Fill(ifl, fitres[4]);
2225  sig_diff->Fill(ifl, fitres[4]-fitres[1]);
2226  sig_width->Fill(ifl, fitres[3]);
2227  sig_sigma->Fill(ifl, fitres[6]);
2228  sig_meanerr->Fill(ifl, parserr[4]);
2229  if (fitres[4]-fitres[1] !=0) sig_meanerrp->Fill(ifl, 100*parserr[4]/(fitres[4]-fitres[1]));
2230  if (gaupr[2]!=0) sig_signf->Fill(ifl,(fitres[4]-fitres[1])/gaupr[2]);
2231 
2232  ped_statmean->Fill(ifl,pedstll[izone]->GetMean());
2233  sig_statmean->Fill(ifl,signall[izone]->GetMean());
2234  ped_rms->Fill(ifl,pedstll[izone]->GetRMS());
2235  sig_rms->Fill(ifl,signall[izone]->GetRMS());
2236  }
2237 
2238  if ((iijj==2) || (iijj==3) || (iijj==1)) {
2239  if (signall[izone]->GetEntries() >5 && fitres[4]>0.1) {
2240  //GMA need to put this==1 in future
2241  float fact=0.812;
2242  if (abs(kl)<=4) fact=0.895;
2243  if (!m_digiInput) fact *=0.19; //conversion factor for GeV/fC
2244 
2245  float fact2 = 0;
2246  if (iijj==2) fact2 = invang[jk][nphimx];
2247  if (iijj==3) fact2 = invang[jk][ij];
2248  if (iijj==1) fact2 = com_invang[jk][ij];
2249 
2250  float calibc = fact*fact2/(fitres[4]*signall[izone]->GetEntries());
2251  float caliberr= TMath::Abs(calibc*parserr[4]/max(0.001,fitres[4]));
2252 
2253  if (iijj==2) {
2254  int ieta = (jk<15) ? jk+1 : 14-jk;
2255  mean_phi_hst->Fill(ieta, calibc);
2256  mean_phi_hst->SetBinError(mean_phi_hst->FindBin(ieta), caliberr);
2257  file_out<<"intieta "<<jk<<" "<<ij<<" "<<ieta<<" "<<mean_phi_hst->FindBin(double(ieta))<<" "<<calibc<<" "<<caliberr<<endl;
2258  } else if (iijj==3) {
2259  const_eta[jk]->Fill(ij+1,calibc);
2260  const_eta[jk]->SetBinError(const_eta[jk]->FindBin(ij+1), caliberr);
2261 
2262  peak_eta[jk]->Fill(ij+1,fitres[4]);
2263  peak_eta[jk]->SetBinError(peak_eta[jk]->FindBin(ij+1),parserr[4]);
2264 
2265  int ieta = (jk<15) ? jk+1 : 14-jk;
2266  const_eta_phi->Fill(ieta, ij+1,calibc);
2267  file_out<<"intietax "<<jk<<" "<<ij<<" "<<ieta<<" "<<const_eta_phi->FindBin(ieta, ij+1)<<endl;
2268  if (caliberr >0) {
2269  const_eta_phi->SetBinError(const_eta_phi->FindBin(ieta, ij+1),caliberr);
2270 
2271  mean_eta[ij] +=calibc/(caliberr*caliberr);
2272  mean_phi[jk] +=calibc/(caliberr*caliberr);
2273 
2274  rms_eta[ij] +=1./(caliberr*caliberr);
2275  rms_phi[jk] +=1./(caliberr*caliberr);
2276 
2277  } else {
2278  const_eta_phi->SetBinError(const_eta_phi->FindBin(ieta, ij+1), 0.0);
2279  }
2280  } else if (iijj==1) {
2281  const_hpdrm[jk]->Fill(ij+1,calibc);
2282  const_hpdrm[jk]->SetBinError(const_hpdrm[jk]->FindBin(ij+1), caliberr);
2283 
2284  peak_hpdrm[jk]->Fill(ij+1,fitres[4]);
2285  peak_hpdrm[jk]->SetBinError(peak_hpdrm[jk]->FindBin(ij+1),parserr[4]);
2286  }
2287 
2288  file_out<<"HO 4 "<<iijj<<" "<< std::setw(3)<<kl<<" "<<std::setw(3)<<ij+1<<" "
2289  <<std::setw(7)<<calibc<<" "<<std::setw(7)<<caliberr<<endl;
2290  }
2291  }
2292 
2293  } else { //if (signall[izone]->GetEntries() >10) {
2294  signall[izone]->Draw();
2295  float varx = 0.000;
2296  int kl = (jk<15) ? jk+1 : 14-jk;
2297  file_out<<"histinfo"<<iijj<<" nof "
2298  <<std::setw(3)<< kl<<" "
2299  <<std::setw(3)<< ij+1<<" "
2300  <<std::setw(5)<<pedstll[izone]->GetEntries()<<" "
2301  <<std::setw(6)<<pedstll[izone]->GetMean()<<" "
2302  <<std::setw(6)<<pedstll[izone]->GetRMS()<<" "
2303  <<std::setw(5)<<signall[izone]->GetEntries()<<" "
2304  <<std::setw(6)<<signall[izone]->GetMean()<<" "
2305  <<std::setw(6)<<signall[izone]->GetRMS()<<" "
2306  <<std::setw(6)<< varx<<" "
2307  <<std::setw(3)<< varx<<endl;
2308 
2309  file_out <<"fitres x"<<iijj<<" "<<kl<<" "<<ij+1<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<endl;
2310  file_out <<"parserr"<<iijj<<" "<<kl<<" "<<ij+1<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<" "<< varx<<endl;
2311 
2312  }
2313  iiter++;
2314  if (iiter%nsample==0) {
2315  c0->Update();
2316 
2317  for (int kl=0; kl<nsample; kl++) {
2318  if (gx0[kl]) {delete gx0[kl];gx0[kl] = 0;}
2319  if (ped0fun[kl]) {delete ped0fun[kl];ped0fun[kl] = 0;}
2320  if (signal[kl]) {delete signal[kl];signal[kl] = 0;}
2321  if (pedfun[kl]) {delete pedfun[kl];pedfun[kl] = 0;}
2322  if (sigfun[kl]) {delete sigfun[kl];sigfun[kl] = 0;}
2323  if (signalx[kl]) {delete signalx[kl];signalx[kl] = 0;}
2324  if (signall[kl]) {delete signall[kl];signall[kl] = 0;}
2325  if (pedstll[kl]) {delete pedstll[kl];pedstll[kl] = 0;}
2326  }
2327 
2328  }
2329  } //for (int jk=0; jk<netamx; jk++) {
2330  } //for (int ij=0; ij<nphimx; ij++) {
2331 
2332  // if (iijj==0) {
2333  // sprintf(out_file, "comb_hosig_allring_%i.jpg", irunold);
2334  // c0->SaveAs(out_file);
2335  // iiter = 0;
2336  // } else {
2337  // // c0->Update();
2338  // }
2339 
2340  // iiter = 0;
2341  } //end of iijj
2342  if (iiter%nsample!=0) {
2343  c0->Update();
2344  for (int kl=0; kl<nsample; kl++) {
2345  if (gx0[kl]) {delete gx0[kl];gx0[kl] = 0;}
2346  if (ped0fun[kl]) {delete ped0fun[kl];ped0fun[kl] = 0;}
2347  if (signal[kl]) {delete signal[kl];signal[kl] = 0;}
2348  if (pedfun[kl]) {delete pedfun[kl];pedfun[kl] = 0;}
2349  if (sigfun[kl]) {delete sigfun[kl];sigfun[kl] = 0;}
2350  if (signalx[kl]) {delete signalx[kl];signalx[kl] = 0;}
2351  if (signall[kl]) {delete signall[kl];signall[kl] = 0;}
2352  if (pedstll[kl]) {delete pedstll[kl];pedstll[kl] = 0;}
2353  }
2354  }
2355 
2356  delete c0;
2357 
2358  xsiz = 600; //int xsiz = 600;
2359  ysiz = 800; //int ysiz = 800;
2360 
2361  gStyle->SetTitleFontSize(0.05);
2362  gStyle->SetTitleSize(0.025,"XYZ");
2363  gStyle->SetLabelSize(0.025,"XYZ");
2364  gStyle->SetStatFontSize(.045);
2365 
2366  gStyle->SetOptStat(0);
2367  ps.NewPage(); TCanvas *c1 = new TCanvas("c1", " Pedestal vs signal", xsiz, ysiz);
2368  ped_evt->Draw(); c1->Update();
2369 
2370  ps.NewPage();
2371  ped_statmean->Draw(); c1->Update();
2372 
2373  ps.NewPage();
2374  ped_rms->Draw(); c1->Update();
2375 
2376  ps.NewPage();
2377  ped_mean->Draw(); c1->Update();
2378 
2379  ps.NewPage();
2380  ped_width->Draw(); c1->Update();
2381 
2382  ps.NewPage();
2383  sig_evt->Draw(); c1->Update();
2384 
2385  ps.NewPage();
2386  sig_statmean->Draw(); c1->Update();
2387 
2388  ps.NewPage();
2389  sig_rms->Draw(); c1->Update();
2390 
2391  ps.NewPage();
2392  fit_chi->Draw(); c1->Update();
2393 
2394  ps.NewPage();
2395  fit_sigevt->Draw(); c1->Update();
2396 
2397  ps.NewPage();
2398  fit_bkgevt->Draw(); c1->Update();
2399 
2400  ps.NewPage();
2401  sig_mean->Draw(); c1->Update();
2402 
2403  ps.NewPage();
2404  sig_width->Draw(); c1->Update();
2405 
2406  ps.NewPage();
2407  sig_sigma->Draw(); c1->Update();
2408 
2409  ps.NewPage();
2410  sig_meanerr->Draw(); c1->Update();
2411 
2412  ps.NewPage();
2413  sig_meanerrp->Draw(); c1->Update();
2414 
2415  ps.NewPage();
2416  sig_signf->Draw(); c1->Update();
2417 
2418  ps.Close();
2419  delete c1;
2420 
2421  file_out.close();
2422 
2423  if (m_figure) {
2424  xsiz = 700;
2425  ysiz = 450;
2426 
2427  gStyle->SetTitleFontSize(0.09);
2428  gStyle->SetPadBottomMargin(0.17);
2429  gStyle->SetPadLeftMargin(0.18);
2430  gStyle->SetPadRightMargin(0.01);
2431  gStyle->SetOptLogy(0);
2432  gStyle->SetOptStat(0);
2433 
2434  TCanvas *c2 = new TCanvas("c2", "runfile", xsiz, ysiz);
2435  c2->Divide(5,3);
2436 
2437  for (int side=0; side <2; side++) {
2438  gStyle->SetNdivisions(303,"XY");
2439  gStyle->SetPadRightMargin(0.01);
2440  int nmn = 0;
2441  int nmx = netamx/2;
2442  if (side==1) {
2443  nmn = netamx/2;
2444  nmx = netamx;
2445  }
2446 
2447  int nzone = 0;
2448 
2449  for (int ij=nmn; ij<nmx; ij++) {
2450 
2451  c2->cd(nzone+1);
2452  const_eta[ij]->GetXaxis()->SetTitle("#phi index");
2453  const_eta[ij]->GetXaxis()->SetTitleSize(.08);
2454  const_eta[ij]->GetXaxis()->CenterTitle();
2455  const_eta[ij]->GetXaxis()->SetTitleOffset(0.9);
2456  const_eta[ij]->GetXaxis()->SetLabelSize(.085);
2457  const_eta[ij]->GetXaxis()->SetLabelOffset(.01);
2458 
2459  const_eta[ij]->GetYaxis()->SetLabelSize(.08);
2460  const_eta[ij]->GetYaxis()->SetLabelOffset(.01);
2461  if (m_digiInput) {
2462  const_eta[ij]->GetYaxis()->SetTitle("GeV/fC");
2463  } else {
2464  const_eta[ij]->GetYaxis()->SetTitle("GeV/MIP-GeV!!");
2465  }
2466 
2467  const_eta[ij]->GetYaxis()->SetTitleSize(.085);
2468  const_eta[ij]->GetYaxis()->CenterTitle();
2469  const_eta[ij]->GetYaxis()->SetTitleOffset(1.3);
2470  const_eta[ij]->SetMarkerSize(0.60);
2471  const_eta[ij]->SetMarkerColor(2);
2472  const_eta[ij]->SetMarkerStyle(20);
2473 
2474 
2475  const_eta[ij]->Draw();
2476  nzone++;
2477  }
2478 
2479  sprintf(out_file, "calibho_%i_side%i.eps", irunold, side);
2480  c2->SaveAs(out_file);
2481 
2482  sprintf(out_file, "calibho_%i_side%i.jpg", irunold, side);
2483  c2->SaveAs(out_file);
2484 
2485  nzone = 0;
2486  for (int ij=nmn; ij<nmx; ij++) {
2487  c2->cd(nzone+1);
2488  peak_eta[ij]->GetXaxis()->SetTitle("#phi index");
2489  peak_eta[ij]->GetXaxis()->SetTitleSize(.08);
2490  peak_eta[ij]->GetXaxis()->CenterTitle();
2491  peak_eta[ij]->GetXaxis()->SetTitleOffset(0.90);
2492  peak_eta[ij]->GetXaxis()->SetLabelSize(.08);
2493  peak_eta[ij]->GetXaxis()->SetLabelOffset(.01);
2494 
2495  peak_eta[ij]->GetYaxis()->SetLabelSize(.08);
2496  peak_eta[ij]->GetYaxis()->SetLabelOffset(.01);
2497  if (m_digiInput) {
2498  peak_eta[ij]->GetYaxis()->SetTitle("fC");
2499  } else {
2500  peak_eta[ij]->GetYaxis()->SetTitle("GeV");
2501  }
2502 
2503  peak_eta[ij]->GetYaxis()->SetTitleSize(.085);
2504  peak_eta[ij]->GetYaxis()->CenterTitle();
2505  peak_eta[ij]->GetYaxis()->SetTitleOffset(1.3);
2506 
2507  peak_eta[ij]->SetMarkerSize(0.60);
2508  peak_eta[ij]->SetMarkerColor(2);
2509  peak_eta[ij]->SetMarkerStyle(20);
2510 
2511  peak_eta[ij]->Draw();
2512  nzone++;
2513  }
2514 
2515  sprintf(out_file, "peakho_%i_side%i.eps", irunold, side);
2516  c2->SaveAs(out_file);
2517 
2518  sprintf(out_file, "peakho_%i_side%i.jpg", irunold, side);
2519  c2->SaveAs(out_file);
2520  }
2521  delete c2;
2522 
2523  // if (m_combined) {
2524  gStyle->SetTitleFontSize(0.045);
2525  gStyle->SetPadRightMargin(0.13);
2526  gStyle->SetPadBottomMargin(0.15);
2527  gStyle->SetPadLeftMargin(0.1);
2528  gStyle->SetOptStat(0);
2529  xsiz = 700;
2530  ysiz = 600;
2531  TCanvas *c1 = new TCanvas("c1", "Fitted const in each tower", xsiz, ysiz);
2532  const_eta_phi->GetXaxis()->SetTitle("#eta");
2533  const_eta_phi->GetXaxis()->SetTitleSize(0.065);
2534  const_eta_phi->GetXaxis()->SetTitleOffset(0.85); //6);
2535  const_eta_phi->GetXaxis()->CenterTitle();
2536  const_eta_phi->GetXaxis()->SetLabelSize(0.045);
2537  const_eta_phi->GetXaxis()->SetLabelOffset(0.01);
2538 
2539  const_eta_phi->GetYaxis()->SetTitle("#phi");
2540  const_eta_phi->GetYaxis()->SetTitleSize(0.075);
2541  const_eta_phi->GetYaxis()->SetTitleOffset(0.5);
2542  const_eta_phi->GetYaxis()->CenterTitle();
2543  const_eta_phi->GetYaxis()->SetLabelSize(0.045);
2544  const_eta_phi->GetYaxis()->SetLabelOffset(0.01);
2545 
2546  const_eta_phi->Draw("colz");
2547  sprintf(out_file, "high_hoconst_eta_phi_%i.jpg",irunold);
2548  c1->SaveAs(out_file);
2549 
2550  delete c1;
2551 
2552  for (int jk=0; jk<netamx; jk++) {
2553  int ieta = (jk<15) ? jk+1 : 14-jk;
2554  if (rms_phi[jk]>0) {
2555  mean_phi_ave->Fill(ieta, mean_phi[jk]/rms_phi[jk]);
2556  mean_phi_ave->SetBinError(mean_phi_ave->FindBin(ieta), pow(double(rms_phi[jk]), -0.5));
2557  }
2558  }
2559 
2560  for (int ij=0; ij<nphimx; ij++) {
2561  if (rms_eta[ij] >0) {
2562  mean_eta_ave->Fill(ij+1, mean_eta[ij]/rms_eta[ij]);
2563  mean_eta_ave->SetBinError(mean_eta_ave->FindBin(ij+1), pow(double(rms_eta[ij]), -0.5));
2564  }
2565  }
2566 
2567  ysiz =450;
2568  gStyle->SetPadLeftMargin(0.13);
2569  gStyle->SetPadRightMargin(0.03);
2570 
2571 
2572  TCanvas *c2y = new TCanvas("c2", "Avearge signal in eta and phi", xsiz, ysiz);
2573  c2y->Divide(2,1);
2574  mean_eta_ave->GetXaxis()->SetTitle("#phi");
2575  mean_eta_ave->GetXaxis()->SetTitleSize(0.085);
2576  mean_eta_ave->GetXaxis()->SetTitleOffset(0.65);
2577  mean_eta_ave->GetXaxis()->CenterTitle();
2578  mean_eta_ave->GetXaxis()->SetLabelSize(0.05);
2579  mean_eta_ave->GetXaxis()->SetLabelOffset(0.001);
2580 
2581  mean_eta_ave->GetYaxis()->SetTitle("Signal (GeV)/MIP");
2582  mean_eta_ave->GetYaxis()->SetTitleSize(0.055);
2583  mean_eta_ave->GetYaxis()->SetTitleOffset(1.3);
2584  mean_eta_ave->GetYaxis()->CenterTitle();
2585  mean_eta_ave->GetYaxis()->SetLabelSize(0.045);
2586  mean_eta_ave->GetYaxis()->SetLabelOffset(0.01);
2587  mean_eta_ave->SetMarkerSize(0.60);
2588  mean_eta_ave->SetMarkerColor(2);
2589  mean_eta_ave->SetMarkerStyle(20);
2590 
2591  c2y->cd(1); mean_eta_ave->Draw();
2592 
2593  mean_phi_ave->GetXaxis()->SetTitle("#eta");
2594  mean_phi_ave->GetXaxis()->SetTitleSize(0.085);
2595  mean_phi_ave->GetXaxis()->SetTitleOffset(0.65); //55);
2596  mean_phi_ave->GetXaxis()->CenterTitle();
2597  mean_phi_ave->GetXaxis()->SetLabelSize(0.05);
2598  mean_phi_ave->GetXaxis()->SetLabelOffset(0.001);
2599 
2600  mean_phi_ave->GetYaxis()->SetTitle("Signal (GeV)/MIP");
2601  mean_phi_ave->GetYaxis()->SetTitleSize(0.055);
2602  mean_phi_ave->GetYaxis()->SetTitleOffset(1.3);
2603  mean_phi_ave->GetYaxis()->CenterTitle();
2604  mean_phi_ave->GetYaxis()->SetLabelSize(0.045);
2605  mean_phi_ave->GetYaxis()->SetLabelOffset(0.01);
2606  mean_phi_ave->SetMarkerSize(0.60);
2607  mean_phi_ave->SetMarkerColor(2);
2608  mean_phi_ave->SetMarkerStyle(20);
2609 
2610  c2y->cd(2); mean_phi_ave->Draw();
2611 
2612  sprintf(out_file, "high_hoaverage_eta_phi_%i.jpg",irunold);
2613  c2y->SaveAs(out_file);
2614 
2615  delete c2y;
2616  // } else { //m_combined
2617 
2618  xsiz = 800;
2619  ysiz = 450;
2620  TCanvas *c3 = new TCanvas("c3", "Avearge signal in eta and phi", xsiz, ysiz);
2621  c3->Divide(2,1);
2622  mean_phi_hst->GetXaxis()->SetTitle("#eta");
2623  mean_phi_hst->GetXaxis()->SetTitleSize(0.065);
2624  mean_phi_hst->GetXaxis()->SetTitleOffset(0.9);
2625  mean_phi_hst->GetXaxis()->CenterTitle();
2626  mean_phi_hst->GetXaxis()->SetLabelSize(0.065);
2627  mean_phi_hst->GetXaxis()->SetLabelOffset(0.001);
2628 
2629  mean_phi_hst->GetYaxis()->SetTitle("GeV/MIP");
2630  mean_phi_hst->GetYaxis()->SetTitleSize(0.055);
2631  mean_phi_hst->GetYaxis()->SetTitleOffset(0.9);
2632  mean_phi_hst->GetYaxis()->CenterTitle();
2633  mean_phi_hst->GetYaxis()->SetLabelSize(0.065);
2634  mean_phi_hst->GetYaxis()->SetLabelOffset(0.01);
2635 
2636  mean_phi_hst->SetMarkerColor(4);
2637  mean_phi_hst->SetMarkerSize(0.8);
2638  mean_phi_hst->SetMarkerStyle(20);
2639  mean_phi_hst->Draw();
2640 
2641  sprintf(out_file, "low_mean_phi_hst_%i.jpg",irunold);
2642  c3->SaveAs(out_file);
2643 
2644  delete c3;
2645 
2646  // } //m_combined
2647 
2648 
2649  gStyle->SetOptLogy(1);
2650  gStyle->SetPadTopMargin(.1);
2651  gStyle->SetPadLeftMargin(.15);
2652  xsiz = 800;
2653  ysiz = 500;
2654  TCanvas *c0x = new TCanvas("c0x", "Signal in each ring", xsiz, ysiz);
2655 
2656  c0x->Divide(3,2);
2657  for (int ij=0; ij<ringmx; ij++) {
2658  int iread = (ij==2) ? routmx : rout12mx;
2659  if (m_digiInput) {
2660  com_sigrsg[ij][iread]->GetXaxis()->SetTitle("Signal/ped (fC)");
2661  } else {
2662  com_sigrsg[ij][iread]->GetXaxis()->SetTitle("Signal/ped (GeV)");
2663  }
2664  com_sigrsg[ij][iread]->GetXaxis()->SetTitleSize(0.060);
2665  com_sigrsg[ij][iread]->GetXaxis()->SetTitleOffset(1.05);
2666  com_sigrsg[ij][iread]->GetXaxis()->CenterTitle();
2667  com_sigrsg[ij][iread]->GetXaxis()->SetLabelSize(0.065);
2668  com_sigrsg[ij][iread]->GetXaxis()->SetLabelOffset(0.01);
2669 
2670  com_sigrsg[ij][iread]->GetYaxis()->SetLabelSize(0.065);
2671  com_sigrsg[ij][iread]->GetYaxis()->SetLabelOffset(0.01);
2672 
2673 
2674  com_sigrsg[ij][iread]->SetLineWidth(3);
2675  com_sigrsg[ij][iread]->SetLineColor(4);
2676 
2677  c0x->cd(ij+1); com_sigrsg[ij][iread]->Draw();
2678 
2679  com_crossg[ij][iread]->SetLineWidth(2);
2680  com_crossg[ij][iread]->SetLineColor(2);
2681  com_crossg[ij][iread]->Draw("same");
2682  }
2683  sprintf(out_file, "hosig_ring_%i.jpg",irunold);
2684  c0x->SaveAs(out_file);
2685  delete c0x;
2686 
2687  gStyle->SetTitleFontSize(0.06);
2688  gStyle->SetOptStat(0);
2689  gStyle->SetOptLogy(0);
2690 
2691  TCanvas *c0 = new TCanvas("c0", "Signal in each ring", xsiz, ysiz);
2692 
2693  c0->Divide(3,2);
2694  for (int jk=0; jk<ringmx; jk++) {
2695  peak_hpdrm[jk]->GetXaxis()->SetTitle("RM #");
2696  peak_hpdrm[jk]->GetXaxis()->SetTitleSize(0.070);
2697  peak_hpdrm[jk]->GetXaxis()->SetTitleOffset(1.0);
2698  peak_hpdrm[jk]->GetXaxis()->CenterTitle();
2699  peak_hpdrm[jk]->GetXaxis()->SetLabelSize(0.065);
2700  peak_hpdrm[jk]->GetXaxis()->SetLabelOffset(0.01);
2701 
2702  peak_hpdrm[jk]->GetYaxis()->SetTitle("Peak(GeV)/MIP");
2703 
2704  peak_hpdrm[jk]->GetYaxis()->SetTitleSize(0.07);
2705  peak_hpdrm[jk]->GetYaxis()->SetTitleOffset(1.3);
2706  peak_hpdrm[jk]->GetYaxis()->CenterTitle();
2707  peak_hpdrm[jk]->GetYaxis()->SetLabelSize(0.065);
2708  peak_hpdrm[jk]->GetYaxis()->SetLabelOffset(0.01);
2709  // peak_hpdrm[jk]->SetLineWidth(3);
2710  // peak_hpdrm[jk]->SetLineColor(4);
2711  peak_hpdrm[jk]->SetMarkerSize(0.60);
2712  peak_hpdrm[jk]->SetMarkerColor(2);
2713  peak_hpdrm[jk]->SetMarkerStyle(20);
2714 
2715 
2716  c0->cd(jk+1); peak_hpdrm[jk]->Draw();
2717  }
2718  sprintf(out_file, "comb_peak_hpdrm_%i.jpg",irunold);
2719  c0->SaveAs(out_file);
2720 
2721  delete c0;
2722 
2723  TCanvas *c1y = new TCanvas("c1y", "Signal in each ring", xsiz, ysiz);
2724 
2725  c1y->Divide(3,2);
2726  for (int jk=0; jk<ringmx; jk++) {
2727  const_hpdrm[jk]->GetXaxis()->SetTitle("RM #");
2728  const_hpdrm[jk]->GetXaxis()->SetTitleSize(0.070);
2729  const_hpdrm[jk]->GetXaxis()->SetTitleOffset(1.3);
2730  const_hpdrm[jk]->GetXaxis()->CenterTitle();
2731  const_hpdrm[jk]->GetXaxis()->SetLabelSize(0.065);
2732  const_hpdrm[jk]->GetXaxis()->SetLabelOffset(0.01);
2733 
2734  if (m_digiInput) {
2735  const_hpdrm[jk]->GetYaxis()->SetTitle("Peak(fC)");
2736  } else {
2737  const_hpdrm[jk]->GetYaxis()->SetTitle("Peak(GeV)");
2738  }
2739  const_hpdrm[jk]->GetYaxis()->SetTitleSize(0.065);
2740  const_hpdrm[jk]->GetYaxis()->SetTitleOffset(1.0);
2741  const_hpdrm[jk]->GetYaxis()->CenterTitle();
2742  const_hpdrm[jk]->GetYaxis()->SetLabelSize(0.065);
2743  const_hpdrm[jk]->GetYaxis()->SetLabelOffset(0.01);
2744  // const_hpdrm[jk]->SetLineWidth(3);
2745  // const_hpdrm[jk]->SetLineColor(4);
2746  const_hpdrm[jk]->SetMarkerSize(0.60);
2747  const_hpdrm[jk]->SetMarkerColor(2);
2748  const_hpdrm[jk]->SetMarkerStyle(20);
2749 
2750  c1y->cd(jk+1); const_hpdrm[jk]->Draw();
2751  }
2752 
2753  sprintf(out_file, "comb_const_hpdrm_%i.jpg",irunold);
2754  c1y->SaveAs(out_file);
2755 
2756  delete c1y;
2757 
2758  } //if (m_figure) {
2759 
2760  // ps.Close();
2761  // file_out.close();
2762 
2763  }// if (m_constant){
2764 
2765 
2766  if (m_figure) {
2767  for (int ij=0; ij<nphimx; ij++) {
2768  for (int jk=0; jk<netamx; jk++) {
2769  stat_eta[jk]->Fill(ij+1,sigrsg[jk][ij]->GetEntries());
2770  statmn_eta[jk]->Fill(ij+1,sigrsg[jk][ij]->GetMean());
2771  }
2772  }
2773 
2774  xsiz = 700;
2775  ysiz = 450;
2776  gStyle->SetTitleFontSize(0.09);
2777  gStyle->SetPadBottomMargin(0.14);
2778  gStyle->SetPadLeftMargin(0.17);
2779  gStyle->SetPadRightMargin(0.01);
2780  gStyle->SetNdivisions(303,"XY");
2781  gStyle->SetOptLogy(1);
2782 
2783  TCanvas *c2x = new TCanvas("c2x", "runfile", xsiz, ysiz);
2784  c2x->Divide(5,3);
2785  for (int side=0; side <2; side++) {
2786  int nmn = 0;
2787  int nmx = netamx/2;
2788  if (side==1) {
2789  nmn = netamx/2;
2790  nmx = netamx;
2791  }
2792  int nzone = 0;
2793  char name[200];
2794 
2795  for (int ij=nmn; ij<nmx; ij++) {
2796  int ieta = (ij<15) ? ij+1 : 14-ij;
2797  c2x->cd(nzone+1);
2798  if (m_digiInput) {
2799  sprintf(name,"fC(#eta=%i)",ieta);
2800  } else {
2801  sprintf(name,"GeV(#eta=%i)",ieta);
2802  }
2803  sigrsg[ij][nphimx]->GetXaxis()->SetTitle(name);
2804  sigrsg[ij][nphimx]->GetXaxis()->SetTitleSize(.08);
2805  sigrsg[ij][nphimx]->GetXaxis()->CenterTitle();
2806  sigrsg[ij][nphimx]->GetXaxis()->SetTitleOffset(0.90);
2807  sigrsg[ij][nphimx]->GetXaxis()->SetLabelSize(.08);
2808  sigrsg[ij][nphimx]->GetXaxis()->SetLabelOffset(.01);
2809 
2810  sigrsg[ij][nphimx]->GetYaxis()->SetLabelSize(.08);
2811  sigrsg[ij][nphimx]->GetYaxis()->SetLabelOffset(.01);
2812  sigrsg[ij][nphimx]->SetLineWidth(2);
2813  sigrsg[ij][nphimx]->SetLineColor(4);
2814  sigrsg[ij][nphimx]->Draw();
2815  crossg[ij][nphimx]->SetLineWidth(2);
2816  crossg[ij][nphimx]->SetLineColor(2);
2817  crossg[ij][nphimx]->Draw("same");
2818  nzone++;
2819  }
2820 
2821  sprintf(out_file, "sig_ho_%i_side%i.eps", irunold, side);
2822  c2x->SaveAs(out_file);
2823 
2824  sprintf(out_file, "sig_ho_%i_side%i.jpg", irunold, side);
2825  c2x->SaveAs(out_file);
2826  }
2827 
2828  gStyle->SetOptLogy(0);
2829  c2x = new TCanvas("c2x", "runfile", xsiz, ysiz);
2830  c2x->Divide(5,3);
2831  for (int side=0; side <2; side++) {
2832  int nmn = 0;
2833  int nmx = netamx/2;
2834  if (side==1) {
2835  nmn = netamx/2;
2836  nmx = netamx;
2837  }
2838  int nzone = 0;
2839 
2840  nzone = 0;
2841  for (int ij=nmn; ij<nmx; ij++) {
2842  c2x->cd(nzone+1);
2843  statmn_eta[ij]->SetLineWidth(2);
2844  statmn_eta[ij]->SetLineColor(4);
2845  statmn_eta[ij]->GetXaxis()->SetTitle("#phi index");
2846  statmn_eta[ij]->GetXaxis()->SetTitleSize(.08);
2847  statmn_eta[ij]->GetXaxis()->CenterTitle();
2848  statmn_eta[ij]->GetXaxis()->SetTitleOffset(0.9);
2849  statmn_eta[ij]->GetYaxis()->SetLabelSize(.08);
2850  statmn_eta[ij]->GetYaxis()->SetLabelOffset(.01);
2851  statmn_eta[ij]->GetXaxis()->SetLabelSize(.08);
2852  statmn_eta[ij]->GetXaxis()->SetLabelOffset(.01);
2853  if (m_digiInput) {
2854  statmn_eta[ij]->GetYaxis()->SetTitle("fC");
2855  } else {
2856  statmn_eta[ij]->GetYaxis()->SetTitle("GeV");
2857  }
2858  statmn_eta[ij]->GetYaxis()->SetTitleSize(.075);
2859  statmn_eta[ij]->GetYaxis()->CenterTitle();
2860  statmn_eta[ij]->GetYaxis()->SetTitleOffset(1.30);
2861 
2862  statmn_eta[ij]->Draw();
2863  nzone++;
2864  }
2865 
2866  sprintf(out_file, "statmnho_%i_side%i.eps", irunold, side);
2867  c2x->SaveAs(out_file);
2868 
2869  sprintf(out_file, "statmnho_%i_side%i.jpg", irunold, side);
2870  c2x->SaveAs(out_file);
2871 
2872  gStyle->SetOptLogy(1);
2873  gStyle->SetNdivisions(203,"XY");
2874 
2875  nzone = 0;
2876  for (int ij=nmn; ij<nmx; ij++) {
2877  c2x->cd(nzone+1);
2878  stat_eta[ij]->SetLineWidth(2);
2879  stat_eta[ij]->SetLineColor(4);
2880  stat_eta[ij]->GetXaxis()->SetTitle("#phi index");
2881  stat_eta[ij]->GetXaxis()->SetTitleSize(.08);
2882  stat_eta[ij]->GetXaxis()->CenterTitle();
2883  stat_eta[ij]->GetXaxis()->SetTitleOffset(0.80);
2884  stat_eta[ij]->GetXaxis()->SetLabelSize(.08);
2885  stat_eta[ij]->GetXaxis()->SetLabelOffset(.01);
2886  stat_eta[ij]->GetYaxis()->SetLabelSize(.08);
2887  stat_eta[ij]->GetYaxis()->SetLabelOffset(.01);
2888 
2889  stat_eta[ij]->Draw();
2890  nzone++;
2891  }
2892 
2893  sprintf(out_file, "statho_%i_side%i.eps", irunold, side);
2894  c2x->SaveAs(out_file);
2895 
2896  sprintf(out_file, "statho_%i_side%i.jpg", irunold, side);
2897  c2x->SaveAs(out_file);
2898  }
2899  delete c2x;
2900 
2901  } //if (m_figure) {
2902 
2903  if (!m_constant) { //m_constant
2904  for (int j=0; j<netamx; j++) {
2905  for (int i=0; i<nphimx; i++) {
2906  if (crossg[j][i]) { delete crossg[j][i];}
2907  if (sigrsg[j][i]) { delete sigrsg[j][i];}
2908  }
2909  }
2910  }
2911 }
2912 
2913 //define this as a plug-in
2915 
2916 
2917 
RunNumber_t run() const
Definition: EventID.h:42
tuple arglist
Definition: fitWZ.py:38
TH1F * corrsgrb[netamx][nphimx]
TH2F * sig_effi[neffip]
EventNumber_t event() const
Definition: EventID.h:44
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
long int flag
Definition: mlp_lapack.h:47
TH1F * const_hpdrm[ringmx]
static const int ringmx
TH1F * com_crossg[ringmx][routmx+1]
std::string theRootFileName
TH1F * corrsgr[netamx][nphimx]
virtual void beginJob()
static const int npixrigh[21]
static const int mapx2[6][3]
static const int npixleft[21]
TH1F * corrsglu[netamx][nphimx]
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
HOCalibAnalyzer(const edm::ParameterSet &)
static const int rout12mx
TH1F * crossg[netamx][nphimx+1]
static const int phimap[4][21]
#define abs(x)
Definition: mlp_lapack.h:159
TH1F * com_corrsgall[ringmx][sectmx]
float invang[netamx][nphimx+1]
TH1F * corrsgc[netamx][nphimx]
#define min(a, b)
Definition: mlp_lapack.h:161
TH1F * peak_hpdrm[ringmx]
static const int mapx1[6][3]
TH1F * corrsglb[netamx][nphimx]
TH1F * sigrsg[netamx][nphimx+1]
TH1F * com_corrsglu[ringmx][sectmx]
static const int npixriup[21]
static const int neffip
edm::InputTag hoCalibVariableCollectionTag
void set_sigma(double &x, bool mdigi)
void beginJob()
Definition: Breakpoints.cc:15
TProfile * hopedtime[netamx][nphimx]
TProfile * hbtime[netamx][nphimx]
std::string theoutputtxtFile
TH1F * com_sigrsg[ringmx][routmx+1]
TProfile * com_hbtime[ringmx][sectmx]
float com_invang[ringmx][routmx+1]
int iEvent
Definition: GenABIO.cc:243
int np
Definition: AMPTWrapper.h:33
virtual void analyze(const edm::Event &, const edm::EventSetup &)
const T & max(const T &a, const T &b)
virtual void endJob()
Double_t totalfunc(Double_t *x, Double_t *par)
T sqrt(T t)
Definition: SSEVec.h:48
vector< float > sig_reg[netamx][nphimx+1]
std::string theoutputpsFile
static const int npixribt[21]
TH1F * corrsgl[netamx][nphimx]
TH1F * com_corrsgc[ringmx][sectmx]
int iphifit
TH1F * com_corrsgrb[ringmx][sectmx]
int j
Definition: DBlmapReader.cc:9
TH1F * com_corrsgr[ringmx][sectmx]
void fcnsg(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t flag)
double f[11][100]
const int nbgpr
double fitprm[nsgpr][netamx]
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:361
TProfile * hotime[netamx][nphimx]
static const int npixleup[21]
Double_t gausX(Double_t *x, Double_t *par)
int k[5][pyjets_maxn]
const int nsgpr
TH1F * corrsgru[netamx][nphimx]
static const int routmx
TH1F * com_corrsgl[ringmx][sectmx]
static const int mapx0m[9][2]
TH1F * com_corrsglb[ringmx][sectmx]
const T & get() const
Definition: EventSetup.h:55
void set_mean(double &x, bool mdigi)
TH1F * peak_eta[netamx]
static const int sectmx
TProfile * com_hopedtime[ringmx][sectmx]
vector< float > cro_ssg[netamx][nphimx+1]
edm::EventID id() const
Definition: EventBase.h:56
const int nphimx
const int netamx
static const int mapx0p[9][2]
tuple gMinuit
Definition: fitWZ.py:35
TH1F * const_eta[netamx]
T * make() const
make new ROOT object
TH1F * com_corrsgru[ringmx][sectmx]
int ietafit
tuple cout
Definition: gather_cfg.py:121
TProfile * com_hotime[ringmx][sectmx]
void fcnbg(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t flag)
Definition: DDAxes.h:10
Double_t langaufun(Double_t *x, Double_t *par)
static const int etamap[4][21]
static const int npixlebt[21]
static const int mypow_2_ncut
TH1F * statmn_eta[netamx]
tuple size
Write out results.
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
TProfile * sigvsevt[15][ncut]
TH1F * corrsgall[netamx][nphimx]
TH1F * stat_eta[netamx]