TFile *pfile = 0;
Int_t rebin = 2;

struct chdef {
  char *name;
  char *name2;
  int color;
  int style;
  char *label;
} decays[] = { 	       
  {"tot",0,1,0,"Total"},
  {"B0KpiTot",0,2,3003,"B^{0}#rightarrow K #pi+cc"},
  {"B0Kpi",0,2,3004,"B^{0}#rightarrow K #pi"},
  {"antiB0Kpi",0,2,3005,"#bar{B}^{0} #rightarrow K #pi"},
  {"BsKK",0,8,1001,"B_{s} #rightarrow KK"},
  {"BsKpiTot",0,6,3003,"B_{s} #rightarrow K #pi+cc"},
  {"BsKpi",0,6,3004,"B_{s} #rightarrow K #pi"},
  {"antiBsKpi",0,6,3005,"#bar{B}^{0} #rightarrow K#pi"},
  {"B0pipi",0,1,1001,"B^{0}#rightarrow #pi #pi"},
  {"B0KK",0,4,1001,"B^{0} #rightarrow KK"},
  {"Bspipi",0,5,1001,"B_{s} #rightarrow #pi #pi"},
  {"LbppiTot",0,9,3003,"#Lambda_{b}^{0}#rightarrow p #pi+cc"},
  {"Lbppi",0,9,3004,"#Lambda_{b}^{0}#rightarrow p #pi"},
  {"antiLbppi",0,9,3005,"#bar{#Lambda}_{b}^{0}#rightarrow p #pi"},
  {"LbpKTot",0,11,3003,"#Lambda_{b}^{0}#rightarrow p K+cc"},
  {"LbpK",0,11,3004,"#Lambda_{b}^{0}#rightarrow p K"},
  {"antiLbpK",0,11,3005,"#bar{#Lambda}_{b}^{0}#rightarrow p K"}
  {"ARGUS",0,14,1001,"Bkg ARGUS"},
  {"exp","bkg",13,1001,"Bkg Exp"},
};
  
Int_t nchannels = 19;

Int_t ExtractProjections(const char *name, TH1 **histo_proj) {
  Int_t rchannel = 0;
  for (int i=0;i<nchannels;++i) {
    histo_proj[i] = (TH1*) pfile->Get(Form("projs/histo_%s_%s",name,decays[i].name));
    if (!histo_proj[i]) {
      //cerr << name << " projection for " << decays[i].name << " don't found" << endl;
      continue;
    }

    rchannel += 1;
    histo_proj[i]->UseCurrentStyle();
    histo_proj[i]->SetLineColor(decays[i].color);
    histo_proj[i]->SetFillColor(decays[i].color);
    histo_proj[i]->SetFillStyle(decays[i].style);

    if (rebin>1) histo_proj[i]->Rebin(rebin);
  }

  return rchannel;
}

void DoDifferences(TPad *&pad_diff, const char *name,
		   TH1 *histo_data, TH1 *histo_projtot, 
		   TH1 * &diff_histo) {
  pad_diff = new TPad(Form("pad_diff_%s",name),"",0,0,1.,.2);
  pad_diff->Draw();
  pad_diff->SetBottomMargin(.05);
  pad_diff->SetTopMargin(.0);
  pad_diff->cd();
  
  diff_histo = (TH1*) histo_data->Clone(Form("diff_%s",name));
  diff_histo->SetTitle("");
  diff_histo->SetYTitle("#Delta/#sigma");
  diff_histo->GetYaxis()->SetTitleSize(.3);
  diff_histo->GetYaxis()->SetTitleOffset(.15);
  diff_histo->GetYaxis()->CenterTitle();
  diff_histo->GetYaxis()->SetLabelSize(0);
  diff_histo->SetXTitle("");
  diff_histo->GetXaxis()->SetLabelSize(0);
  diff_histo->SetFillStyle(1001);
  diff_histo->SetFillColor(kRed);
  diff_histo->SetMinimum(-5);
  diff_histo->SetMaximum(5);
  for (Int_t i=1;i<=diff_histo->GetNbinsX();++i) {
    if (diff_histo->GetBinContent(i)==0) continue;

    Int_t oi = histo_projtot->GetXaxis()->FindBin(diff_histo->GetBinCenter(i));
    diff_histo->SetBinContent(i,(diff_histo->GetBinContent(i)-histo_projtot->GetBinContent(oi))
			      /diff_histo->GetBinError(i));
  }
  diff_histo->Draw();
  diff_histo->GetYaxis()->SetNdivisions(10);
  diff_histo->GetXaxis()->SetNdivisions(3);
  pad_diff->SetGrid(1,1);
}

