PWG4/macros/electrons/makeCombinedData.C

   1 /////////////////////////////////////////
   2 // A set of plotting utilities for the
   3 // electron chapter of the PPR
   4 //
   5 // J.L.Klay (Cal Poly)
   6 // 28-Oct-2009
   7 ////////////////////////////////////////
   8
   9 TH1F *alltte,    *alltrk,    *allemc;
  10 TH1F *sumtte,    *sumtrk,    *sumemc;  //all but misid'ed hadrons
  11 TH1F *btte,      *btrk,      *bemc;
  12 TH1F *ctte,      *ctrk,      *cemc;
  13 TH1F *cbtte,     *cbtrk,     *cbemc;
  14 TH1F *convtte,   *convtrk,   *convemc;
  15 TH1F *daltte,    *daltrk,    *dalemc;
  16 TH1F *wztte,     *wztrk,     *wzemc;
  17 TH1F *htte,      *htrk,      *hemc;
  18
  19 TH1F *allMC;
  20 TH1F *bMC, *cMC, *cbMC;
  21 TH1F *convMC, *dalMC, *wzMC;
  22 TH1F *mchad;
  23 TH1F *allheratio, *behratio;
  24
  25 TF1* fpow;
  26
  27 void makeData(char* jjfname = "data/scaled25Oct09/TOTALhistosscaled-LHC09b2-0.root",
  28               char* bfname = "data/scaled25Oct09/histosscaledLHC09b4AODc.root",
  29               char* wfname = "data/scaled25Oct09/histosWboson.root") {
  30
  31   //TO AVOID DOUBLE-COUNTING ELECTRONS:
  32   //NOTE:  Jet-Jet events are "minimum bias" which means that there
  33   //are B-Jets included and for the higher pThard bins the B-jets will
  34   //play a more important role.  So to AVOID DOUBLE-COUNTING
  35   //B-electrons, we will suppress the electrons in the Jet-Jet events
  36   //that come from B->e or B->C->e when we combine the productions
  37   //into the total histograms
  38
  39   //TO AVOID DOUBLE-COUNTING HADRONS:
  40   //For the hadron yields, we need to use only the yields from the
  41   //Jet-Jet events + W events and NOT USE the hadrons from the B-Jet events
  42
  43   //For bjet and jet-jet events
  44   double pyscale = (1.E6)*0.5*208*208*100/360; //seconds*lumi*Pb*Pb*acceptance
  45   double bscale = pyscale; //Do we need to scale by Branching ratio
  46                            //for forced semi-leptonic decays?
  47   double wscale = pyscale;
  48
  49   TFile* jjfile = new TFile(jjfname);
  50   if(!jjfile) { printf("NO JET-JET FILE\n"); return; }
  51   TH2F* jjtte = (TH2F*)jjfile->Get("AnaElectron_hPtNPEleTTEScaled");
  52   TH2F* jjemc = (TH2F*)jjfile->Get("AnaElectron_hPtNPEleEMCALScaled");
  53   TH2F* jjtrk = (TH2F*)jjfile->Get("AnaElectron_hPtNPEleTPCTRDScaled");
  54   TH1F* jjmult = (TH1F*)jjfile->Get("AnaElectron_hRefMultScaled");
  55   Int_t nEvtJJ = jjmult->GetEntries();
  56   jjtte->Scale(pyscale);
  57   jjemc->Scale(pyscale);
  58   jjtrk->Scale(pyscale);
  59   TH2F* jjmcele = (TH2F*)jjfile->Get("AnaElectron_hPtMCElectronScaled");
  60   TH1F* jjmchad = (TH1F*)jjfile->Get("AnaElectron_hPtMCHadronScaled");
  61   jjmcele->Scale(pyscale);
  62   jjmchad->Scale(pyscale);
  63
  64   TFile* bfile = new TFile(bfname);
  65   if(!bfile) { printf("NO B-JET FILE\n"); return; }
  66   TH2F* bjtte = (TH2F*)bfile->Get("AnaElectron_hPtNPEleTTEScaled");
  67   TH2F* bjemc = (TH2F*)bfile->Get("AnaElectron_hPtNPEleEMCALScaled");
  68   TH2F* bjtrk = (TH2F*)bfile->Get("AnaElectron_hPtNPEleTPCTRDScaled");
  69   TH1F* bjmult = (TH1F*)bfile->Get("AnaElectron_hRefMultScaled");
  70   Int_t nEvtB = bjmult->GetEntries();
  71   bjtte->Scale(bscale);
  72   bjemc->Scale(bscale);
  73   bjtrk->Scale(bscale);
  74   TH2F* bjmcele = (TH2F*)bfile->Get("AnaElectron_hPtMCElectronScaled");
  75   TH1F* bjmchad = (TH1F*)bfile->Get("AnaElectron_hPtMCHadronScaled");
  76   bjmcele->Scale(bscale);
  77   bjmchad->Scale(bscale);
  78
  79   TFile* wfile = new TFile(wfname);
  80   if(!wfile) { printf("NO W-BOSON FILE\n"); return; }
  81   TH2F* wjtte = (TH2F*)wfile->Get("AnaElectron_hPtNPEleTTE");
  82   TH2F* wjemc = (TH2F*)wfile->Get("AnaElectron_hPtNPEleEMCAL");
  83   TH2F* wjtrk = (TH2F*)wfile->Get("AnaElectron_hPtNPEleTPCTRD");
  84   TH1F* wjmult = (TH1F*)wfile->Get("AnaElectron_hRefMult");
  85   Int_t nEvtW = wjmult->GetEntries();
  86   wjtte->Scale(wscale); //already scaled by evts
  87   wjemc->Scale(wscale); //already scaled by evts
  88   wjtrk->Scale(wscale); //already scaled by evts
  89   TH2F* wjmcele = (TH2F*)wfile->Get("AnaElectron_hPtMCElectron");
  90   TH1F* wjmchad = (TH1F*)wfile->Get("AnaElectron_hPtMCHadron");
  91   wjmcele->Scale(wscale);
  92   wjmchad->Scale(wscale);
  93
  94   printf("Event statistics: %d (JET-JET)  %d (B-JET)  %d (W-Boson)\n",nEvtJJ,nEvtB,nEvtW);
  95
  96   makeMCElectrons(jjmcele,bjmcele,wjmcele,jjmchad,wjmchad);
  97   makeTTEElectrons(jjtte,bjtte,wjtte);
  98   makeEMCElectrons(jjemc,bjemc,wjemc);
  99   makeTRKElectrons(jjtrk,bjtrk,wjtrk);
 100
 101   TFile* fout = new TFile("CombinedCocktailHistograms.root","RECREATE");
 102   fout->cd();
 103   alltte->Write();
 104   alltrk->Write();
 105   allemc->Write();
 106
 107   sumtte->Write();
 108   sumtrk->Write();
 109   sumemc->Write();
 110
 111   btte->Write();
 112   btrk->Write();
 113   bemc->Write();
 114
 115   ctte->Write();
 116   ctrk->Write();
 117   cemc->Write();
 118
 119   cbtte->Write();
 120   cbtrk->Write();
 121   cbemc->Write();
 122
 123   convtte->Write();
 124   convtrk->Write();
 125   convemc->Write();
 126
 127   daltte->Write();
 128   daltrk->Write();
 129   dalemc->Write();
 130
 131   wztte->Write();
 132   wztrk->Write();
 133   wzemc->Write();
 134
 135   htte->Write();
 136   htrk->Write();
 137   hemc->Write();
 138
 139   allMC->Write();
 140   bMC->Write();
 141   cMC->Write();
 142   cbMC->Write();
 143   convMC->Write();
 144   dalMC->Write();
 145   wzMC->Write();
 146   mchad->Write();
 147
 148   allheratio->Write();
 149   behratio->Write();
 150
 151   fout->Close();
 152
 153 }
 154
 155
 156 void ScaleAndConfigure(TH1F* hist,Double_t scale, Int_t color,Bool_t keepErr)
 157 {
 158   hist->Scale(scale);
 159   hist->SetLineColor(color);
 160   hist->SetMarkerColor(color);
 161   hist->SetLineWidth(2);
 162   if(keepErr == kFALSE) {
 163     //remove the error bars - useful for MC rates
 164     for(Int_t i = 1; i <= hist->GetNbinsX(); i++) {
 165       if(hist->GetBinContent(i) > 0.) {
 166         if(hist->GetBinError(i)/hist->GetBinContent(i) > 0.5) {
 167           Double_t avg = 0.;
 168           if(i > 1 && i < hist->GetNbinsX())
 169             avg = (hist->GetBinContent(i-1) + hist->GetBinContent(i+1))/2.;
 170           hist->SetBinContent(i,avg);
 171         }
 172       }
 173       hist->SetBinError(i,0.);
 174     }
 175   } else {
 176     //Set the error bars to statistics of the bin content
 177     for(Int_t i = 1; i <= hist->GetNbinsX(); i++) {
 178       if(hist->GetBinContent(i) > 0.) {
 179         if(hist->GetBinError(i)/hist->GetBinContent(i) > 0.5) {
 180           Double_t avg = 0;
 181           if(i > 1 && i < hist->GetNbinsX())
 182             avg = (hist->GetBinContent(i-1) + hist->GetBinContent(i+1))/2.;
 183           hist->SetBinContent(i,avg);
 184         }
 185       }
 186       hist->SetBinError(i,TMath::Sqrt(hist->GetBinContent(i)));
 187     }
 188   }
 189 }
 190
 191 TH1F* GetPtCutHisto(TH1F* input)
 192 {
 193   //Given a rate histogram vs pt, return the histogram with yield
 194   //above a given pTcut
 195
 196   TH1F* result = (TH1F*)input->Clone();
 197   char name[100];
 198   sprintf(name,"%s_ptCut",result->GetName());
 199   //  printf("####### %s #########\n",name);
 200   result->SetNameTitle(name,name);
 201   for(Int_t i = 1; i <= result->GetNbinsX(); i++) {
 202     Double_t val = input->Integral(i,result->GetNbinsX());
 203     //printf("<Bin %d> value %.2f integral above %.2f\n",i,input->GetBinContent(i),val);
 204     result->SetBinContent(i,val);
 205     result->SetBinError(i,0.);
 206   }
 207   //printf("################\n");
 208   return result;
 209
 210 }
 211
 212 void makeMCElectrons(TH2F* jjmcele, TH2F* bjmcele, TH2F* wjmcele, TH1F* jjmchad, TH1F* wjmchad) {
 213
 214   //Jet-Jet MC electrons
 215   TH1F* jjallmc = (TH1F*)jjmcele->ProjectionX("jjallmc",1,1);
 216   TH1F* jjbmc = (TH1F*)jjmcele->ProjectionX("jjbmc",2,2);
 217   TH1F* jjcmc = (TH1F*)jjmcele->ProjectionX("jjcmc",3,3);
 218   TH1F* jjcandbmc = (TH1F*)jjmcele->ProjectionX("jjcandbmc",4,4);
 219   TH1F* jjconvmc = (TH1F*)jjmcele->ProjectionX("jjconvmc",5,5);
 220   TH1F* jjdalmc = (TH1F*)jjmcele->ProjectionX("jjdalmc",6,6);
 221
 222   //Bottom-Jet MC electrons
 223   TH1F* bjallmc = (TH1F*)bjmcele->ProjectionX("bjallmc",1,1);
 224   TH1F* bjbmc = (TH1F*)bjmcele->ProjectionX("bjbmc",2,2);
 225   TH1F* bjcmc = (TH1F*)bjmcele->ProjectionX("bjcmc",3,3);
 226   TH1F* bjcandbmc = (TH1F*)bjmcele->ProjectionX("bjcandbmc",4,4);
 227   TH1F* bjconvmc = (TH1F*)bjmcele->ProjectionX("bjconvmc",5,5);
 228   TH1F* bjdalmc = (TH1F*)bjmcele->ProjectionX("bjdalmc",6,6);
 229
 230   //W-Jet MC electrons
 231   TH1F* wjallmc = (TH1F*)wjmcele->ProjectionX("wjallmc",1,1);
 232   TH1F* wjbmc = (TH1F*)wjmcele->ProjectionX("wjbmc",2,2);
 233   TH1F* wjcmc = (TH1F*)wjmcele->ProjectionX("wjcmc",3,3);
 234   TH1F* wjcandbmc = (TH1F*)wjmcele->ProjectionX("wjcandbmc",4,4);
 235   TH1F* wjconvmc = (TH1F*)wjmcele->ProjectionX("wjconvmc",5,5);
 236   TH1F* wjdalmc = (TH1F*)wjmcele->ProjectionX("wjdalmc",6,6);
 237   TH1F* wjwzmc = (TH1F*)wjmcele->ProjectionX("wjwzmc",7,7);
 238
 239   //MC Hadrons (from jj events only)
 240   TCanvas *ctemp = new TCanvas("ctemp");
 241   ctemp->Divide(2,3);
 242   ctemp->cd(1); gPad->SetLogy(); jjmchad->Draw();
 243   mchad = (TH1F*)jjmchad->Clone(); mchad->SetName("mchad");
 244   TH1F* temphad = (TH1F*)mchad->Clone(); temphad->SetName("temphad");
 245   smoothWithFit(temphad,1e5,-3,10,60,100);
 246   for(Int_t i = 10; i<= mchad->GetNbinsX(); i++) {
 247     Double_t pt = mchad->GetBinCenter(i);
 248     mchad->SetBinContent(i,temphad->GetFunction("fpow")->Eval(pt));
 249   }
 250   jjmchad->Draw();
 251   temphad->Draw("same");
 252   mchad->Draw("same");
 253   mchad->SetTitle("MC hadrons in Pb+Pb in EMCAL acceptance");
 254   mchad->SetName("mchad");
 255   mchad->SetXTitle("p_T (GeV/c)");
 256
 257   bMC = (TH1F*)bjbmc->Clone(); bMC->SetName("bMC");  //B-Jet + W-jet
 258   bMC->Add(wjbmc);
 259   ctemp->cd(2); gPad->SetLogy(); bMC->Draw();
 260   TH1F* foob = (TH1F*)bMC->Clone(); foob->SetName("foob");
 261   smoothWithFit(foob,1e6,-3,15,50,100);
 262   for(Int_t i = 10; i<= bMC->GetNbinsX(); i++) bMC->SetBinContent(i,foob->GetBinContent(i));
 263   bMC->Draw();
 264   foob->Draw("same");
 265
 266   cMC = (TH1F*)jjcmc->Clone(); cMC->SetName("cMC"); //Jet-Jet + W-jet
 267   cMC->Add(wjcmc);
 268   ctemp->cd(3); gPad->SetLogy(); cMC->Draw();
 269   TH1F* fooc = (TH1F*)cMC->Clone(); fooc->SetName("fooc");
 270   smoothWithFit(fooc,1e6,-3,5,30,100);
 271   for(Int_t i = 5; i<= cMC->GetNbinsX(); i++) cMC->SetBinContent(i,fooc->GetBinContent(i));
 272   cMC->Draw();
 273   fooc->Draw("same");
 274
 275   cbMC = (TH1F*)bjcandbmc->Clone(); cbMC->SetName("cbMC"); //B-Jet + W-jet
 276   cbMC->Add(wjcandbmc);
 277   ctemp->cd(4); gPad->SetLogy(); cbMC->Draw();
 278   TH1F* foocb = (TH1F*)cbMC->Clone(); foocb->SetName("foocb");
 279   smoothWithFit(foocb,1e6,-3,8,35,100);
 280   for(Int_t i = 8; i<= cbMC->GetNbinsX(); i++) cbMC->SetBinContent(i,foocb->GetBinContent(i));
 281   cbMC->Draw();
 282   foocb->Draw("same");
 283
 284   convMC = (TH1F*)jjconvmc->Clone(); convMC->SetName("convMC"); //Jet-Jet + W-jet
 285   convMC->Add(wjconvmc);
 286   ctemp->cd(5); gPad->SetLogy(); convMC->Draw();
 287   TH1F* fooconv = (TH1F*)convMC->Clone(); fooconv->SetName("fooconv");
 288   smoothWithFit(fooconv,1e6,-3,6,40,100);
 289   //  for(Int_t i = 6; i<= convMC->GetNbinsX(); i++) convMC->SetBinContent(i,fooconv->GetBinContent(i));
 290   convMC->Draw();
 291   fooconv->Draw("same");
 292
 293   dalMC = (TH1F*)jjdalmc->Clone(); dalMC->SetName("dalMC"); //Jet-Jet + W-jet
 294   dalMC->Add(wjdalmc);
 295   ctemp->cd(6); gPad->SetLogy(); dalMC->Draw();
 296   TH1F* foodal = (TH1F*)dalMC->Clone(); foodal->SetName("foodal");
 297   smoothWithFit(foodal,1e6,-3,3,20,100);
 298   for(Int_t i = 3; i<= dalMC->GetNbinsX(); i++) dalMC->SetBinContent(i,foodal->GetBinContent(i));
 299   dalMC->Draw();
 300   foodal->Draw("same");
 301
 302   wzMC = (TH1F*)wjwzmc->Clone(); wzMC->SetName("wzMC"); //W-jet only
 303   TCanvas* cw= new TCanvas("cw");
 304   cw->cd(); gPad->SetLogy(); wzMC->GetYaxis()->SetRangeUser(1,1e3); wzMC->Draw();
 305   TH1F* foowz = (TH1F*)wzMC->Clone(); foowz->SetName("foowz");
 306   TF1* fws = new TF1("fws","[0]*(1+exp((x-[1])/[2]))^-1",39,100);
 307   fws->SetParameters(100,30,5);
 308   foowz->Fit(fws,"R");
 309   TF1* fwzexp = new TF1("fwzexp","[0]+[1]*log(x/[2])^2",4,40);
 310   fwzexp->SetParameters(100,10,3);
 311   foowz->Fit(fwzexp,"R");
 312   for(Int_t i = 8; i<= wzMC->GetNbinsX(); i++) {
 313     Double_t pt = wzMC->GetBinCenter(i);
 314     if(pt < 40) wzMC->SetBinContent(i,fwzexp->Eval(pt));
 315     if(pt > 40) wzMC->SetBinContent(i,fws->Eval(pt));
 316   }
 317   wzMC->GetYaxis()->SetRangeUser(1,1e3);
 318   wzMC->Draw();
 319   fws->Draw("same");
 320   foowz->Draw("same");
 321
 322   //All mc electrons is the sum of
 323   //Jet-Jet: conversions + direct charm + dalitz
 324   //Bottom-Jet: direct bottom + indirect bottom
 325   //W-Jet: all (because these events are exclusive of the others)
 326   allMC = (TH1F*)wzMC->Clone(); allMC->SetNameTitle("allMC","All MC Electrons");
 327   allMC->Add(convMC); allMC->Add(cMC); allMC->Add(dalMC);
 328   allMC->Add(bMC); allMC->Add(cbMC);
 329
 330   //Hadron/electron ratios
 331   allheratio = (TH1F*)allMC->Clone(); allheratio->SetName("allheratio");
 332   behratio = (TH1F*)bMC->Clone(); behratio->SetName("behratio");
 333   allheratio->SetTitle("MC hadrons and electrons in Pb+Pb, 5.5 TeV");
 334   allheratio->SetXTitle("p_{T} (GeV/c)");
 335   allheratio->SetYTitle("Hadrons/Electrons");
 336   for(Int_t i = 1; i <= allheratio->GetNbinsX(); i++) {
 337     Double_t vale = allMC->GetBinContent(i);
 338     Double_t valb = bMC->GetBinContent(i);
 339     Double_t valh = mchad->GetBinContent(i);
 340     //printf("pT %.2f, Hadron %.1f, Electron %.1f, B-electron
 341     //%.1f\n",all->GetBinCenter(i),valh,vale,valb);
 342     if(vale>0) allheratio->SetBinContent(i,valh/vale);
 343     else allheratio->SetBinContent(i,0.);
 344
 345     if(valb>0) behratio->SetBinContent(i,valh/valb);
 346     else behratio->SetBinContent(i,0.);
 347
 348     allheratio->SetBinError(i,0.);
 349     behratio->SetBinError(i,0.);
 350   }
 351
 352   double myscale = 1.; //we already scaled them
 353   ScaleAndConfigure(allMC,myscale,kBlack,kFALSE);
 354   ScaleAndConfigure(bMC,myscale,kRed,kFALSE);
 355   ScaleAndConfigure(cMC,myscale,kBlue,kFALSE);
 356   ScaleAndConfigure(cbMC,myscale,kViolet,kFALSE);
 357   ScaleAndConfigure(convMC,myscale,kOrange-3,kFALSE);
 358   ScaleAndConfigure(dalMC,myscale,kGreen-3,kFALSE);
 359   ScaleAndConfigure(wzMC,myscale,kOrange-7,kFALSE);
 360   ScaleAndConfigure(mchad,myscale,kGreen+2,kFALSE);
 361
 362   return;
 363 }
 364
 365 void makeTTEElectrons(TH2F* jjtte, TH2F* bjtte, TH2F* wjtte) {
 366
 367   //Jet-Jet TTE Electrons
 368   TH1F* jjalltte = (TH1F*)jjtte->ProjectionX("jjalltte",1,1);
 369   TH1F* jjbtte = (TH1F*)jjtte->ProjectionX("jjbtte",2,2);
 370   TH1F* jjctte = (TH1F*)jjtte->ProjectionX("jjctte",3,3);
 371   TH1F* jjcbtte = (TH1F*)jjtte->ProjectionX("jjcbtte",4,4);
 372   TH1F* jjconvtte = (TH1F*)jjtte->ProjectionX("jjconvtte",5,5);
 373   TH1F* jjdaltte = (TH1F*)jjtte->ProjectionX("jjdaltte",6,6);
 374   TH1F* jjwztte = (TH1F*)jjtte->ProjectionX("jjwztte",7,7);
 375   TH1F* jjhtte = (TH1F*)jjtte->ProjectionX("jjhtte",9,9);
 376
 377   //B-Jet TTE Electrons
 378   TH1F* bjalltte = (TH1F*)bjtte->ProjectionX("bjalltte",1,1);
 379   TH1F* bjbtte = (TH1F*)bjtte->ProjectionX("bjbtte",2,2);
 380   TH1F* bjctte = (TH1F*)bjtte->ProjectionX("bjctte",3,3);
 381   TH1F* bjcbtte = (TH1F*)bjtte->ProjectionX("bjcbtte",4,4);
 382   TH1F* bjconvtte = (TH1F*)bjtte->ProjectionX("bjconvtte",5,5);
 383   TH1F* bjdaltte = (TH1F*)bjtte->ProjectionX("bjdaltte",6,6);
 384   TH1F* bjwztte = (TH1F*)bjtte->ProjectionX("bjwztte",7,7);
 385   TH1F* bjhtte = (TH1F*)bjtte->ProjectionX("bjhtte",9,9);
 386
 387   //W-Jet TTE Electrons
 388   TH1F* wjalltte = (TH1F*)wjtte->ProjectionX("wjalltte",1,1);
 389   TH1F* wjbtte = (TH1F*)wjtte->ProjectionX("wjbtte",2,2);
 390   TH1F* wjctte = (TH1F*)wjtte->ProjectionX("wjctte",3,3);
 391   TH1F* wjcbtte = (TH1F*)wjtte->ProjectionX("wjcbtte",4,4);
 392   TH1F* wjconvtte = (TH1F*)wjtte->ProjectionX("wjconvtte",5,5);
 393   TH1F* wjdaltte = (TH1F*)wjtte->ProjectionX("wjdaltte",6,6);
 394   TH1F* wjwztte = (TH1F*)wjtte->ProjectionX("wjwztte",7,7);
 395   TH1F* wjhtte = (TH1F*)wjtte->ProjectionX("wjhtte",9,9);
 396
 397   btte = (TH1F*)bjbtte->Clone(); btte->SetName("btte");  //B-Jet + W-jet
 398   btte->Add(wjbtte)
 399   TCanvas * ctemptte = new TCanvas("ctemptte");
 400   ctemptte->Divide(2,3);
 401   ctemptte->cd(1); gPad->SetLogy(); btte->Draw();
 402   TH1F* foobtte = (TH1F*)btte->Clone(); foobtte->SetName("foobtte");
 403   smoothWithFit(foobtte,1e5,-3,8,40,100);
 404   btte->Draw(); foobtte->Draw("same");
 405
 406   ctte = (TH1F*)jjctte->Clone(); ctte->SetName("ctte"); //Jet-Jet + W-jet
 407   ctte->Add(wjctte);
 408   ctemptte->cd(2); gPad->SetLogy(); ctte->Draw();
 409   TH1F* fooctte = (TH1F*)ctte->Clone(); fooctte->SetName("fooctte");
 410   smoothWithFit(fooctte,1e5,-3,3,12,100);
 411   ctte->Draw(); fooctte->Draw("same");
 412
 413   cbtte = (TH1F*)bjcbtte->Clone(); cbtte->SetName("cbtte"); //B-Jet + W-jet
 414   cbtte->Add(wjcbtte);
 415   ctemptte->cd(3); gPad->SetLogy(); cbtte->Draw();
 416   TH1F* foocbtte = (TH1F*)cbtte->Clone(); foocbtte->SetName("foocbtte");
 417   smoothWithFit(foocbtte,1e5,-3,5,20,100);
 418   cbtte->Draw(); foocbtte->Draw("same");
 419
 420   convtte = (TH1F*)jjconvtte->Clone(); convtte->SetName("convtte"); //Jet-Jet + W-jet
 421   convtte->Add(wjconvtte);
 422   ctemptte->cd(4); gPad->SetLogy(); convtte->Draw();
 423   TH1F* fooconvtte = (TH1F*)convtte->Clone(); fooconvtte->SetName("fooconvtte");
 424   smoothWithFit(fooconvtte,1e5,-3,6,13,100);
 425   convtte->Draw(); fooconvtte->Draw("same");
 426
 427   daltte = (TH1F*)jjdaltte->Clone(); daltte->SetName("daltte"); //Jet-Jet + W-jet
 428   daltte->Add(wjdaltte);
 429   ctemptte->cd(5); gPad->SetLogy(); daltte->Draw();
 430   TH1F* foodaltte = (TH1F*)daltte->Clone(); foodaltte->SetName("foodaltte");
 431   smoothWithFit(foodaltte,1e2,-5,2,10,100);
 432   daltte->Draw(); foodaltte->Draw("same");
 433
 434   htte = (TH1F*)jjhtte->Clone(); htte->SetName("htte");
 435   htte->Add(wjhtte);
 436   ctemptte->cd(6); gPad->SetLogy(); htte->Draw();
 437   TH1F* foohtte = (TH1F*)htte->Clone(); foohtte->SetName("foohtte");
 438   smoothWithFit(foohtte,1e5,-3,10,40,100);
 439   htte->Draw(); foohtte->Draw("same");
 440
 441   wztte = (TH1F*)wjwztte->Clone(); wztte->SetName("wztte"); //W-jet only
 442   TCanvas* ctempwztte = new TCanvas("ctempwztte");
 443   ctempwztte->cd(); gPad->SetLogy(); wztte->Draw();
 444   TH1F* foowztte = (TH1F*)wztte->Clone(); foowztte->SetName("foowztte");
 445   smoothWithFit(foowztte,1e5,-3,10,40,100);
 446   wztte->Draw(); foowztte->Draw("same");
 447
 448   //All TTE electrons is the sum of
 449   //Jet-Jet: conversions + direct charm + dalitz + misid
 450   //Bottom-Jet: direct bottom + indirect bottom
 451   //W-Jet: all (because these events are exclusive of the others)
 452   alltte = (TH1F*)wjalltte->Clone(); alltte->SetName("alltte");
 453   alltte->Add(jjconvtte); alltte->Add(jjctte); alltte->Add(jjdaltte); alltte->Add(jjhtte);
 454   alltte->Add(bjbtte); alltte->Add(bjcbtte);
 455   sumtte = (TH1F*)wjalltte->Clone(); sumtte->SetName("sumtte");
 456   sumtte->Add(jjconvtte); sumtte->Add(jjctte); sumtte->Add(jjdaltte);
 457   sumtte->Add(bjbtte); sumtte->Add(bjcbtte);
 458
 459   double myscale = 1.; //we already scaled them
 460   ScaleAndConfigure(alltte,myscale,kBlack,kFALSE);
 461   ScaleAndConfigure(sumtte,myscale,kBlack,kFALSE);
 462   ScaleAndConfigure(btte,myscale,kRed,kFALSE);
 463   ScaleAndConfigure(ctte,myscale,kBlue,kFALSE);
 464   ScaleAndConfigure(cbtte,myscale,kViolet,kFALSE);
 465   ScaleAndConfigure(convtte,myscale,kOrange-3,kFALSE);
 466   ScaleAndConfigure(daltte,myscale,kGreen-3,kFALSE);
 467   ScaleAndConfigure(htte,myscale,kGreen+2,kFALSE);
 468   ScaleAndConfigure(wztte,myscale,kOrange-7,kFALSE);
 469
 470   return;
 471 }
 472
 473 void makeEMCElectrons(TH2F* jjemc, TH2F* bjemc, TH2F* wjemc) {
 474
 475   //Jet-Jet EMC Electrons
 476   TH1F* jjallemc = (TH1F*)jjemc->ProjectionX("jjallemc",1,1);
 477   TH1F* jjbemc = (TH1F*)jjemc->ProjectionX("jjbemc",2,2);
 478   TH1F* jjcemc = (TH1F*)jjemc->ProjectionX("jjcemc",3,3);
 479   TH1F* jjcbemc = (TH1F*)jjemc->ProjectionX("jjcbemc",4,4);
 480   TH1F* jjconvemc = (TH1F*)jjemc->ProjectionX("jjconvemc",5,5);
 481   TH1F* jjdalemc = (TH1F*)jjemc->ProjectionX("jjdalemc",6,6);
 482   TH1F* jjwzemc = (TH1F*)jjemc->ProjectionX("jjwzemc",7,7);
 483   TH1F* jjhemc = (TH1F*)jjemc->ProjectionX("jjhemc",9,9);
 484
 485   //B-Jet EMC Electrons
 486   TH1F* bjallemc = (TH1F*)bjemc->ProjectionX("bjallemc",1,1);
 487   TH1F* bjbemc = (TH1F*)bjemc->ProjectionX("bjbemc",2,2);
 488   TH1F* bjcemc = (TH1F*)bjemc->ProjectionX("bjcemc",3,3);
 489   TH1F* bjcbemc = (TH1F*)bjemc->ProjectionX("bjcbemc",4,4);
 490   TH1F* bjconvemc = (TH1F*)bjemc->ProjectionX("bjconvemc",5,5);
 491   TH1F* bjdalemc = (TH1F*)bjemc->ProjectionX("bjdalemc",6,6);
 492   TH1F* bjwzemc = (TH1F*)bjemc->ProjectionX("bjwzemc",7,7);
 493   TH1F* bjhemc = (TH1F*)bjemc->ProjectionX("bjhemc",9,9);
 494
 495   //W-Jet EMC Electrons
 496   TH1F* wjallemc = (TH1F*)wjemc->ProjectionX("wjallemc",1,1);
 497   TH1F* wjbemc = (TH1F*)wjemc->ProjectionX("wjbemc",2,2);
 498   TH1F* wjcemc = (TH1F*)wjemc->ProjectionX("wjcemc",3,3);
 499   TH1F* wjcbemc = (TH1F*)wjemc->ProjectionX("wjcbemc",4,4);
 500   TH1F* wjconvemc = (TH1F*)wjemc->ProjectionX("wjconvemc",5,5);
 501   TH1F* wjdalemc = (TH1F*)wjemc->ProjectionX("wjdalemc",6,6);
 502   TH1F* wjwzemc = (TH1F*)wjemc->ProjectionX("wjwzemc",7,7);
 503   TH1F* wjhemc = (TH1F*)wjemc->ProjectionX("wjhemc",9,9);
 504
 505   bemc = (TH1F*)bjbemc->Clone(); bemc->SetName("bemc");  //B-Jet + W-jet
 506   bemc->Add(wjbemc);
 507   TCanvas * ctempemc = new TCanvas("ctempemc","",0,0,800,800);
 508   ctempemc->Divide(2,3);
 509   ctempemc->cd(1); gPad->SetLogy(); bemc->Draw();
 510   TH1F* foobemc = (TH1F*)bemc->Clone(); foobemc->SetName("foobemc");
 511   smoothWithFit(foobemc,1e5,-3,8,40,100);
 512   TRandom rand;
 513   for(Int_t i = 8; i<= bemc->GetNbinsX(); i++) {
 514     Double_t dither = rand.Gaus(0.,foobemc->GetBinContent(i)/2.);
 515     if(dither + foobemc->GetBinContent(i) < foobemc->GetBinContent(i)/100.) dither = 0.;
 516     bemc->SetBinContent(i,foobemc->GetBinContent(i)+dither);
 517     bemc->SetBinError(i,sqrt(foobemc->GetBinContent(i)));
 518   }
 519   bemc->Draw(); foobemc->Draw("same");
 520
 521   cemc = (TH1F*)jjcemc->Clone(); cemc->SetName("cemc"); //Jet-Jet + W-jet
 522   cemc->Add(wjcemc);
 523   ctempemc->cd(2); gPad->SetLogy(); cemc->Draw();
 524   TH1F* foocemc = (TH1F*)cemc->Clone(); foocemc->SetName("foocemc");
 525   smoothWithFit(foocemc,1e5,-3,5,14,100);
 526   for(Int_t i = 5; i<= cemc->GetNbinsX(); i++) {
 527     Double_t dither = rand.Gaus(0.,foocemc->GetBinContent(i)/2.);
 528     if(dither + foocemc->GetBinContent(i) < foocemc->GetBinContent(i)/100.) dither = 0.;
 529     cemc->SetBinContent(i,foocemc->GetBinContent(i)+dither);
 530     cemc->SetBinError(i,sqrt(foocemc->GetBinContent(i)));
 531   }
 532   cemc->Draw(); foocemc->Draw("same");
 533
 534   cbemc = (TH1F*)bjcbemc->Clone(); cbemc->SetName("cbemc"); //B-Jet + W-jet
 535   cbemc->Add(wjcbemc);
 536   ctempemc->cd(3); gPad->SetLogy(); cbemc->Draw();
 537   TH1F* foocbemc = (TH1F*)cbemc->Clone(); foocbemc->SetName("foocbemc");
 538   smoothWithFit(foocbemc,1e5,-3,8,20,100);
 539   for(Int_t i = 8; i<= cbemc->GetNbinsX(); i++) {
 540     Double_t dither = rand.Gaus(0.,foocbemc->GetBinContent(i)/2.);
 541     if(dither + foocbemc->GetBinContent(i) < foocbemc->GetBinContent(i)/100.) dither = 0.;
 542     cbemc->SetBinContent(i,foocbemc->GetBinContent(i)+dither);
 543     cbemc->SetBinError(i,sqrt(foocbemc->GetBinContent(i)));
 544   }
 545   cbemc->Draw(); foocbemc->Draw("same");
 546
 547   convemc = (TH1F*)jjconvemc->Clone(); convemc->SetName("convemc"); //Jet-Jet + W-jet
 548   convemc->Add(wjconvemc);
 549   ctempemc->cd(4); gPad->SetLogy(); convemc->Draw();
 550   TH1F* fooconvemc = (TH1F*)convemc->Clone(); fooconvemc->SetName("fooconvemc");
 551   smoothWithFit(fooconvemc,1e6,-3,5,15,100);
 552   for(Int_t i = 5; i<= convemc->GetNbinsX(); i++) {
 553     Double_t dither = rand.Gaus(0.,fooconvemc->GetBinContent(i)/2.);
 554     if(dither + fooconvemc->GetBinContent(i) < fooconvemc->GetBinContent(i)/100.) dither = 0.;
 555     convemc->SetBinContent(i,fooconvemc->GetBinContent(i)+dither);
 556     convemc->SetBinError(i,sqrt(fooconvemc->GetBinContent(i)));
 557   }
 558   convemc->Draw(); fooconvemc->Draw("same");
 559
 560   dalemc = (TH1F*)jjdalemc->Clone(); dalemc->SetName("dalemc"); //Jet-Jet + W-jet
 561   dalemc->Add(wjdalemc);
 562   ctempemc->cd(5); gPad->SetLogy(); dalemc->Draw();
 563   TH1F* foodalemc = (TH1F*)dalemc->Clone(); foodalemc->SetName("foodalemc");
 564   for(Int_t i = 18; i <= dalemc->GetNbinsX(); i++) {
 565     dalemc->SetBinContent(i,0);
 566     dalemc->SetBinError(i,0);
 567     foodalemc->SetBinContent(i,0);
 568     foodalemc->SetBinError(i,0);
 569   }
 570   //  smoothWithFit(foodalemc,1e6,-3,8.,12.,100);
 571   dalemc->Draw(); //foodalemc->Draw("same");
 572
 573   hemc = (TH1F*)jjhemc->Clone(); hemc->SetName("hemc");
 574   hemc->Add(wjhemc);
 575   ctempemc->cd(6); gPad->SetLogy(); hemc->Draw();
 576   TH1F* foohemc = (TH1F*)hemc->Clone(); foohemc->SetName("foohemc");
 577   smoothWithFit(foohemc,1e5,-3,10,40,100);
 578   for(Int_t i = 10; i<= hemc->GetNbinsX(); i++) {
 579     Double_t dither = rand.Gaus(0.,foohemc->GetBinContent(i)/2.);
 580     if(dither + foohemc->GetBinContent(i) < foohemc->GetBinContent(i)/100.) dither = 0.;
 581     hemc->SetBinContent(i,foohemc->GetBinContent(i)+dither);
 582     hemc->SetBinError(i,sqrt(foohemc->GetBinContent(i)));
 583   }
 584   hemc->Draw(); foohemc->Draw("same");
 585
 586   wzemc = (TH1F*)wjwzemc->Clone(); wzemc->SetName("wzemc"); //W-jet only
 587   TCanvas* ctempwzemc = new TCanvas("ctempwzemc");
 588   ctempwzemc->cd(); gPad->SetLogy(); wzemc->Draw();
 589   TH1F* foowzemc = (TH1F*)wzemc->Clone(); foowzemc->SetName("foowzemc");
 590   TF1* fwsemc = new TF1("fwsemc","[0]*(1+exp((x-[1])/[2]))^-1",30,50);
 591   fwsemc->SetParameters(10,30,3);
 592   foowzemc->Fit(fwsemc,"R");
 593   TF1* fwzexpemc = new TF1("fwzexpemc","[0]+[1]*log(x/[2])^2",5,20);
 594   fwzexpemc->SetParameters(10,10,3);
 595   //  foowzemc->Fit(fwzexpemc,"R");
 596   /*
 597   for(Int_t i = 8; i<= wzMC->GetNbinsX(); i++) {
 598     Double_t pt = wzMC->GetBinCenter(i);
 599     if(pt < 40) wzMC->SetBinContent(i,fwzexp->Eval(pt));
 600     if(pt > 40) wzMC->SetBinContent(i,fws->Eval(pt));
 601   }
 602   */
 603   wzemc->Draw(); foowzemc->Draw("same");
 604   fwzexpemc->Draw("same");
 605
 606   //All EMC electrons is the sum of
 607   //Jet-Jet: conversions + direct charm + dalitz + misid
 608   //Bottom-Jet: direct bottom + indirect bottom
 609   //W-Jet: all (because these events are exclusive of the others)
 610   allemc = (TH1F*)wjallemc->Clone(); allemc->SetName("allemc");
 611   allemc->Add(jjconvemc); allemc->Add(jjcemc); allemc->Add(jjdalemc); allemc->Add(jjhemc);
 612   allemc->Add(bjbemc); allemc->Add(bjcbemc);
 613   sumemc = (TH1F*)wjallemc->Clone(); sumemc->SetName("sumemc");
 614   sumemc->Add(jjconvemc); sumemc->Add(jjcemc); sumemc->Add(jjdalemc);
 615   sumemc->Add(bjbemc); sumemc->Add(bjcbemc);
 616
 617   double myscale = 1.; //we already scaled them
 618   ScaleAndConfigure(allemc,myscale,kBlack,kFALSE);
 619   ScaleAndConfigure(sumemc,myscale,kBlack,kFALSE);
 620   ScaleAndConfigure(bemc,myscale,kRed,kFALSE);
 621   ScaleAndConfigure(cemc,myscale,kBlue,kFALSE);
 622   ScaleAndConfigure(cbemc,myscale,kViolet,kFALSE);
 623   ScaleAndConfigure(convemc,myscale,kOrange-3,kFALSE);
 624   ScaleAndConfigure(dalemc,myscale,kGreen-3,kFALSE);
 625   ScaleAndConfigure(hemc,myscale,kGreen+2,kFALSE);
 626   ScaleAndConfigure(wzemc,myscale,kOrange-7,kFALSE);
 627
 628   return;
 629 }
 630
 631 void makeTRKElectrons(TH2F* jjtrk, TH2F* bjtrk, TH2F* wjtrk) {
 632
 633   //Jet-Jet TRK Electrons
 634   TH1F* jjalltrk = (TH1F*)jjtrk->ProjectionX("jjalltrk",1,1);
 635   TH1F* jjbtrk = (TH1F*)jjtrk->ProjectionX("jjbtrk",2,2);
 636   TH1F* jjctrk = (TH1F*)jjtrk->ProjectionX("jjctrk",3,3);
 637   TH1F* jjcbtrk = (TH1F*)jjtrk->ProjectionX("jjcbtrk",4,4);
 638   TH1F* jjconvtrk = (TH1F*)jjtrk->ProjectionX("jjconvtrk",5,5);
 639   TH1F* jjdaltrk = (TH1F*)jjtrk->ProjectionX("jjdaltrk",6,6);
 640   TH1F* jjwztrk = (TH1F*)jjtrk->ProjectionX("jjwztrk",7,7);
 641   TH1F* jjhtrk = (TH1F*)jjtrk->ProjectionX("jjhtrk",9,9);
 642
 643   //B-Jet TRK Electrons
 644   TH1F* bjalltrk = (TH1F*)bjtrk->ProjectionX("bjalltrk",1,1);
 645   TH1F* bjbtrk = (TH1F*)bjtrk->ProjectionX("bjbtrk",2,2);
 646   TH1F* bjctrk = (TH1F*)bjtrk->ProjectionX("bjctrk",3,3);
 647   TH1F* bjcbtrk = (TH1F*)bjtrk->ProjectionX("bjcbtrk",4,4);
 648   TH1F* bjconvtrk = (TH1F*)bjtrk->ProjectionX("bjconvtrk",5,5);
 649   TH1F* bjdaltrk = (TH1F*)bjtrk->ProjectionX("bjdaltrk",6,6);
 650   TH1F* bjwztrk = (TH1F*)bjtrk->ProjectionX("bjwztrk",7,7);
 651   TH1F* bjhtrk = (TH1F*)bjtrk->ProjectionX("bjhtrk",9,9);
 652
 653   //W-Jet TRK Electrons
 654   TH1F* wjalltrk = (TH1F*)wjtrk->ProjectionX("wjalltrk",1,1);
 655   TH1F* wjbtrk = (TH1F*)wjtrk->ProjectionX("wjbtrk",2,2);
 656   TH1F* wjctrk = (TH1F*)wjtrk->ProjectionX("wjctrk",3,3);
 657   TH1F* wjcbtrk = (TH1F*)wjtrk->ProjectionX("wjcbtrk",4,4);
 658   TH1F* wjconvtrk = (TH1F*)wjtrk->ProjectionX("wjconvtrk",5,5);
 659   TH1F* wjdaltrk = (TH1F*)wjtrk->ProjectionX("wjdaltrk",6,6);
 660   TH1F* wjwztrk = (TH1F*)wjtrk->ProjectionX("wjwztrk",7,7);
 661   TH1F* wjhtrk = (TH1F*)wjtrk->ProjectionX("wjhtrk",9,9);
 662
 663   btrk = (TH1F*)bjbtrk->Clone(); btrk->SetName("btrk");  //B-Jet + W-jet
 664   btrk->Add(wjbtrk);
 665   TCanvas * ctemptrk = new TCanvas("ctemptrk","",0,0,800,800);
 666   ctemptrk->Divide(2,3);
 667   ctemptrk->cd(1); gPad->SetLogy(); btrk->Draw();
 668   TH1F* foobtrk = (TH1F*)btrk->Clone(); foobtrk->SetName("foobtrk");
 669   smoothWithFit(foobtrk,1e5,-3,8,40,100);
 670   TRandom rand2;
 671   for(Int_t i = 8; i<= btrk->GetNbinsX(); i++) {
 672     Double_t dither = rand2.Gaus(0.,foobtrk->GetBinContent(i)/2.);
 673     if(dither + foobtrk->GetBinContent(i) < foobtrk->GetBinContent(i)/100.) dither = 0.;
 674     btrk->SetBinContent(i,foobtrk->GetBinContent(i)+dither);
 675     btrk->SetBinError(i,sqrt(foobtrk->GetBinContent(i)));
 676   }
 677   btrk->Draw(); foobtrk->Draw("same");
 678
 679   ctrk = (TH1F*)jjctrk->Clone(); ctrk->SetName("ctrk"); //Jet-Jet + W-jet
 680   ctrk->Add(wjctrk);
 681   ctemptrk->cd(2); gPad->SetLogy(); ctrk->Draw();
 682   TH1F* fooctrk = (TH1F*)ctrk->Clone(); fooctrk->SetName("fooctrk");
 683   smoothWithFit(fooctrk,1e5,-3,5,20,100);
 684   for(Int_t i = 5; i<= ctrk->GetNbinsX(); i++) {
 685     Double_t dither = rand2.Gaus(0.,fooctrk->GetBinContent(i)/2.);
 686     if(dither + fooctrk->GetBinContent(i) < fooctrk->GetBinContent(i)/100.) dither = 0.;
 687     ctrk->SetBinContent(i,fooctrk->GetBinContent(i)+dither);
 688     ctrk->SetBinError(i,sqrt(fooctrk->GetBinContent(i)));
 689   }
 690   ctrk->Draw(); fooctrk->Draw("same");
 691
 692   cbtrk = (TH1F*)bjcbtrk->Clone(); cbtrk->SetName("cbtrk"); //B-Jet + W-jet
 693   cbtrk->Add(wjcbtrk);
 694   ctemptrk->cd(3); gPad->SetLogy(); cbtrk->Draw();
 695   TH1F* foocbtrk = (TH1F*)cbtrk->Clone(); foocbtrk->SetName("foocbtrk");
 696   smoothWithFit(foocbtrk,1e5,-3,5,20,100);
 697   for(Int_t i = 5; i<= cbtrk->GetNbinsX(); i++) {
 698     Double_t dither = rand2.Gaus(0.,foocbtrk->GetBinContent(i)/2.);
 699     if(dither + foocbtrk->GetBinContent(i) < foocbtrk->GetBinContent(i)/100.) dither = 0.;
 700     cbtrk->SetBinContent(i,foocbtrk->GetBinContent(i)+dither);
 701     cbtrk->SetBinError(i,sqrt(foocbtrk->GetBinContent(i)));
 702   }
 703   cbtrk->Draw(); foocbtrk->Draw("same");
 704
 705   convtrk = (TH1F*)jjconvtrk->Clone(); convtrk->SetName("convtrk"); //Jet-Jet + W-jet
 706   convtrk->Add(wjconvtrk);
 707   ctemptrk->cd(4); gPad->SetLogy(); convtrk->Draw();
 708   TH1F* fooconvtrk = (TH1F*)convtrk->Clone(); fooconvtrk->SetName("fooconvtrk");
 709   smoothWithFit(fooconvtrk,1e6,-3,5,20,100);
 710   for(Int_t i = 5; i<= convtrk->GetNbinsX(); i++) {
 711     Double_t dither = rand2.Gaus(0.,fooconvtrk->GetBinContent(i)/2.);
 712     if(dither + fooconvtrk->GetBinContent(i) < fooconvtrk->GetBinContent(i)/100.) dither = 0.;
 713     convtrk->SetBinContent(i,fooconvtrk->GetBinContent(i)+dither);
 714     convtrk->SetBinError(i,sqrt(fooconvtrk->GetBinContent(i)));
 715   }
 716   convtrk->Draw(); fooconvtrk->Draw("same");
 717
 718   daltrk = (TH1F*)jjdaltrk->Clone(); daltrk->SetName("daltrk"); //Jet-Jet + W-jet
 719   daltrk->Add(wjdaltrk);
 720   ctemptrk->cd(5); gPad->SetLogy(); daltrk->Draw();
 721   TH1F* foodaltrk = (TH1F*)daltrk->Clone(); foodaltrk->SetName("foodaltrk");
 722   smoothWithFit(foodaltrk,1e6,-3,5.,40.,100);
 723   for(Int_t i = 5; i<= daltrk->GetNbinsX(); i++) {
 724     Double_t dither = rand2.Gaus(0.,foodaltrk->GetBinContent(i)/2.);
 725     if(dither + foodaltrk->GetBinContent(i) < foodaltrk->GetBinContent(i)/100.) dither = 0.;
 726     daltrk->SetBinContent(i,foodaltrk->GetBinContent(i)+dither);
 727     daltrk->SetBinError(i,sqrt(foodaltrk->GetBinContent(i)));
 728   }
 729   daltrk->Draw(); foodaltrk->Draw("same");
 730
 731   htrk = (TH1F*)jjhtrk->Clone(); htrk->SetName("htrk");
 732   htrk->Add(wjhtrk);
 733   ctemptrk->cd(6); gPad->SetLogy(); htrk->Draw();
 734   TH1F* foohtrk = (TH1F*)htrk->Clone(); foohtrk->SetName("foohtrk");
 735   smoothWithFit(foohtrk,1e5,-3,5,80,100);
 736   for(Int_t i = 5; i<= htrk->GetNbinsX(); i++) {
 737     Double_t dither = rand2.Gaus(0.,foohtrk->GetBinContent(i)/2.);
 738     if(dither + foohtrk->GetBinContent(i) < foohtrk->GetBinContent(i)/100.) dither = 0.;
 739     htrk->SetBinContent(i,foohtrk->GetBinContent(i)+dither);
 740     htrk->SetBinError(i,sqrt(foohtrk->GetBinContent(i)));
 741   }
 742   htrk->Draw(); foohtrk->Draw("same");
 743
 744   wztrk = (TH1F*)wjwztrk->Clone(); wztrk->SetName("wztrk"); //W-jet only
 745   TCanvas* ctempwztrk = new TCanvas("ctempwztrk");
 746   ctempwztrk->cd(); gPad->SetLogy(); wztrk->Draw();
 747   TH1F* foowztrk = (TH1F*)wztrk->Clone(); foowztrk->SetName("foowztrk");
 748   TF1* fwstrk = new TF1("fwstrk","[0]*(1+exp((x-[1])/[2]))^-1",30,50);
 749   fwstrk->SetParameters(10,30,3);
 750   foowztrk->Fit(fwstrk,"R");
 751   TF1* fwzexptrk = new TF1("fwzexptrk","[0]+[1]*log(x/[2])^2",8,40);
 752   fwzexptrk->SetParameters(10,10,3);
 753   //  foowztrk->Fit(fwzexptrk,"R");
 754   /*
 755   for(Int_t i = 8; i<= wzMC->GetNbinsX(); i++) {
 756     Double_t pt = wzMC->GetBinCenter(i);
 757     if(pt < 40) wzMC->SetBinContent(i,fwzexp->Eval(pt));
 758     if(pt > 40) wzMC->SetBinContent(i,fws->Eval(pt));
 759   }
 760   */
 761   wztrk->Draw(); foowztrk->Draw("same");
 762
 763   //All TRK electrons is the sum of
 764   //Jet-Jet: conversions + direct charm + dalitz + misid
 765   //Bottom-Jet: direct bottom + indirect bottom
 766   //W-Jet: all (because these events are exclusive of the others)
 767   alltrk = (TH1F*)wjalltrk->Clone(); alltrk->SetName("alltrk");
 768   alltrk->Add(jjconvtrk); alltrk->Add(jjctrk); alltrk->Add(jjdaltrk); alltrk->Add(jjhtrk);
 769   alltrk->Add(bjbtrk); alltrk->Add(bjcbtrk);
 770   sumtrk = (TH1F*)wjalltrk->Clone(); sumtrk->SetName("sumtrk");
 771   sumtrk->Add(jjconvtrk); sumtrk->Add(jjctrk); sumtrk->Add(jjdaltrk);
 772   sumtrk->Add(bjbtrk); sumtrk->Add(bjcbtrk);
 773
 774   double myscale = 1.; //we already scaled them
 775   ScaleAndConfigure(alltrk,myscale,kBlack,kFALSE);
 776   ScaleAndConfigure(sumtrk,myscale,kBlack,kFALSE);
 777   ScaleAndConfigure(btrk,myscale,kRed,kFALSE);
 778   ScaleAndConfigure(ctrk,myscale,kBlue,kFALSE);
 779   ScaleAndConfigure(cbtrk,myscale,kViolet,kFALSE);
 780   ScaleAndConfigure(convtrk,myscale,kOrange-3,kFALSE);
 781   ScaleAndConfigure(daltrk,myscale,kGreen-3,kFALSE);
 782   ScaleAndConfigure(htrk,myscale,kGreen+2,kFALSE);
 783   ScaleAndConfigure(wztrk,myscale,kOrange-7,kFALSE);
 784
 785   return;
 786 }
 787
 788 void smoothWithFit(TH1F* hist, Double_t p0, Double_t p1, Double_t min, Double_t max,Double_t remax) {
 789
 790   fpow = new TF1("fpow","[0]*pow(x,[1])",min,max);
 791   fpow->SetParameters(p0,p1);
 792   hist->Fit(fpow,"R");
 793   for(Int_t i = (Int_t)min; i <= (Int_t)remax; i++) {
 794     Double_t pt = hist->GetBinCenter(i);
 795     Double_t val = fpow->Eval(pt);
 796     hist->SetBinContent(i,val);
 797     hist->SetBinError(i,sqrt(val));
 798   }
 799
 800   return;
 801 }
 802
 803