// #define ZOOMED

void DrawMass(TString nSt="") {

  TStyle *plain = new TStyle("Plain","Plain Style(no colors/fill areas)");
  plain->SetCanvasBorderMode(0);
  plain->SetPadBorderMode(0);
  plain->SetPadColor(0);
  plain->SetCanvasColor(0);
  plain->SetTitleColor(0);
  plain->SetStatColor(0);
  plain->SetOptStat(0);
//  plain->SetOptFit(11);
  plain->SetOptFit(-1);
  plain->SetStatW(0.19);
  plain->SetStatH(0.7);
  plain->SetOptTitle(0);
  plain->SetPadLeftMargin(0.15);
  plain->SetPadRightMargin(0.04);
  plain->SetPadTopMargin(0.08);
  plain->SetPadBottomMargin(0.18);
  plain->cd();

  TString name1 = "Masses.root";

  TFile *f1 = new TFile(name1.Data(),"read");
  TDirectory *dir1 = (TDirectory*) f1;//->GetDirectory("KFParticlesFinder/Particles/K0s/Parameters/");

  const int NParameters = 6;

  TH1D *hPart[NParameters];
  TF1  *fPart[NParameters];

  string sPart[NParameters];
  string sPartFit[NParameters];

  string ParNames[NParameters] = {"ks",
                                  "lambda",
                                  "ksi-",
                                  "omega-",
                                  "sigma*-",
                                  "sigma*+"};

  string AxisName[NParameters];

  AxisName[0] = "K_{s}^{0} m_{inv} [GeV/c^{2}]";
  AxisName[1] = "#Lambda m_{inv} [GeV/c^{2}]";
  AxisName[2] = "#Xi^{-} m_{inv} [GeV/c^{2}]";
  AxisName[3] = "#Omega^{-} m_{inv} [GeV/c^{2}]";
  AxisName[4] = "#Sigma^{*-} m_{inv} [GeV/c^{2}]";
  AxisName[5] = "#Sigma^{*+} m_{inv} [GeV/c^{2}]";
  
  for(int i=0; i<NParameters; i++)
  {
    sPart[i]    = ParNames[i];
    sPartFit[i] = ParNames[i]+"GaussFit";
  }

  struct Bins 
  {
    float binStart, binEnd;
  };

  Bins fit_range[NParameters] = { {0.43, 0.6}, {1.09, 1.14}, {1.28, 1.38}, {1.5, 1.7}, {1.32, 1.45}, {1.32, 1.45} };
#ifdef ZOOMED
  Bins zoom_range[NParameters] = { {0.48, 0.52}, {1.1, 1.13}, {1.28, 1.38}, {1.5, 1.7}, {1.32, 1.45}, {1.32, 1.45} };
#endif
  
  for(int i=0; i<NParameters; i++)
  {
    hPart[i] = (TH1D*) f1->Get(sPart[i].data());
    hPart[i] -> SetLineColor(1);
    fPart[i] = new TF1(sPartFit[i].data(),"gaus(0)+pol2(3)",fit_range[i].binStart,fit_range[i].binEnd);
    fPart[i] -> SetNpx(10000);
    fPart[i] -> SetParNames("Const","Mass","Sigma","a0","a1","a2");
    fPart[i] -> SetLineColor(2);
    fPart[i] -> SetParameter(1,(fit_range[i].binStart + fit_range[i].binEnd)/2);
    fPart[i] -> SetParameter(2,(fit_range[i].binEnd - fit_range[i].binStart)/4);
    fPart[i] -> SetParLimits(2,0,1);
    if(i==0)
    {
      fPart[i] -> SetParameter(0, 100000);
      fPart[i] -> SetParameter(1, 0.5);
      fPart[i] -> SetParameter(2, 0.0027);
//      fPart[i] -> SetParameter(3, -3.132e5);
//      fPart[i] -> SetParameter(4, 5.534e5);
//      fPart[i] -> SetParameter(5, -2.438e5);
    }
    if(i==1)
    {
      fPart[i] -> SetParameter(0, 100000);
      fPart[i] -> SetParameter(1, 1.115);
      fPart[i] -> SetParameter(2, 0.002);
    }
    if(i==0)
    {
      fPart[i] -> SetParameter(0, 1000);
      fPart[i] -> SetParameter(1, 0.5);
      fPart[i] -> SetParameter(2, 0.0027);
//      fPart[i] -> SetParameter(3, -3.132e5);
//      fPart[i] -> SetParameter(4, 5.534e5);
//      fPart[i] -> SetParameter(5, -2.438e5);
    }
    if(i==1)
    {
      fPart[i] -> SetParameter(0, 1000);
      fPart[i] -> SetParameter(1, 1.115);
      fPart[i] -> SetParameter(2, 0.002);
    }
    if(i==4)
    {
      fPart[i] -> SetParameter(0, 20);
      fPart[i] -> SetParameter(1, 1.385);
      fPart[i] -> SetParameter(2, 0.005);
//       fPart[i] -> SetParameter(3, 2.26517e+04);
//       fPart[i] -> SetParameter(4, -3.26356e+04);
//       fPart[i] -> SetParameter(5,  1.18000e+04);
    }
    if(i==5)
    {
      fPart[i] -> SetParameter(0, 20);
      fPart[i] -> SetParameter(1, 1.385);
      fPart[i] -> SetParameter(2, 0.005);
    }
    
    hPart[i] -> Fit(sPartFit[i].data(),"","",fit_range[i].binStart,fit_range[i].binEnd);
    for(int j=3; j<6; j++)
      fPart[i] -> FixParameter(j,fPart[i]->GetParameter(j));
    hPart[i] -> Fit(sPartFit[i].data(),"","",fit_range[i].binStart,fit_range[i].binEnd);
#ifdef ZOOMED
    hPart[i]->GetXaxis()->SetRangeUser(zoom_range[i].binStart,zoom_range[i].binEnd);
#endif

  }

  TCanvas *c1;
  c1 = new TCanvas("c1","Particles",0,0,1200,1000);
  c1 -> UseCurrentStyle();
  c1 -> Divide(2,3);

  int NPad=0;

  for(int i=0; i<NParameters; i++)
  {
    NPad++;
    c1->cd(NPad);
    //hPart[i] -> UseCurrentStyle();
    hPart[i] -> GetXaxis()-> SetTitle(AxisName[i].data());
    hPart[i] -> GetXaxis()-> SetTitleSize(0.08);
    hPart[i] -> GetXaxis()-> SetLabelSize(0.07);
    hPart[i] -> GetXaxis()-> SetTitleOffset(1);
    hPart[i] -> GetXaxis()-> SetNdivisions(5,5,0);
    hPart[i] -> GetYaxis()-> SetTitle("Entries");
    hPart[i] -> GetYaxis()-> SetTitleSize(0.08);
    hPart[i] -> GetYaxis()-> SetLabelSize(0.07);
    hPart[i] -> GetYaxis()-> SetTitleOffset(1);
    hPart[i] -> GetYaxis()-> SetNdivisions(3,5,0);
    hPart[i] -> Draw();
    if(i<3)fPart[i] -> Draw("same");
  }

#ifdef ZOOMED
  c1->SaveAs("Mass_zoomed.pdf");
  c1->SaveAs("Mass_zoomed.png");
#else
  c1->SaveAs("Mass.pdf");
  c1->SaveAs("Mass.png");
#endif
  
  TF1 *fSInt[NParameters];
  TF1 *fBInt[NParameters];
  for(int i=0; i<NParameters; i++)
  {
    float start = fPart[i]->GetParameter(1) - 2*fPart[i]->GetParameter(2);
    float end   = fPart[i]->GetParameter(1) + 2*fPart[i]->GetParameter(2);

    sPartFit[i] += "int";
    fSInt[i] = new TF1(sPartFit[i].data(),"gaus",fit_range[i].binStart,fit_range[i].binEnd);
    for(int iPar=0; iPar<3; iPar++)
      fSInt[i] -> SetParameter(iPar, fPart[i] -> GetParameter(iPar));

    sPartFit[i] += "BG";
    fBInt[i] = new TF1(sPartFit[i].data(),"pol2",fit_range[i].binStart,fit_range[i].binEnd);
    for(int iPar=3; iPar<6; iPar++)
      fBInt[i] -> SetParameter(iPar-3, fPart[i] -> GetParameter(iPar));

    float S = fSInt[i] -> Integral(start,end);
    float B = fBInt[i] -> Integral(start,end);
    const float inorm = ( hPart[i]->GetXaxis()->GetXmax() - hPart[i]->GetXaxis()->GetXmin() )*hPart[i]->GetXaxis()->GetNbins();

    cout << ParNames[i] << endl;
    std::cout << "S    " << i << "  " << S*inorm << std::endl;
    std::cout << "B    " << i << "  " << B*inorm << std::endl;
    std::cout << "S/B  " << i << "  " << S/B << std::endl;
  }
}