// mass projections
void mass_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do Mass projections" << endl;

  TCanvas *canvas_mass = new TCanvas("canvas_mass","Mass projection",1);
  TPad *pad_main = new TPad("pad_main_mass","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_masspipi");
  histo_data->UseCurrentStyle();
  histo_data->SetXTitle("Mass #pi#pi [GeV/c^{2}]");
  histo_data->SetYTitle(Form("Events per %2.2f Mev/c^{2}",histo_data->GetBinWidth(1)*1e3));
  histo_data->Draw("e");

  // various decay channels
  //TH1 *histo_mass[nchannels];
  TH1 **histo_mass = new TH1*[nchannels];
  int rchannel = ExtractProjections("Mass",histo_mass);

  TLegend *leg = new TLegend(.6,.4,.75,.5+.1*rchannel);
  leg->AddEntry(histo_data,"Data","pl");
  

  Double_t scale_factor = histo_data->Integral()/histo_mass[0]->Integral()*
    //histo_mass[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_mass[0]->GetBinWidth(1);

  for (int i=0;i<nchannels;++i) {
    if (!histo_mass[i]) continue;
    cout << "Mass projection for " << decays[i].name << endl;

    // scale the histogram
    histo_mass[i]->Scale(scale_factor);
    histo_mass[i]->Draw("c same");

    leg->AddEntry(histo_mass[i],decays[i].label,"lf");
  }

  leg->Draw();

  canvas_mass->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"IDPtot",histo_data,histo_mass[0],diff_histo);
}


// ID1 projections
void ID1_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do ID1 projections" << endl;

  TCanvas *canvas_mass = new TCanvas("canvas_ID1","ID1 projection",1);
  TPad *pad_main = new TPad("pad_main_ID1","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_id1");
  histo_data->SetYTitle(Form("Events per %4.2f Interval",
			     histo_data->GetBinWidth(1)));
  histo_data->SetXTitle("ID1");
  histo_data->UseCurrentStyle();
  histo_data->Draw("e");

  
  // various decay channels
  //TH1 *histo_proj[nchannels];
  TH1 **histo_proj = new TH1*[nchannels];
  int rchannel = ExtractProjections("ID1",histo_proj);

  TLegend *leg = new TLegend(.2,.75-.1*rchannel,.42,.85);
  leg->AddEntry(histo_data,"Data","pl");

  Double_t scale_factor = histo_data->Integral()/histo_proj[0]->Integral()*
    //histo_proj[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_proj[0]->GetBinWidth(1);

  for (int i=0;i<nchannels;++i) {
    if (!histo_proj[i]) continue;

    cout << "ID1 projection for " << decays[i].name << endl;

    // scale the histogram
    histo_proj[i]->Scale(scale_factor);
    histo_proj[i]->Draw("c same");

    leg->AddEntry(histo_proj[i],decays[i].label,"lf");
  }

  leg->Draw();

  canvas_mass->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"ID1",histo_data,histo_proj[0],diff_histo);
}


// ID2 projections
void ID2_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do ID2 projections" << endl;

  TCanvas *canvas_mass = new TCanvas("canvas_ID2","ID2 projection",1);
  TPad *pad_main = new TPad("pad_main_ID2","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_id2");
  histo_data->SetYTitle(Form("Events per %4.2f Interval",
			     histo_data->GetBinWidth(1)));
  histo_data->SetXTitle("ID2");
  histo_data->UseCurrentStyle();
  histo_data->Draw("e");

  
  // various decay channels
  //TH1 *histo_proj[nchannels];
  TH1 **histo_proj = new TH1*[nchannels];
  int rchannel = ExtractProjections("ID2",histo_proj);

  TLegend *leg = new TLegend(.2,.75-.1*rchannel,.42,.85);
  leg->AddEntry(histo_data,"Data","pl");

  // scale the histogram
  Double_t scale_factor = histo_data->Integral()/histo_proj[0]->Integral()*
    //histo_proj[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_proj[0]->GetBinWidth(1);

  for (int i=0;i<nchannels;++i) {
    if (!histo_proj[i]) continue;

    cout << "ID2 projection for " << decays[i].name << endl;

    histo_proj[i]->Scale(scale_factor);
    histo_proj[i]->Draw("c same");

    leg->AddEntry(histo_proj[i],decays[i].label,"lf");
  }

  leg->Draw();

  canvas_mass->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"ID2",histo_data,histo_proj[0],diff_histo);
}


