PWGGA/PHOSTasks/PHOS_PbPbQA/macros/MakeMmixPi0.C

   1 const Double_t kMean=0.135 ; //Approximate peak position to facilitate error estimate
   2
   3 //-----------------------------------------------------------------------------
   4 void MakeMmixPi0(const Int_t centrality=0,
   5                  const char* cModule="")
   6 {
   7   //---------------------------------------------------------------------------
   8   // This macro processes PWGPP QA output of the analysis task PHOSPbPbQA
   9   // (see analysis code in the class AliAnalysisTaskPHOSPbPbQA).
  10   // It fits Real/Mixed ratio, normalize Mixed and subtract it from Real
  11   // Arguments:
  12   // centrality: 0 or 1 for centralities 0-20% or 20-100%
  13   // cModule: string with the PHOS module to analyze: "", "SM1", "SM2" or "SM3"
  14   //-
  15   // Yuri Kharlov. 19.11.2011
  16   //---------------------------------------------------------------------------
  17 /* $Id$ */
  18
  19   TFile * f = new TFile("PHOSPbPb_all.root") ;
  20   TList *histESD = (TList*) f->Get("PHOSPbPbQAResults");
  21   char key[125] ;
  22
  23   const char* pid="All";
  24   const char* txtCent[] = {"0-20%","20-100%"};
  25   TH2F * hCentrality  = (TH2F*)f->Get("hCenPHOS") ;
  26   TH1D * hCentrality1 = hCentrality->ProjectionX();
  27   TString inputKey;
  28   TString outputKey = Form("%s10_cent%d",pid,centrality);
  29
  30   TH2F *h , *hAdd;
  31   TH2F *hm, *hmAdd;
  32   if (centrality == 0 || centrality == 1) {// centrality 0-20%, 20-100%
  33     inputKey = Form("hPi0%s%s_cen%d"  ,pid,cModule,centrality);
  34     TH2F *h = (TH2F*)f->Get(inputKey) ;
  35     inputKey = Form("hMiPi0%s%s_cen%d",pid,cModule,centrality);
  36     TH2F *hm= (TH2F*)f->Get(inputKey) ;
  37     if (h==0) {
  38       printf("Missing histogram %s\n",inputKey);
  39       return;
  40     }
  41   }
  42   else {
  43     printf("Wrong centrality %d. Allowed values are 0,1,2,3,10.\n",centrality);
  44     return;
  45   }
  46
  47   // Int_t nPadX=2,nPadY=2;
  48   // Int_t nbin=4 ;
  49   // Double_t xa[]={0.8, 2.0, 3.5, 6.0, 20.} ;
  50
  51   Int_t nPadX=1,nPadY=1;
  52   Int_t nbin=1 ;
  53   Double_t xa[]={2.0, 20.} ;
  54
  55   PPRstyle();
  56   gStyle->SetOptFit(111);
  57   gStyle->SetPadLeftMargin(0.14);
  58   gStyle->SetPadRightMargin(0.01);
  59   gStyle->SetPadTopMargin(0.05);
  60   gStyle->SetPadBottomMargin(0.08);
  61
  62   TString txtModule;
  63   if      (strcmp(cModule,"")   ==0) txtModule="All PHOS modules";
  64   else if (strcmp(cModule,"SM1")==0) txtModule="PHOS module 4";
  65   else if (strcmp(cModule,"SM2")==0) txtModule="PHOS module 3";
  66   else if (strcmp(cModule,"SM3")==0) txtModule="PHOS module 2";
  67
  68   TPaveText *label = new TPaveText(0.15,0.80,0.40,0.90,"NDC");
  69   label->SetFillColor(kWhite);
  70   label->SetBorderSize(1);
  71   label->AddText(Form("Centrality %s" ,txtCent[centrality]));
  72   label->AddText(Form("%s",txtModule.Data()));
  73
  74   //Fit real only
  75   //Linear Bg
  76   char kkey[55];
  77   sprintf(kkey,outputKey.Data()) ;
  78   char key2[155];
  79   sprintf(key,"Mix%s",kkey) ;
  80   sprintf(key2,"%s_mr1",key) ;
  81   TH1D * mr1 = new TH1D(key2,"Mass",nbin,xa) ;
  82   sprintf(key2,"%s_sr1",key) ;
  83   TH1D * sr1 = new TH1D(key2,"Width",nbin,xa) ;
  84   sprintf(key2,"%s_ar1",key) ;
  85   TH1D * ar1 = new TH1D(key2,"a",nbin,xa) ;
  86   sprintf(key2,"%s_br1",key) ;
  87   TH1D * br1 = new TH1D(key2,"a",nbin,xa) ;
  88   sprintf(key2,"%s_yr1",key) ;
  89   TH1D * nr1 = new TH1D(key2,"Raw yield",nbin,xa) ;
  90   sprintf(key2,"%s_yr1int",key) ;
  91   TH1D * nr1int = new TH1D(key2,"Raw yield, integrated",nbin,xa) ;
  92
  93   //Quadratic Bg
  94   sprintf(key2,"%s_mr2",key) ;
  95   TH1D * mr2 = new TH1D(key2,"Mass",nbin,xa) ;
  96   sprintf(key2,"%s_sr2",key) ;
  97   TH1D * sr2 = new TH1D(key2,"Width",nbin,xa) ;
  98   sprintf(key2,"%s_ar2",key) ;
  99   TH1D * ar2 = new TH1D(key2,"a",nbin,xa) ;
 100   sprintf(key2,"%s_br2",key) ;
 101   TH1D * br2 = new TH1D(key2,"a",nbin,xa) ;
 102   sprintf(key2,"%s_cr2",key) ;
 103   TH1D * cr2 = new TH1D(key2,"a",nbin,xa) ;
 104   sprintf(key2,"%s_yr2",key) ;
 105   TH1D * nr2 = new TH1D(key2,"Raw yield",nbin,xa) ;
 106   sprintf(key2,"%s_yr2int",key) ;
 107   TH1D * nr2int = new TH1D(key2,"Raw yield, integrated",nbin,xa) ;
 108
 109   TF1 * fit1 = new TF1("fit",CB,0.,1.,6) ;
 110   fit1->SetParName(0,"A") ;
 111   fit1->SetParName(1,"m_{0}") ;
 112   fit1->SetParName(2,"#sigma") ;
 113   fit1->SetParName(3,"a_{0}") ;
 114   fit1->SetParName(4,"a_{1}") ;
 115   fit1->SetParName(5,"a_{2}") ;
 116   fit1->SetLineWidth(2) ;
 117   fit1->SetLineColor(2) ;
 118   TF1 * fgs = new TF1("gs",CBs,0.,1.,3) ;
 119   fgs->SetParName(0,"A") ;
 120   fgs->SetParName(1,"m_{0}") ;
 121   fgs->SetParName(2,"#sigma") ;
 122   fgs->SetLineColor(2) ;
 123   fgs->SetLineWidth(2) ;
 124
 125   TF1 * fit2 = new TF1("fit2",CB2,0.,1.,7) ;
 126   fit2->SetParName(0,"A") ;
 127   fit2->SetParName(1,"m_{0}") ;
 128   fit2->SetParName(2,"#sigma") ;
 129   fit2->SetParName(3,"a_{0}") ;
 130   fit2->SetParName(4,"a_{1}") ;
 131   fit2->SetParName(5,"a_{2}") ;
 132   fit2->SetParName(6,"a_{3}") ;
 133   fit2->SetLineWidth(2) ;
 134   fit2->SetLineColor(4) ;
 135   fit2->SetLineStyle(2) ;
 136
 137   TF1 * fbg1 = new TF1("bg1",BG1,0.,1.,3) ;
 138   TF1 * fbg2 = new TF1("bg2",BG2,0.,1.,4) ;
 139
 140   TCanvas * c3 = new TCanvas("mggFit1_Signal","mggFitCB",10,10,1200,800) ;
 141   c3->Divide(nPadX,nPadY) ;
 142
 143   TCanvas * c1 = new TCanvas("mggFit1","mggFit1",10,10,1200,800) ;
 144   c1->Divide(nPadX,nPadY) ;
 145   c1->cd(0) ;
 146
 147   TCanvas *c4, *c2;
 148   TAxis * pta=h->GetYaxis() ;
 149   TAxis * ma=h->GetXaxis() ;
 150   for(Int_t i=1;i<=nbin;i++){
 151     if(i<=15)
 152        c1->cd(i) ;
 153     else{
 154       if(c2==0){
 155         c2 = new TCanvas("mggFit2","mggFit2",10,10,1200,800) ;
 156         c2->Divide(nPadX,nPadY) ;
 157       }
 158       c2->cd(i-16) ;
 159     }
 160     printf("\t%.1f<pt<%.1f GeV/c\n",xa[i-1],xa[i]);
 161     Int_t imin=pta->FindBin(xa[i-1]+0.0001);
 162     Int_t imax=pta->FindBin(xa[i]-0.0001) ;
 163     Double_t pt=(xa[i]+xa[i-1])/2. ;
 164     TH1D * hp = h->ProjectionX("re",imin,imax) ;
 165     hp->Sumw2() ;
 166     TH1D * hpm= hm->ProjectionX("mi",imin,imax) ;
 167     hpm->Sumw2() ;
 168     // if(i<=11){
 169     if(i<=10){
 170       hp ->Rebin(4) ;
 171       hpm->Rebin(4) ;
 172     }
 173     else{
 174       hp ->Rebin(5) ;
 175       hpm->Rebin(5) ;
 176     }
 177     for(Int_t ib=1; ib<=hp->GetNbinsX();ib++){if(hp ->GetBinContent(ib)==0)hp ->SetBinError(ib,1.);}
 178     for(Int_t ib=1; ib<=hp->GetNbinsX();ib++){if(hpm->GetBinContent(ib)==0)hpm->SetBinError(ib,1.);}
 179     TH1D * hpm2   = (TH1D*)hpm->Clone("Bg1") ;
 180     TH1D * hpcopy = (TH1D*)hp ->Clone("hpcopy") ;
 181     TH1D * hp2    = (TH1D*)hp ->Clone("hp2") ;
 182     hpcopy->SetXTitle("M_{#gamma#gamma} (GeV/c^{2})");
 183     hp2   ->SetXTitle("M_{#gamma#gamma} (GeV/c^{2})");
 184     hpcopy->Divide(hpm) ;
 185     sprintf(key,"%3.1f<p_{T}<%3.1f GeV/c",xa[i-1],xa[i]) ;
 186     hpcopy->SetTitle(key) ;
 187     hpcopy->SetMarkerStyle(20) ;
 188     hpcopy->SetMarkerSize(0.7) ;
 189
 190     fit1->SetParameters(0.001,0.136,0.005,0.0002,-0.002,0.0) ;
 191     fit1->SetParLimits(0,0.000,1.000) ;
 192     fit1->SetParLimits(1,0.125,0.145) ;
 193     fit1->SetParLimits(2,0.002,0.015) ;
 194
 195     Double_t rangeMin=0.06;
 196     Double_t rangeMax=0.26;
 197
 198     hpcopy->Fit(fit1,"Q" ,"",rangeMin,rangeMax) ;
 199     hpcopy->Fit(fit1,"MQ","",rangeMin,rangeMax) ;
 200
 201     ar1->SetBinContent(i,fit1->GetParameter(3)) ;
 202     ar1->SetBinError  (i,fit1->GetParError(3)) ;
 203     br1->SetBinContent(i,fit1->GetParameter(4)) ;
 204     br1->SetBinError  (i,fit1->GetParError(4)) ;
 205
 206     fit2->SetParameters(fit1->GetParameters()) ;
 207     fit2->SetParLimits(0,0.000,1.000) ;
 208     fit2->SetParLimits(1,0.125,0.145) ;
 209     fit2->SetParLimits(2,0.002,0.015) ;
 210
 211     hpcopy->Fit(fit2,"+NQ","",rangeMin,rangeMax) ;
 212     hpcopy->Fit(fit2,"+MQ","",rangeMin,rangeMax) ;
 213
 214     ar2->SetBinContent(i,fit2->GetParameter(3)) ;
 215     ar2->SetBinError  (i,fit2->GetParError(3)) ;
 216     br2->SetBinContent(i,fit2->GetParameter(4)) ;
 217     br2->SetBinError  (i,fit2->GetParError(4)) ;
 218     cr2->SetBinContent(i,fit2->GetParameter(5)) ;
 219     cr2->SetBinError  (i,fit2->GetParError(5)) ;
 220     hpcopy->SetAxisRange(0.04,0.3,"X") ;
 221     hpcopy->Draw() ;
 222     label ->Draw();
 223     if(c2)
 224       c2->Update() ;
 225     else
 226       c1->Update() ;
 227 //    if(getchar()=='q')return ;
 228
 229     fbg1->SetParameters(fit1->GetParameter(3),fit1->GetParameter(4),fit1->GetParameter(5));
 230     fbg2->SetParameters(fit2->GetParameter(3),fit2->GetParameter(4),fit2->GetParameter(5),
 231                         fit2->GetParameter(6));
 232
 233     Double_t intRangeMin = PeakPosition(pt,centrality)-3.*PeakWidth(pt,centrality) ;
 234     Double_t intRangeMax = PeakPosition(pt,centrality)+3.*PeakWidth(pt,centrality) ;
 235     // printf("pt=%f, %f < M < %f\n",pt,intRangeMin,intRangeMax);
 236     Int_t    intBinMin   = hp->GetXaxis()->FindBin(intRangeMin) ;
 237     Int_t    intBinMax   = hp->GetXaxis()->FindBin(intRangeMax) ;
 238     Double_t errStat     = hpm->Integral(intBinMin,intBinMax);
 239
 240     hpm ->Multiply(fbg1) ;
 241     hpm2->Multiply(fbg2) ;
 242     hp  ->Add(hpm ,-1.) ;
 243     hp2 ->Add(hpm2,-1.) ;
 244
 245     if(i<=15)
 246        c3->cd(i) ;
 247     else{
 248       if(c4==0){
 249         c4 = new TCanvas("mggFit2_Signal","mggFit2",10,10,1200,800) ;
 250         c4->Divide(nPadX,nPadY) ;
 251       }
 252       c4->cd(i-15) ;
 253     }
 254
 255     Int_t binPi0 = hp->FindBin(kMean);
 256     fgs->SetParameters(hp->Integral(binPi0-1,binPi0+1)/3.,fit1->GetParameter(1),fit1->GetParameter(2)) ;
 257     fgs->SetParLimits(0,0.000,1.e+5) ;
 258     fgs->SetParLimits(1,0.120,0.145) ;
 259     fgs->SetParLimits(2,0.002,0.015) ;
 260     hp->Fit(fgs,"Q","",rangeMin,rangeMax) ;
 261     hp->SetMaximum(hp2->GetMaximum()*1.4) ;
 262     hp->SetMinimum(hp2->GetMinimum()*1.1) ;
 263     hp->SetMarkerStyle(20) ;
 264     hp->SetMarkerSize(0.7) ;
 265     mr1->SetBinContent(i,fgs->GetParameter(1)) ;
 266     mr1->SetBinError  (i,fgs->GetParError(1) ) ;
 267     sr1->SetBinContent(i,TMath::Abs(fgs->GetParameter(2))) ;
 268     sr1->SetBinError  (i,fgs->GetParError(2) ) ;
 269
 270     Double_t y=fgs->Integral(intRangeMin,intRangeMax)/hp->GetXaxis()->GetBinWidth(1) ;
 271     nr1->SetBinContent(i,y) ;
 272     Double_t ey=0 ;
 273     if(fgs->GetParameter(0)!=0. && fgs->GetParameter(2)!=0.){
 274       Double_t en=fgs->GetParError(0)/fgs->GetParameter(0) ;
 275       Double_t es=fgs->GetParError(2)/fgs->GetParameter(2) ;
 276       ey=y*TMath::Sqrt(en*en+es*es) ;
 277     }
 278     nr1->SetBinError(i,ey) ;
 279
 280     Double_t npiInt = hp->Integral(intBinMin,intBinMax) ;
 281     Double_t norm   = fbg1->GetParameter(0) ;
 282     Double_t normErr= fbg1->GetParError(0) ;
 283     if(npiInt>0.){
 284       nr1int->SetBinContent(i,npiInt) ;
 285       nr1int->SetBinError(i,TMath::Sqrt(npiInt + norm*errStat + normErr*normErr*errStat*errStat + norm*norm*errStat)) ;
 286     }
 287     hp2->SetAxisRange(0.04,0.3,"X") ;
 288     hp2->SetMaximum(hp2->GetMaximum()*1.4) ;
 289     hp2->SetMinimum(hp2->GetMinimum()*1.1) ;
 290     hp2->SetMarkerStyle(20) ;
 291     hp2->SetMarkerSize(0.7) ;
 292     hp2->Fit(fgs,"Q","",rangeMin,rangeMax) ;
 293     mr2->SetBinContent(i,fgs->GetParameter(1)) ;
 294     mr2->SetBinError  (i,fgs->GetParError(1)) ;
 295     sr2->SetBinContent(i,TMath::Abs(fgs->GetParameter(2))) ;
 296     sr2->SetBinError  (i,fgs->GetParError(2)) ;
 297     y=fgs->Integral(intRangeMin,intRangeMax)/hp->GetXaxis()->GetBinWidth(1) ;
 298     nr2->SetBinContent(i,y) ;
 299     ey=0 ;
 300     if(fgs->GetParameter(0)!=0. && fgs->GetParameter(2)!=0.){
 301       Double_t en=fgs->GetParError(0)/fgs->GetParameter(0) ;
 302       Double_t es=fgs->GetParError(2)/fgs->GetParameter(2) ;
 303       ey=y*TMath::Sqrt(en*en+es*es) ;
 304     }
 305     nr2->SetBinError(i,ey) ;
 306     npiInt=hp2->Integral(intBinMin,intBinMax) ;
 307     norm=fbg2->GetParameter(0) ;
 308     normErr=fbg2->GetParError(0) ;
 309     if(npiInt>0.){
 310       nr2int->SetBinContent(i,npiInt) ;
 311       nr2int->SetBinError(i,TMath::Sqrt(npiInt + norm*errStat + normErr*normErr*errStat*errStat + norm*norm*errStat)) ;
 312     }
 313     hp2->SetTitle(key) ;
 314     hp2->Draw() ;
 315     label ->Draw();
 316     if(c4)
 317       c4->Update() ;
 318     else
 319       c3->Update() ;
 320
 321     delete hp ;
 322 //    delete hp2 ;
 323 //    delete hpcopy ;
 324     delete hpm ;
 325     delete hpm2 ;
 326   }
 327   char name[55] ;
 328
 329   sprintf(name,"Pi0_ratio_cent%d%s.eps" ,centrality,cModule) ;
 330   c1->Print(name) ;
 331   sprintf(name,"Pi0_signal_cent%d%s.eps",centrality,cModule) ;
 332   c3->Print(name) ;
 333
 334   //Normalize by the number of events
 335   Int_t cMin,cMax;
 336   if      (centrality == 0) {
 337     cMin=1;
 338     cMax=20;
 339   }
 340   else if (centrality == 1) {
 341     cMin=21;
 342     cMax=80;
 343   }
 344   Double_t nevents = hCentrality1->Integral(cMin,cMax);
 345   nr1   ->Scale(1./nevents) ;
 346   nr1int->Scale(1./nevents) ;
 347   nr2   ->Scale(1./nevents) ;
 348   nr2int->Scale(1./nevents) ;
 349
 350   printf("\t==============================\n");
 351   printf("\t|       N events = %i8    |\n",nevents);
 352   printf("\t==============================\n");
 353
 354   TFile fout("LHC10h_Pi0_FitResult.root","update");
 355   mr1->Write(0,TObject::kOverwrite) ;
 356   sr1->Write(0,TObject::kOverwrite) ;
 357   ar1->Write(0,TObject::kOverwrite) ;
 358   br1->Write(0,TObject::kOverwrite) ;
 359   nr1->Write(0,TObject::kOverwrite) ;
 360   nr1int->Write(0,TObject::kOverwrite) ;
 361   ar2->Write(0,TObject::kOverwrite) ;
 362   br2->Write(0,TObject::kOverwrite) ;
 363   cr2->Write(0,TObject::kOverwrite) ;
 364   mr2->Write(0,TObject::kOverwrite) ;
 365   sr2->Write(0,TObject::kOverwrite) ;
 366   nr2->Write(0,TObject::kOverwrite) ;
 367   nr2int->Write(0,TObject::kOverwrite) ;
 368   fout.Close() ;
 369
 370 }
 371
 372 //-----------------------------------------------------------------------------
 373 Double_t PeakPosition(Double_t pt, Int_t centrality = 0){
 374   //Fit to the measured peak position
 375   Double_t pi0Peak = 0.135;
 376   if      (centrality == 0) pi0Peak=0.1413;
 377   else if (centrality == 1) pi0Peak=0.1380;
 378   else if (centrality == 2) pi0Peak=0.1367;
 379   else if (centrality == 3) pi0Peak=0.1359;
 380   return pi0Peak;
 381   // return 4.99292e-003*exp(-pt*9.32300e-001)+1.34944e-001 ;
 382 }
 383 //-----------------------------------------------------------------------------
 384 Double_t PeakWidth(Double_t pt, Int_t centrality = 0){
 385   //fit to the measured peak width
 386   Double_t pi0Sigma = 0.07;
 387   if      (centrality == 0) pi0Sigma=0.0084;
 388   else if (centrality == 1) pi0Sigma=0.0071;
 389   else if (centrality == 2) pi0Sigma=0.0068;
 390   else if (centrality == 3) pi0Sigma=0.0064;
 391   return pi0Sigma;
 392
 393   // Double_t a=0.0068 ;
 394   // Double_t b=0.0025 ;
 395   // Double_t c=0.000319 ;
 396   // return TMath::Sqrt(a*a+b*b/pt/pt+c*c*pt*pt) ;
 397 }
 398
 399 //-----------------------------------------------------------------------------
 400 Double_t CB(Double_t * x, Double_t * par){
 401   //Parameterization of Real/Mixed ratio
 402   Double_t m=par[1] ;
 403   Double_t s=par[2] ;
 404   Double_t dx=(x[0]-m)/s ;
 405   Double_t dMpi0 = x[0]-kMean;
 406   Double_t peak = par[0] * TMath::Exp(-dx*dx/2.);
 407   Double_t bg   = par[3] + par[4]*dMpi0 + par[5]*dMpi0*dMpi0;
 408   return peak+bg ;
 409 }
 410
 411 //-----------------------------------------------------------------------------
 412 Double_t CB2(Double_t * x, Double_t * par){
 413   //Another parameterization of Real/Mixed ratio7TeV/README
 414   Double_t m=par[1] ;
 415   Double_t s=par[2] ;
 416   Double_t dx=(x[0]-m)/s ;
 417   Double_t dMpi0 = x[0]-kMean;
 418   Double_t peak = par[0] * TMath::Exp(-dx*dx/2.);
 419   Double_t bg   = par[3] + par[4]*dMpi0 + par[5]*dMpi0*dMpi0 + par[6]*dMpi0*dMpi0*dMpi0;
 420   return peak+bg ;
 421 }
 422 //-----------------------------------------------------------------------------
 423 Double_t CBs(Double_t * x, Double_t * par){
 424   //Parameterizatin of signal
 425   Double_t m=par[1] ;
 426   Double_t s=par[2] ;
 427   Double_t dx=(x[0]-m)/s ;
 428   Double_t peak = par[0] * TMath::Exp(-dx*dx/2.);
 429   return peak ;
 430 }
 431 //-----------------------------------------------------------------------------
 432 Double_t BG1(Double_t * x, Double_t * par){
 433   //Normalizatino of Mixed
 434   return par[0]+par[1]*(x[0]-kMean)+par[2]*(x[0]-kMean)*(x[0]-kMean) ;
 435 }
 436 //-----------------------------------------------------------------------------
 437 Double_t BG2(Double_t * x, Double_t * par){
 438   //Another normalization of  Mixed
 439   return par[0]+par[1]*(x[0]-kMean)+par[2]*(x[0]-kMean)*(x[0]-kMean)+par[3]*(x[0]-kMean)*(x[0]-kMean)*(x[0]-kMean) ;
 440 }
 441
 442
 443 //-----------------------------------------------------------------------------
 444 PPRstyle()
 445 {
 446
 447   //////////////////////////////////////////////////////////////////////
 448   //
 449   // ROOT style macro for the TRD TDR
 450   //
 451   //////////////////////////////////////////////////////////////////////
 452
 453   gStyle->SetPalette(1);
 454   gStyle->SetCanvasBorderMode(-1);
 455   gStyle->SetCanvasBorderSize(1);
 456   gStyle->SetCanvasColor(10);
 457
 458   gStyle->SetFrameFillColor(10);
 459   gStyle->SetFrameBorderSize(1);
 460   gStyle->SetFrameBorderMode(-1);
 461   gStyle->SetFrameLineWidth(1.2);
 462   gStyle->SetFrameLineColor(1);
 463
 464   gStyle->SetHistFillColor(0);
 465   gStyle->SetHistLineWidth(1);
 466   gStyle->SetHistLineColor(1);
 467
 468   gStyle->SetPadColor(10);
 469   gStyle->SetPadBorderSize(1);
 470   gStyle->SetPadBorderMode(-1);
 471
 472   gStyle->SetStatColor(10);
 473   gStyle->SetTitleColor(kBlack,"X");
 474   gStyle->SetTitleColor(kBlack,"Y");
 475
 476   gStyle->SetLabelSize(0.04,"X");
 477   gStyle->SetLabelSize(0.04,"Y");
 478   gStyle->SetLabelSize(0.04,"Z");
 479   gStyle->SetTitleSize(0.04,"X");
 480   gStyle->SetTitleSize(0.04,"Y");
 481   gStyle->SetTitleSize(0.04,"Z");
 482   gStyle->SetTitleFont(42,"X");
 483   gStyle->SetTitleFont(42,"Y");
 484   gStyle->SetTitleFont(42,"X");
 485   gStyle->SetLabelFont(42,"X");
 486   gStyle->SetLabelFont(42,"Y");
 487   gStyle->SetLabelFont(42,"Z");
 488   gStyle->SetStatFont(42);
 489
 490   gStyle->SetTitleOffset(1.0,"X");
 491   gStyle->SetTitleOffset(1.4,"Y");
 492
 493   gStyle->SetFillColor(kWhite);
 494   gStyle->SetTitleFillColor(kWhite);
 495
 496   gStyle->SetOptDate(0);
 497   gStyle->SetOptTitle(1);
 498   gStyle->SetOptStat(0);
 499   gStyle->SetOptFit(0);
 500
 501 }
 502