Project

General

Profile

makeEffHist_v2.C

Version 2 - Muhammad Elnimr, 02/04/2015 12:31 PM

 
1
#define makeEffHist_cxx
2
#include "makeEffHist.h"
3
#include <TH2.h>
4
#include <TStyle.h>
5
#include <TCanvas.h>
6

    
7
void makeEffHist::Loop(int nEvents=0,char * filePrefex="")
8
{
9
//   In a ROOT session, you can do:
10
//      Root > .L makeEffHist.C
11
//      Root > makeEffHist t
12
//      Root > t.GetEntry(12); // Fill t data members with entry number 12
13
//      Root > t.Show();       // Show values of entry 12
14
//      Root > t.Show(16);     // Read and show values of entry 16
15
//      Root > t.Loop();       // Loop on all entries
16
//
17
  prettyPlots();
18
//     This is the loop skeleton where:
19
//    jentry is the global entry number in the chain
20
//    ientry is the entry number in the current Tree
21
//  Note that the argument to GetEntry must be:
22
//    jentry for TChain::GetEntry
23
//    ientry for TTree::GetEntry and TBranch::GetEntry
24
//
25
//       To read only selected branches, Insert statements like:
26
// METHOD1:
27
//    fChain->SetBranchStatus("*",0);  // disable all branches
28
//    fChain->SetBranchStatus("branchname",1);  // activate branchname
29
// METHOD2: replace line
30
//    fChain->GetEntry(jentry);       //read all branches
31
//by  b_branchname->GetEntry(ientry); //read only this branch
32
   if (fChain == 0) return;
33

    
34
  double fMinMCKE=0.05;           // Minimum MC particle kinetic energy (GeV).
35
  double fMinMCLen=20.;          // Minimum MC particle length in tpc (cm).
36
  double fMatchColinearity=0.97;  // Minimum matching colinearity.
37
  double fMatchDisp=2.0;         // Maximum matching displacement.
38
  double fWMatchDisp=15;        // Maximum matching displacement in the w direction.
39
  double pdg=-13;
40

    
41
   Long64_t nentries = fChain->GetEntriesFast();
42
   TH1F * ntracks=new TH1F("ntracks","No. of Tracks",100,0,10);
43
   TH2F * trkytrkz=new TH2F("trkytrkz","Track Spts. y vs. z",1000,-10,300,1000,-150,150);
44
   TH2F * trkytrkx=new TH2F("trkytrkx","Track Spts. y vs. x",1000,-50,250,1000,-150,150);
45
   TH2F * trkztrkx=new TH2F("trkztrkx","Track Spts. z vs. x",1000,-50,250,1000,-100,200);
46
   TH3F * trkztrkxtrky=new TH3F("xyz","xyz",10,-10,300,10,-100,200,10,-10,300);
47

    
48

    
49
   TH2F * hgoodlen=new TH2F("hgoodlen","hgoodlen",300,0,300,300,0,300);
50

    
51
   TH1F * htrkx=new TH1F("htrkx","Track Spacepoints x (cm)",100,-100,200);
52
   TH1F * htrky=new TH1F("htrky","Track Spacepoints y (cm)",100,-100,200);
53
   TH1F * htrkz=new TH1F("htrkz","Track Spacepoints z (cm)",100,-100,200);
54
   TH1F * htrkstartx=new TH1F("htrkstartx","Track Start x (cm)",100,-100,300);
55
   TH1F * htrkstarty=new TH1F("htrkstarty","Track Start y (cm)",100,-100,200);
56
   TH1F * htrkstartz=new TH1F("htrkstartz","Track Start z (cm)",100,-100,200);
57
   TH1F * htrkstartx_diff=new TH1F("htrkstartx_diff","htrkstartx_diff",100,-50,50);
58
   TH1F * htrkstarty_diff=new TH1F("htrkstarty_diff","htrkstarty_diff",100,-50,50);
59
   TH1F * htrkstartz_diff=new TH1F("htrkstartz_diff","htrkstartz_diff",100,-50,50);
60
   TH1F * htrkendx=new TH1F("htrkendx","Track End x (cm)",100,-100,300);
61
   TH1F * htrkendy=new TH1F("htrkendy","Track End y (cm)",100,-100,200);
62
   TH1F * htrkendz=new TH1F("htrkendz","Track End z (cm)",100,-100,200);
63
   TH1F * hntracks_reco=new TH1F("hntracks_reco","No. reco tracks",10,-1,4);
64

    
65
   TH1F * htrkspacepoints=new TH1F("htrkspacepoints","All spacepoints",50,0,100);
66
   TH1F * hgtrkspacepoints=new TH1F("hgtrkspacepoints","Good spacepoints",50,0,100);
67
   TH1F * hetrkspacepoints=new TH1F("hetrkspacepoints","Eff. spacepoints",50,0,100);
68

    
69
   TH1F * htrkclust=new TH1F("htrkclust","All clust",50,0,100);
70
   TH1F * hgtrkclust=new TH1F("hgtrkclust","Good clust",50,0,100);
71
   TH1F * hetrkclust=new TH1F("hetrkclust","Eff. clust",50,0,100);
72

    
73
   TH1F * htrkhits=new TH1F("htrkhits","All hits",25,0,400);
74
   TH1F * hgtrkhits=new TH1F("hgtrkhits","Good hits",25,0,400);
75
   TH1F * hetrkhits=new TH1F("hetrkhits","Eff. hits",25,0,400);
76

    
77

    
78
   TH2F * hnHitsnClust=new TH2F("hnHitsnClust","nhits vs. nclusts",300,-10,290,50,-2,48);
79
   TH2F * hnHitsnSpcpts=new TH2F("hnHitsnSpcpts","nhits vs. nspcpts",300,-10,290,100,-10,90);
80
   TH2F * hnClustSpcpts=new TH2F("hnClustSpcpts","nclusts vs. nspcpts",50,-2,48,100,-10,90);
81

    
82
   Long64_t nbytes = 0, nb = 0;
83
   MCHists mchists("somename");
84
   if(nEvents==0);
85
   else
86
     nentries=nEvents;
87
   int COUNTER=0;
88
   for (Long64_t jentry=0; jentry<nentries;jentry++) {
89
      Long64_t ientry = LoadTree(jentry);
90
      if (ientry < 0) break;
91
      nb = fChain->GetEntry(jentry);   nbytes += nb;
92
      // if (Cut(ientry) < 0) continue;
93
      Long64_t i=0;
94
      double pdg=-13;
95
      int  nTracks=0;
96
      int nspacepoints=0;
97
      int skipEvent=0;
98
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
99
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
100
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
101
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
102
      //
103
      //
104
      //  Reco track loop
105
      //
106
      //
107
      int totalSpacepoints=0;
108
      for(int j=0;j<ntracks_reco;j++)
109
        {
110
          /*          if(trkpdg[j]!=-13 )
111
            {
112
              cout << "----" << trkpdg[j] << endl;
113
              cout << " track stuff " << trkd2[j] << endl;
114
              }*/
115
          //          if(trklen_cut_MC[0]<5  && trklen[j]>60 ) std::cout << " Event No." <<  event << std::endl;
116
          if(!(1/*trkang[j]<0 &&trkd2[j]<0.01 &&trkmatchdisp[j]>0  && abs(trklen[j]-10-trklen_cut_MC[0]+10)<30 */))
117
            {
118
            skipEvent=1;
119
            break;
120
            }
121
          totalSpacepoints+=ntrkhits[j];
122
        }
123
      if(trklen_cut_MC[0]<200)
124
              continue;
125
      if(skipEvent) 
126
        continue;
127
      /*      if(ntracks_reco==0)
128
              continue;  */
129
      if(trkpdg_MC[0]!=-13)
130
        continue;
131
      
132
      //      if(trklen[j]-10-trklen_cut_MC[0]+10<30)
133
      //      if(trklen_cut_MC[0]>180 &&trklen[i]< 100 )
134
      //        cout << "Event number is :" << event << endl;
135

    
136
      /*      mchists.fHmcstartx->Fill(trkstartx_MC[i]);
137
      mchists.fHmcstarty->Fill(trkstarty_MC[i]);
138
      mchists.fHmcstartz->Fill(trkstartz_MC[i]);
139
      mchists.fHmcendx->Fill(trkendx_MC[i]);
140
      mchists.fHmcendy->Fill(trkendy_MC[i]);
141
      mchists.fHmcendz->Fill(trkendz_MC[i]);*/
142
      mchists.fHmctheta->Fill(trktheta_MC[i]);
143
      mchists.fHmcphi->Fill(trkphi_MC[i]);
144
      mchists.fHmctheta_xz->Fill(trktheta_xz_MC[i]);
145
      mchists.fHmctheta_yz->Fill(trktheta_yz_MC[i]);
146
      mchists.fHmcmom->Fill(trkmom_MC[i]);
147
      mchists.fHmclen->Fill(trklen_cut_MC[i]);
148
      mchists.fHmclenshort->Fill(trklen_cut_MC[i]);
149
      hntracks_reco->Fill(ntracks_reco);
150
      htrkspacepoints->Fill(totalSpacepoints);
151
      htrkclust->Fill(nclust);
152

    
153
      int nhits_trk=0;
154
      /*      for(int p=0;p<nhits;p++)
155
        {
156
          if(hit_trkid[p]>0)
157
            nhits_trk++;
158
            }*/
159
      //      htrkhits->Fill(nhits_trk);
160
      htrkhits->Fill(nhits);
161

    
162
      COUNTER++;
163
      
164

    
165
      if(ntracks_reco==0 )
166
        {
167
          if(trklen_cut_MC[0]>200)
168
            {
169
              cout << "=============================================" << endl; 
170
              cout << " Run : " << run << "Subrun :  "  << subrun << " Event : "  << event <<endl;
171
              cout << "Zero Tracks Reconstructed" << endl;
172
              cout << "=============================================" << endl; 
173
            }
174
        }
175
      else
176
        {
177
          for(int j=0;j<ntracks_reco;j++)
178
            {
179
              if(!( trklen[j]/trklen_cut_MC[0] > 0.5))
180
                {
181
                  if(trklen_cut_MC[0]>200)
182
                    {
183
                      cout << "=============================================" << endl; 
184
                      cout << " Run : " << run << "Subrun :  "  << subrun << " Event : "  << event <<endl;
185
                      cout << "Bad event : Reco length " << trklen[j] <<" MC length (within readoutwindow) " << trklen_cut_MC[0]<<endl;
186
                      cout << "=============================================" << endl; 
187
                    }
188
                }
189
            }
190
        }
191
      
192

    
193

    
194

    
195

    
196
      for(int j=0;j<ntracks_reco;j++)
197
        {
198
          /*          if( 0 
199
             //trkd2[j]>0.2
200
               || trklen_cut_MC[i]<fMinMCLen
201
               || trkcolinearity[j]<fMatchColinearity 
202
              // ||trklen[j]/trklen_cut_MC[i] < 0.5
203
               || abs(trkwmatchdisp[j]) <= fWMatchDisp &&  trkmatchdisp[j] <= fMatchDisp 
204
               ||ntrkhits[j]==0
205
               )
206
               continue;*/
207
          //          if(abs(trklen[j]-trklen_cut_MC[j]+100)>10 ) continue;
208
          //if(abs(trklen[j]-trklen_cut_MC[j])>10 ) continue;
209
          //          if(ntrkhits[j]==0) continue;
210
          //          cout << "Event no. is : " << event << endl;
211
          /*          htrkstartx->Fill(trkstartx[j]);
212
          htrkstarty->Fill(trkstarty[j]);
213
          htrkstartz->Fill(trkstartz[j]);*/
214

    
215
        
216
          htrkendx->Fill(trkendx[j]);
217
          htrkendy->Fill(trkendy[j]);
218
          htrkendz->Fill(trkendz[j]);
219

    
220
          
221

    
222
          htrkstartx->Fill(trkstartx[j]);
223
          htrkstarty->Fill(trkstarty[j]);
224
          htrkstartz->Fill(trkstartz[j]);
225

    
226
          htrkstartx_diff->Fill(trkendx[j]-trkstartx[j]);
227
          htrkstarty_diff->Fill(trkendy[j]-trkstarty[j]);
228
          htrkstartz_diff->Fill(trkendz[j]-trkstartz[j]);
229

    
230
          mchists.fHstartdx->Fill(trkstartx_MC[i]-trkstartx[j]);
231
          mchists.fHstartdy->Fill(trkstarty_MC[i]-trkstarty[j]);
232
          mchists.fHstartdz->Fill(trkstartz_MC[i]-trkstartz[j]);
233
          for(int tt=0;tt<ntrkhits[j];tt++)
234
            {
235
              
236
              trkytrkz->Fill(trkz[j][tt],trky[j][tt]);
237
              trkytrkx->Fill(trkx[j][tt],trky[j][tt]);
238
              trkztrkx->Fill(trkx[j][tt],trkz[j][tt]);
239
              trkztrkxtrky->Fill(trkz[j][tt],trky[j][tt],trkx[j][tt]);
240
              
241
              htrkx->Fill(trkx[j][tt]);
242
              htrky->Fill(trky[j][tt]);
243
              htrkz->Fill(trkz[j][tt]);
244
              
245
            }
246
          
247
          nspacepoints+=ntrkhits[j];
248
          nTracks++;
249
          mchists.fHtrkd2->Fill(trkd2[j]);
250
          //          if( trkmatchdisp[j] <= fMatchDisp  && trklen_cut_MC[i]> -999)
251
          //          if(  trklen_cut_MC[i]> -999)
252
          {
253
            mchists.fHlvsl->Fill(trklen_cut_MC[i],trklen[j]);
254
            mchists.fHlvslvstrkang->Fill(trkang[j],trklen_cut_MC[i],trklen[j]);
255
            mchists.fHdl->Fill(trklen[j]-trklen_cut_MC[i]);
256
          }
257

    
258
          bool good = (/*trkang[j]>0 &&ntracks_reco==1&& trklen_cut_MC[i]>10 &&abs(trkwmatchdisp[j]) <= fWMatchDisp &&  trkmatchdisp[j] <= fMatchDisp &&*/
259
                       // trklenratio[j] > 0.5 &&
260
                        trklen[j]/trklen_cut_MC[i] > 0.5
261
                       //                                && pdg==-13 
262
                        //                        && trklen_cut_MC[i]>fMinMCLen
263
                        //&& trkcolinearity[j]>fMatchColinearity
264
                        ) ;
265
          //          bool good=1;
266
          //trklen[j]/trklen_cut_MC[i] > 0.5 && pdg==13 && trklen_cut_MC[i]>fMinMCLen  && trkcolinearity[j]>fMatchColinearity);
267
          if(good) {
268
            /* mchists.fHgstartx->Fill(mcstart.X());
269
               mchists.fHgstarty->Fill(mcstart.Y());
270
            mchists.fHgstartz->Fill(mcstart.Z());
271
            mchists.fHgendx->Fill(mcend.X());
272
            mchists.fHgendy->Fill(mcend.Y());
273
            mchists.fHgendz->Fill(mcend.Z());*/
274
            mchists.fHgtheta->Fill(trktheta_MC[j]);
275
            mchists.fHgphi->Fill(trkphi_MC[j]);
276
            mchists.fHgtheta_xz->Fill(trktheta_xz_MC[j]);
277
            mchists.fHgtheta_yz->Fill(trktheta_yz_MC[j]);
278
            mchists.fHgmom->Fill(trkmom_MC[j]);
279
            hgoodlen->Fill(trklen_cut_MC[i],trklen[j]);
280
            hgtrkspacepoints->Fill(totalSpacepoints);
281
            mchists.fHglen->Fill(trklen_cut_MC[j]);
282
            mchists.fHglenshort->Fill(trklen_cut_MC[j]);
283
            //            hgtrkhits->Fill(nhits_trk);
284
            hgtrkhits->Fill(nhits);
285
            hgtrkclust->Fill(nclust);
286
          }
287
          else{
288
           
289
          }
290
        }
291

    
292
      //      hnHitsnClust->Fill((double)nhits_trk,(double)nclust);
293
      //      hnHitsnSpcpts->Fill((double)nhits_trk,(double)totalSpacepoints);
294
      hnHitsnClust->Fill((double)nhits,(double)nclust);
295
      hnHitsnSpcpts->Fill((double)nhits,(double)totalSpacepoints);
296
      hnClustSpcpts->Fill((double)totalSpacepoints,(double)nclust);
297
      //
298
      //
299
      // End of  Reco track loop
300
      //
301
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
302
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
303
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
304
      //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
305
      //
306
      ntracks->Fill(nTracks);
307
      //      cout << " event no " << jentry<< " no of spacepoints " << nspacepoints << endl;
308
   }
309

    
310
   cout <<"+++++++++++ "<< endl;
311
   cout << "::effcalc" << endl;
312
   cout <<"+++++++++++ "<< endl;
313
   effcalc(mchists.fHgtheta, mchists.fHmctheta, mchists.fHetheta);
314
   cout <<"+++++++++++ "<< endl;
315
   effcalc(mchists.fHgphi, mchists.fHmcphi, mchists.fHephi);
316
   cout <<"+++++++++++ "<< endl;
317
   effcalc(mchists.fHgtheta_xz, mchists.fHmctheta_xz, mchists.fHetheta_xz);
318
   cout <<"+++++++++++ "<< endl;
319
   effcalc(mchists.fHgtheta_yz, mchists.fHmctheta_yz, mchists.fHetheta_yz);
320
   cout <<"+++++++++++ "<< endl;
321
   effcalc(mchists.fHgmom, mchists.fHmcmom, mchists.fHemom);
322
   cout <<"+++++++++++ "<< endl;
323
   effcalc(mchists.fHglen, mchists.fHmclen, mchists.fHelen);
324
   cout <<"+++++++++++ "<< endl;
325
   effcalc(mchists.fHglenshort, mchists.fHmclenshort, mchists.fHelenshort);
326
   cout <<"+++++++++++ "<< endl;
327
   effcalc( hgtrkspacepoints,  htrkspacepoints, hetrkspacepoints);
328
   cout <<"+++++++++++ "<< endl;
329
   effcalc( hgtrkhits,  htrkhits, hetrkhits);
330
   cout <<"+++++++++++ "<< endl;
331
   effcalc( hgtrkclust,  htrkclust, hetrkclust);
332

    
333

    
334
   //   cout << COUNTER << endl;
335
   //
336
   //
337
   // theta
338
   //
339
     
340
   /*   new TCanvas;
341
   mchists.fHmctheta->Draw();
342
   new TCanvas;
343
   mchists.fHgtheta->Draw();*/
344
  
345
   cout<< "+++++++++++" << endl;
346
   cout << " Phi " << endl;
347
   cout<< "+++++++++++" << endl;
348
    new TCanvas;
349
   mchists.fHmcphi->Draw();
350
   mchists.fHmcphi->SetFillColor(kPink);
351
   //   new TCanvas;
352
   mchists.fHgphi->Draw("same");
353
   mchists.fHgphi->SetFillColor(kBlue);
354
   cout<< "+++++++++++" << endl;
355
   cout << " theta " << endl;
356
   cout<< "+++++++++++" << endl;
357
   new TCanvas;
358
   mchists.fHmctheta->Draw();
359
   mchists.fHmctheta->SetFillColor(kPink);
360
   //   new TCanvas;
361
   mchists.fHgtheta->Draw("same");
362
   mchists.fHgtheta->SetFillColor(kBlue);
363
   cout<< "+++++++++++" << endl;
364
   cout << " Momentum " << endl;
365
   cout<< "+++++++++++" << endl;
366
   new TCanvas;
367
   mchists.fHmcmom->Draw();
368
   mchists.fHmcmom->SetFillColor(kPink);
369
   //  new TCanvas;
370
   mchists.fHgmom->Draw("same");
371
   mchists.fHgmom->SetFillColor(kBlue);
372
   cout<< "+++++++++++" << endl;
373
   cout << " Length " << endl;
374
   cout<< "+++++++++++" << endl;
375
   new TCanvas;
376
   mchists.fHmclen->Draw();
377
   mchists.fHmclen->SetMaximum(700);
378
   mchists.fHmclen->SetMinimum(0);
379
   mchists.fHmclen->SetFillColor(kPink);
380
   //  new TCanvas;
381
   mchists.fHglen->Draw("same");
382
   mchists.fHglen->SetFillColor(kBlue);
383
   //   mchists.fHglen->SetFillStyle(3445);
384
   /////////
385
   /////////
386
   /////////
387
   /////////
388
   /////////
389

    
390
  TString tmpName="display"; 
391
  tmpName+="_hist.root";
392
  TFile * file=new TFile(tmpName,"UPDATE");
393

    
394
   TString tmpSt="";
395
   new TCanvas;
396
   mchists.fHephi->Draw("p");
397
   tmpSt=  mchists.fHephi->GetName();
398
   tmpSt+=filePrefex;
399
   mchists.fHephi->SetName(tmpSt);
400
   mchists.fHephi->Write();
401
   new TCanvas;
402
   mchists.fHetheta->Draw("p");
403
   tmpSt=  mchists.fHetheta->GetName();
404
   tmpSt+=filePrefex;
405
   mchists.fHetheta->SetName(tmpSt);
406
   mchists.fHetheta->Write();
407
   new TCanvas;
408
   mchists.fHetheta_xz->Draw("p");
409
    tmpSt=  mchists.fHetheta_xz->GetName();
410
   tmpSt+=filePrefex;
411
   mchists.fHetheta_xz->SetName(tmpSt);
412
   mchists.fHetheta_xz->Write();
413
   new TCanvas;
414
   mchists.fHetheta_yz->Draw("p");
415
  tmpSt=  mchists.fHetheta_yz->GetName();
416
   tmpSt+=filePrefex;
417
   mchists.fHetheta_yz->SetName(tmpSt);
418
   mchists.fHetheta_yz->Write();
419
   new TCanvas;
420
   mchists.fHelen->Draw("p");
421
     tmpSt=  mchists.fHelen->GetName();
422
   tmpSt+=filePrefex;
423
   mchists.fHelen->SetName(tmpSt);
424
   mchists.fHelen->Write();
425
   new TCanvas;
426
   mchists.fHelenshort->Draw("p");
427
   tmpSt=  mchists.fHelenshort->GetName();
428
   tmpSt+=filePrefex;
429
   mchists.fHelenshort->SetName(tmpSt);
430
   mchists.fHelenshort->Write();
431

    
432
   new TCanvas;
433
   mchists.fHlvsl->Draw();
434
  tmpSt=  mchists.fHlvsl->GetName();
435
   tmpSt+=filePrefex;
436
   mchists.fHlvsl->SetName(tmpSt);
437
   mchists.fHlvsl->Write();
438
   //c1->SaveAs("lvsl.png");
439
   new TCanvas;
440
   mchists.fHdl->Draw();
441
  tmpSt=  mchists.fHdl->GetName();
442
   tmpSt+=filePrefex;
443
   mchists.fHdl->SetName(tmpSt);
444
   mchists.fHdl->Write();
445
   //c1_n2->SaveAs("dl.png");
446
   new TCanvas;
447
   mchists.fHlvslvstrkang->Draw();
448
   new TCanvas;
449
   mchists.fHemom->Draw();
450
   tmpSt=  mchists.fHemom->GetName();
451
   tmpSt+=filePrefex;
452
   mchists.fHemom->SetName(tmpSt);
453
   mchists.fHemom->Write();
454

    
455
   
456
   new TCanvas;
457
   mchists.fHmctheta_xz->Draw();
458
   mchists.fHgtheta_xz->Draw("same");
459
   mchists.fHgtheta_xz->SetFillColor(kGray);
460
   mchists.fHgtheta_xz->SetMinimum(0);
461
   new TCanvas;
462
   mchists.fHmctheta_yz->Draw();
463
   mchists.fHgtheta_yz->Draw("same");
464
   mchists.fHgtheta_yz->SetFillColor(kGray);
465
   
466
   new TCanvas;
467
   hgoodlen->Draw("zcol");
468
   
469
   new TCanvas;
470
  htrkstartx_diff->Draw();
471
  htrkstartx_diff->SetFillColor(kGray);
472
  
473
  new TCanvas;
474
  htrkstarty_diff->Draw();
475
  htrkstarty_diff->SetFillColor(kGray);
476
  
477
  new TCanvas;
478
  htrkstartz_diff->Draw();
479
  htrkstartz_diff->SetFillColor(kGray);
480
  
481
  new TCanvas;
482
  hetrkspacepoints->Draw("p");
483
  tmpSt=  hetrkspacepoints->GetName();
484
  tmpSt+=filePrefex;
485
  hetrkspacepoints->SetName(tmpSt);
486
  hetrkspacepoints->Write();
487

    
488
  new TCanvas;
489
  hetrkhits->Draw("p");
490
  tmpSt=  hetrkhits->GetName();
491
  tmpSt+=filePrefex;
492
  hetrkhits->SetName(tmpSt);
493
  hetrkhits->Write();
494

    
495
 new TCanvas;
496
  hetrkclust->Draw("p");
497
  tmpSt=  hetrkclust->GetName();
498
  tmpSt+=filePrefex;
499
  hetrkclust->SetName(tmpSt);
500
  hetrkclust->Write();
501

    
502
  new TCanvas;
503
  hnHitsnClust->Draw("zcol");
504
  hnHitsnClust->Write();
505
  new TCanvas;
506
  hnHitsnSpcpts->Draw("zcol");
507
  hnHitsnSpcpts->Write();
508
  new TCanvas;
509
  hnClustSpcpts->Draw("zcol");
510
  hnClustSpcpts->Write();
511
  file->Close();
512
  //  hgtrkspacepoints->Draw("same");
513
  //htrkspacepoints->Draw();
514
  
515
  //  new TCanvas;
516
  //  mchists.htrkd2->Draw();
517
  //  new TCanvas
518
  //   saveHist(mchists.fHtrkd2);
519
  // new TCanvas;
520
  //   mchists.fHstartdz->Draw();
521
  // saveHist( mchists.fHdl, mchists.fHlvsl);
522
  //  saveHist(filePrefex,trkytrkz,trkytrkx,trkztrkx,htrkstartx,htrkstarty,htrkstartz,mchists.fHdl,mchists.fHlvsl,htrkendx,htrkendy,htrkendz,hntracks_reco);
523
  /*            mchists.fHmcmom->Draw();*/
524

    
525
}
526

    
527
void effcalc(const TH1* hnum, const TH1* hden, TH1* heff)
528
{
529
  int nbins = hnum->GetNbinsX();
530
  if (nbins != hden->GetNbinsX())
531
    cout << "something is wrong" << endl;
532
  if (nbins != heff->GetNbinsX())
533
    cout << "something is wrong" << endl;
534
  // Loop over bins, including underflow and overflow.
535
  cout << " Hist Name:    "<< heff->GetName()  << endl;
536
  for(int ibin = 0; ibin <= nbins+1; ++ibin) {
537
    double num = hnum->GetBinContent(ibin);
538
    double den = hden->GetBinContent(ibin);
539
    if(den == 0.) {
540
      heff->SetBinContent(ibin, 0.);
541
      heff->SetBinError(ibin, 0.);
542
    }
543
    else {
544
      double eff = num / den;
545
      if(eff < 0.)
546
        eff = 0.;
547
      /*      if(eff > 1.)
548
              eff = 1.;*/
549
      double err = std::sqrt(eff * (1.-eff) / den);
550
      //double err = std::sqrt(num*num+den*den);
551
      heff->SetBinContent(ibin, eff);
552
      heff->SetBinError(ibin, err);
553
    }
554
    if( heff->GetBinContent(ibin)>0.8)
555
      {
556
        cout << " Bin No" << ibin << endl;
557

    
558
        cout << "+++ Before" <<  endl;
559
        cout << " Bin Content" << heff->GetBinContent(ibin-3) << endl;
560
        cout << " Bin Error" << heff->GetBinError(ibin-3) << endl;
561
        cout << " num" << hnum->GetBinContent(ibin-3) << endl;
562
        cout << " den" << hden->GetBinContent(ibin-3) << endl;
563

    
564
        cout << " Bin Content" << heff->GetBinContent(ibin-2) << endl;
565
        cout << " Bin Error" << heff->GetBinError(ibin-2) << endl;
566
        cout << " num" << hnum->GetBinContent(ibin-2) << endl;
567
        cout << " den" << hden->GetBinContent(ibin-2) << endl;
568

    
569
        cout << " Bin Content" << heff->GetBinContent(ibin-1) << endl;
570
        cout << " Bin Error" << heff->GetBinError(ibin-1) << endl;
571
        cout << " num" << hnum->GetBinContent(ibin-1) << endl;
572
        cout << " den" << hden->GetBinContent(ibin-1) << endl;
573

    
574
        cout << "++++++" << endl;
575
        cout << " Bin Content" << heff->GetBinContent(ibin) << endl;
576
        cout << " Bin Error" << heff->GetBinError(ibin) << endl;
577
        cout << " num" << hnum->GetBinContent(ibin) << endl;
578
        cout << " den" << hden->GetBinContent(ibin) << endl;
579
      }
580

    
581
  }
582
  heff->SetMinimum(0.);
583
  heff->SetMaximum(1.05);
584
  heff->SetMarkerStyle(20);
585
}
586

    
587
void saveHist(char * filePrefex,TH2F * h1,TH2F *h2,TH2F* h3,TH1F * h4, TH1F * h5, TH1F * h6,TH1F * h7, TH2F *h8,TH1F * h9, TH1F * h10, TH1F * h11, TH1F * h12)
588
{
589
 
590

    
591

    
592
  TCanvas * tmpCan=new TCanvas("tmpCan","tmpCan",1200,800);
593
  TString tmpSt="";
594
  tmpCan->Divide(4,3);
595
  tmpCan->cd(1);
596
  h1->Draw();
597
  tmpSt=h1->GetName();
598
  tmpSt+=filePrefex;
599
  h1->SetName(tmpSt);
600
  h1->Write();
601
  h1->SetFillColor(kGray);
602
  tmpCan->cd(2);
603
  h2->Draw();
604
  tmpSt=h2->GetName();
605
  tmpSt+=filePrefex;
606
  h2->SetName(tmpSt);
607
  h2->SetFillColor(kGray);
608
  h2->Write();
609
  tmpCan->cd(3);
610
  h3->Draw();
611
 tmpSt=h3->GetName();
612
  tmpSt+=filePrefex;
613
  h3->SetName(tmpSt);
614
  h3->SetFillColor(kGray);
615
  h3->Write();
616
  tmpCan->cd(4);
617
  h4->Draw();
618
 tmpSt=h4->GetName();
619
  tmpSt+=filePrefex;
620
  h4->SetName(tmpSt);
621
  h4->SetFillColor(kGray);
622
  h4->Write();
623
  tmpCan->cd(5);
624
  h5->Draw();
625
  h5->SetFillColor(kGray);
626
  h5->Write();
627
  tmpCan->cd(6);
628
  h6->Draw();
629
 tmpSt=h6->GetName();
630
  tmpSt+=filePrefex;
631
  h6->SetName(tmpSt);
632
  h6->SetFillColor(kGray);
633
  h6->Write();
634
  tmpCan->cd(7);
635
  h7->Draw();
636
 tmpSt=h7->GetName();
637
  tmpSt+=filePrefex;
638
  h7->SetName(tmpSt);
639
  h7->SetFillColor(kGray);
640
  h7->Write();
641
  tmpCan->cd(8);
642
  h8->Draw("zcol");
643
 tmpSt=h8->GetName();
644
  tmpSt+=filePrefex;
645
  h8->SetName(tmpSt);
646
  h8->SetFillColor(kGray);
647
  h8->Write();
648
  tmpCan->cd(9);
649
  h9->Draw();
650
 tmpSt=h9->GetName();
651
  tmpSt+=filePrefex;
652
  h9->SetName(tmpSt);
653
  h9->SetFillColor(kGray);
654
  h9->Write();
655
  tmpCan->cd(10);
656
  h10->Draw();
657
   tmpSt=h10->GetName();
658
  tmpSt+=filePrefex;
659
  h10->SetName(tmpSt);
660
  h10->SetFillColor(kGray);
661
  h10->Write();
662
  tmpCan->cd(11);
663
  h11->Draw();
664
   tmpSt=h11->GetName();
665
  tmpSt+=filePrefex;
666
  h11->SetName(tmpSt);
667
  h11->SetFillColor(kGray);
668
  h11->Write();
669
  tmpCan->cd(12);
670
  h12->Draw();
671
 tmpSt=h12->GetName();
672
  tmpSt+=filePrefex;
673
  h12->SetName(tmpSt);
674
  h12->SetFillColor(kGray);
675
  h12->Write();
676
  tmpCan->Draw();
677

    
678

    
679
  //  tmpCan->SaveAs(histName);
680
  //  delete tmpCan;
681
  //  tmpCan->SaveAs(histName);
682
  //  cout << histName << endl;
683
}
684
void saveHist(TH1F * h1, TH2F * h2)
685
{
686
  TCanvas * tmpCan=new TCanvas("tmpCan","tmpCan",600,800);
687
  tmpCan->Divide(1,2);
688
  tmpCan->cd(1);
689
  h1->Draw();
690
  tmpCan->cd(2);
691
  h2->Draw();
692
  tmpCan->Draw();
693
  //  tmpCan->SaveAs(histName);
694
  //  delete tmpCan;
695
  //  tmpCan->SaveAs(histName);
696
  //  cout << histName << endl;
697
}
698

    
699

    
700
void prettyPlots()
701
{
702
  //    gROOT  -> SetStyle("Default");
703
  //    gROOT  -> SetStyle("Plain");
704
  gStyle->SetNdivisions(110,"xyz");
705
  gStyle -> SetStatFont(42);
706
  //  gStyle -> SetOptStat(111);
707
  gStyle -> SetOptStat(111);
708
  gStyle -> SetOptFit(1111);
709
  gROOT  ->ForceStyle();
710
  gStyle -> SetPalette(1);
711
  gStyle->SetStatX(0.8201149);
712
  gStyle->SetStatY(0.8444915);
713
  //  gStyle->SetStatW(0.8807471);
714
  //  gStyle->SetStatH(0.8644068);
715

    
716
  gStyle -> SetLabelFont(42,"XYZ");
717
  gStyle -> SetTitleFont(42, "XYZ");
718
  gStyle -> SetLabelSize(0.04, "XYZ");
719
  gStyle -> SetTitleSize(0.05, "XYZ");
720
  gStyle -> SetTitleOffset(1.0, "X");
721
  gStyle -> SetTitleOffset(1.0, "Y");
722
  gStyle -> SetTitleOffset(1.0, "Z");
723
}