// IDsum projections
void IDsum_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do IDsum projections" << endl;

  TCanvas *canvas_mass = new TCanvas("canvas_IDsum","IDsum projection",1);
  TPad *pad_main = new TPad("pad_main_IDsum","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_id_sum");
  histo_data->SetYTitle(Form("Events per %4.2f Interval",
			     histo_data->GetBinWidth(1)));
  histo_data->SetXTitle("ID1+ID2");
  histo_data->UseCurrentStyle();
  histo_data->Draw("e");

  
  // various decay channels
  //  TH1 *histo_proj[nchannels];
  TH1 **histo_proj = new TH1*[nchannels];
  int rchannel = ExtractProjections("IDsum",histo_proj);

  TLegend *leg = new TLegend(.2,.75-.1*rchannel,.42,.85);
  leg->AddEntry(histo_data,"Data","pl");

  // scale the histogram
  Double_t scale_factor = histo_data->Integral()/histo_proj[0]->Integral()*
    //histo_proj[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_proj[0]->GetBinWidth(1);

  for (int i=0;i<nchannels;++i) {
    if (!histo_proj[i]) continue;

    cout << "IDsum projection for " << decays[i].name << endl;

    histo_proj[i]->Scale(scale_factor);
    histo_proj[i]->Draw("c same");

    leg->AddEntry(histo_proj[i],decays[i].label,"lf");
  }

  leg->Draw();

  canvas_mass->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"IDsum",histo_data,histo_proj[0],diff_histo);
}


// IDdiff projections
void IDdiff_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do IDdiff projections" << endl;

  TCanvas *canvas_mass = new TCanvas("canvas_IDdiff","IDdiff projection",1);
  TPad *pad_main = new TPad("pad_main_IDsum","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_id_diff");
  histo_data->SetTitle("ID1-ID2");
  histo_data->SetYTitle(Form("Events per %4.2f Interval",
			     histo_data->GetBinWidth(1)));
  histo_data->SetXTitle("ID1-ID2");
  histo_data->GetXaxis()->CenterTitle();
  histo_data->UseCurrentStyle();
  histo_data->Draw("e");

  
  // various decay channels
  //TH1 *histo_proj[nchannels];
  TH1 **histo_proj = new TH1*[nchannels];
  int rchannel = ExtractProjections("IDdiff",histo_proj);

  TLegend *leg = new TLegend(.2,.75-.1*rchannel,.42,.85);
  leg->AddEntry(histo_data,"Data","pl");

  // scale the histogram
  Double_t scale_factor = histo_data->Integral()/histo_proj[0]->Integral()*
    //histo_proj[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_proj[0]->GetBinWidth(1);

  for (int i=0;i<nchannels;++i) {
    if (!histo_proj[i]) continue;

    cout << "IDdiff projection for " << decays[i].name << endl;

    histo_proj[i]->Scale(scale_factor);
    histo_proj[i]->Draw("c same");

    leg->AddEntry(histo_proj[i],decays[i].label,"lf");
  }

  leg->Draw();

  canvas_mass->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"IDdiff",histo_data,histo_proj[0],diff_histo);
}

// Alpha projections
void Alpha_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do Alpha projections" << endl;

  TCanvas *canvas_mass = new TCanvas("canvas_Alpha","Alpha projection",1);
  TPad *pad_main = new TPad("pad_main_Alpha","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_alpha");
  histo_data->SetYTitle(Form("Events per %4.2f interval",
			     histo_data->GetBinWidth(1)));
  histo_data->SetXTitle("#alpha");
  histo_data->UseCurrentStyle();
  histo_data->Draw("e");

  TLegend *leg = new TLegend(.6,.4,.85,.8);
  leg->AddEntry(histo_data,"Data","pl");
  
  // various decay channels
  //TH1 *histo_proj[nchannels];
  TH1 **histo_proj = new TH1*[nchannels];
  int rchannel = ExtractProjections("Alpha",histo_proj);

  // scale the histogram
  Double_t scale_factor = histo_data->Integral()/histo_proj[0]->Integral()*
    //histo_proj[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_proj[0]->GetBinWidth(1);

  for (int i=0;i<nchannels;++i) {
    if (!histo_proj[i]) continue;

    cout << "Alpha projection for " << decays[i].name << endl;

    histo_proj[i]->Scale(scale_factor);
    histo_proj[i]->Draw("c same");

    leg->AddEntry(histo_proj[i],decays[i].label,"lf");
  }

  //leg->Draw();

  canvas_mass->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"IDAlpha",histo_data,histo_proj[0],diff_histo);
}


// Ptot projections
void Ptot_proj() {
  if (!pfile) {
    cerr << "No file found" << endl;
    return;
  }
  cout << "Do Ptot projections" << endl;

  TCanvas *canvas_ptot = new TCanvas("canvas_Ptot","Ptot projection",1);
  TPad *pad_main = new TPad("pad_main_Ptot","",0,.2,1,1);
  pad_main->Draw();
  pad_main->cd();

  TH1 *histo_data = (TH1*) pfile->Get("data/histo_ptot");
  histo_data->UseCurrentStyle();
  histo_data->SetYTitle(Form("Events per %4.2f MeV/c",
			     histo_data->GetBinWidth(1)*1e3));
  histo_data->SetXTitle("P_{tot} [GeV/c]");
  histo_data->Draw("e");

  
  // various decay channels
  //TH1 *histo_proj[nchannels];
  TH1 **histo_proj = new TH1*[nchannels];
  int rchannel = ExtractProjections("Ptot",histo_proj);

  TLegend *leg = new TLegend(.65,.75-.1*rchannel,.85,.85);
  leg->AddEntry(histo_data,"Data","pl");

  TH1 **histo_pdf = new TH1*[nchannels];
  int npdfs = 0;
  for (int i=0;i<nchannels;++i) {
    // take pdf for this channel
    TH2 *param2d = (TH2*) pfile->Get(Form("pdfs/histo_aptot_%s",decays[i].name));

    if (!param2d && decays[i].name2) { // try with the alias
      param2d = (TH2*) pfile->Get(Form("pdfs/histo_aptot_%s",decays[i].name2));
    }

    if (!param2d) {
      histo_pdf[i] = 0;
      cerr << "Ptot pdf for " << decays[i].name << " don't found" << endl;
      continue;
    }
    
    histo_pdf[i] = param2d->ProjectionY();
    npdfs += 1;
  }

  TCanvas *canvas_ptot_chn = 0;

  if (npdfs>1)
    canvas_ptot_chn = new TCanvas("canvas_Ptot_chn",
				  "Ptot proj. chn",1);
  int nrow = npdfs/3;
  if (npdfs%3) nrow += 1;
  canvas_ptot_chn->Divide(3,nrow);

  // scale the histogram
  Double_t scale_factor = histo_data->Integral()/histo_proj[0]->Integral()*
    //histo_proj[0]->GetNbinsX()/histo_data->GetNbinsX();
    histo_data->GetBinWidth(1)/histo_proj[0]->GetBinWidth(1);

  nrow = 0;
  for (int i=0;i<nchannels;++i) {
    if (histo_proj[i]) { // draw projection
      pad_main->cd();
      cout << "Ptot projection for " << decays[i].name << endl;

      histo_proj[i]->Scale(scale_factor);
      histo_proj[i]->Draw("c same");

      leg->AddEntry(histo_proj[i],decays[i].label,"lf");
      pad_main->Update();
        
      if (histo_pdf[i]) {
	++nrow;
	Double_t scale = histo_proj[i]->Integral()/histo_pdf[i]->Integral()*
	  //histo_pdf[i]->GetNbinsX()/histo_proj[i]->GetNbinsX();
	  histo_proj[i]->GetBinWidth(1)/histo_pdf[i]->GetBinWidth(1);
	cout << "pdf check for " << decays[i].name << endl;
	
	canvas_ptot_chn->cd(nrow);
	histo_proj[i]->Draw();
	histo_pdf[i]->Scale(scale);
	histo_pdf[i]->Draw("same");
	canvas_ptot_chn->cd(nrow)->Update();
      }
    }
  }
  if (nrow==0) canvas_ptot_chn->Close();

  pad_main->cd();
  histo_data->Draw("e same");
  leg->Draw();

  canvas_ptot->cd();
  TPad *pad_diff;
  TH1* diff_histo;
  DoDifferences(pad_diff,"IDPtot",histo_data,histo_proj[0],diff_histo);
}

void draw_proj(const char *pfname = "projections.root") {
  gROOT->Macro("../bhhstyle.C");

  pfile = new TFile(pfname);

  TControlBar *macrobar = new TControlBar("vertical","Do projections plot",10,10);
  macrobar->AddButton("Mass #pi#pi","mass_proj()","Draw Mass projections");
  macrobar->AddButton("Alpha","Alpha_proj()","Draw Alpha projections");
  macrobar->AddButton("Ptot","Ptot_proj()","Draw Ptot projections");
  macrobar->AddButton("ID1","ID1_proj()","Draw ID1 projections");
  macrobar->AddButton("ID2","ID2_proj()","Draw ID2 projections");
  macrobar->AddButton("IDsum","IDsum_proj()","Draw ID Diff projections");
  macrobar->AddButton("IDdiff","IDdiff_proj()","Draw ID Sum projections");
  macrobar->AddButton("TBrowser","new TBrowser()","Open a TBrowser object");
  macrobar->Show();
}