class  MyFitFun1 {
public:
  double Evaluate(double *xPtr, double *parPtr) {
    double fitVal = 0.0;
    double xVal = *xPtr;


    double mag1   = parPtr[0];
    double center1 = parPtr[1];
    double gamma1  = parPtr[2];
    double bw1 = mag1/(pow(xVal*xVal - center1*center1,2) + center1*center1*gamma1*gamma1);

    double mag2   = parPtr[3];
    double center2 = parPtr[4];
    double gamma2  = parPtr[5];
    double bw2 = mag2/(pow(xVal*xVal - center2*center2,2) + center2*center2*gamma2*gamma2);

    double bg = 1*fabs(parPtr[6]) + 0*parPtr[7]*xVal;
    fitVal =  bw1 + bw2 + 0*bg;

    return fitVal;
  }
};

class  MyFitFun2 {
public:
  double Evaluate(double *xPtr, double *parPtr) {
    double fitVal = 0.0;
    double xVal = *xPtr;


    double mag1   = fabs(parPtr[0]);
    double center1 = parPtr[1];
    double gamma1  = parPtr[2];
    double bw1 = mag1/(pow(xVal*xVal - center1*center1,2) + center1*center1*gamma1*gamma1);

    double mag2   = fabs(parPtr[3]);
    double center2 = parPtr[4];
    double gamma2  = parPtr[5];
    double bw2 = mag2/(pow(xVal*xVal - center2*center2,2) + center2*center2*gamma2*gamma2);

    double mag3   = fabs(parPtr[6]);
    double center3 = fabs(parPtr[7]);
    double gamma3  = fabs(parPtr[8]);
    double bw3 = mag3/(pow(xVal*xVal - center3*center3,2) + center3*center3*gamma3*gamma3);

    double bg = 1*fabs(parPtr[9]) + 1*parPtr[10]*xVal;
    fitVal =  bw1 + bw2 + bw3 + 0*bg;

    return fitVal;
  }
};




void makePics(int iShow, int pdgLim){

  gROOT->Reset();
  gROOT->SetStyle("Plain");
  gStyle->SetOptStat(0);

  TFile *fFF = new TFile("./outFitToRun/ff.root");
  //TFile *fA  = new TFile("./A_V1.root");
  TFile *fA  = new TFile("./AC_V1.root");

  TFile *fT  = new TFile("./AB_V1_mcThrownStage3.root");
  TH1D *hT; fT->GetObject("hMassKpKmPi0Course3",hT);

  TFile *fS  = new TFile("./AC_V1_mcSeenStage3.root");
  TH1D *hS; fS->GetObject("hMassKpKmPi0Course3",hS);
  hS->Sumw2();
  hT->Sumw2();
  TH1D * hEff = (TH1D *) hT->Clone("hEff");
  hEff->Divide(hS,hT);
  
  if (iShow == 91) hT->Draw();
  if (iShow == 92) hS->Draw();
  if (iShow == 93) hEff->Draw("e");
  //return;

  TH1D *hReal1; fA->GetObject("hMassKpKmPi0Course1",hReal1);
  TH1D *hReal3; fA->GetObject("hMassKpKmPi0Course3",hReal3);

  TH2D *hFitFrac; fFF->GetObject("hFitFrac",hFitFrac);

  TH1D * hJ0 = (TH1D*)hFitFrac->ProjectionX("hJ0",1,2);
  TH1D * hJ1 = (TH1D*)hFitFrac->ProjectionX("hJ1",3,10);
  //TH1D * hJ1 = (TH1D*)hFitFrac->ProjectionX("hJ1",2,4);
  //TH1D * hJ1 = (TH1D*)hFitFrac->ProjectionX("hJ1",5,10);

  //hReal3->SetBinError(5,0);
  //hReal3->SetBinError(6,0);
  //hReal3->SetBinContent(5,0);
  //hReal3->SetBinContent(6,0);

  //hJ0->SetBinError(5,0);
  //hJ0->SetBinError(6,0);
  //hJ0->SetBinContent(5,0);
  //hJ0->SetBinContent(6,0);

  hReal3->Divide(hReal3,hEff);

  TH1D * hRealJ0 = (TH1D *)hReal3->Clone("hRealJ0");
  TH1D * hRealJ1 = (TH1D *) hReal3->Clone("hRealJ1");

  TH1D * hJ1A = (TH1D *) hReal3->Clone("hJ1A");
  TH1D * hJ1B = (TH1D *) hReal3->Clone("hJ1B");
  TH1D * hJ1D = (TH1D *) hReal3->Clone("hJ1C");
  TH1D * hJ1C = (TH1D *) hReal3->Clone("hJ1D");
  
  TH1D * hJ0A = (TH1D *) hReal3->Clone("hJ0A");
  TH1D * hJ0B = (TH1D *) hReal3->Clone("hJ0B");
  TH1D * hJ0C = (TH1D *) hReal3->Clone("hJ0C");

  hJ0A->SetLineStyle(9);
  hJ0B->SetLineStyle(9);
  hJ0C->SetLineStyle(9);

  hJ0A->SetLineColor(2);
  hJ0B->SetLineColor(2);
  hJ0C->SetLineColor(2);

  hJ0A->SetLineWidth(3);
  hJ0B->SetLineWidth(3);
  hJ0C->SetLineWidth(3);

  hJ1A->SetLineStyle(9);
  hJ1B->SetLineStyle(9);
  hJ1C->SetLineStyle(9);
  hJ1D->SetLineStyle(9);

  hJ1A->SetLineColor(3);
  hJ1B->SetLineColor(3);
  hJ1C->SetLineColor(3);
  hJ1D->SetLineColor(3);

  hJ1A->SetLineWidth(3);
  hJ1B->SetLineWidth(3);
  hJ1C->SetLineWidth(3);
  hJ1D->SetLineWidth(3);

  hRealJ0->Multiply(hRealJ0,hJ0);
  hRealJ0->SetLineWidth(2);
  hRealJ0->SetLineColor(2);
  hRealJ0->SetFillStyle(1001);
  hRealJ0->SetFillColor(2);

  hRealJ1->Multiply(hRealJ1,hJ1);
  hRealJ1->SetLineWidth(2);
  hRealJ1->SetLineColor(3);
  hRealJ1->SetFillStyle(1001);
  hRealJ1->SetFillColor(3);

  MyFitFun1 * fptr1 = new MyFitFun1();
  int nPar = 8;
  TF1 *f1 = new TF1("f1",fptr1,&MyFitFun1::Evaluate,0.0,4.0,nPar,
		    "MyFitFun1","Evaluate");

  MyFitFun2 * fptr2 = new MyFitFun2();
  nPar = 11;

  TF1 *f2 = new TF1("f2",fptr2,&MyFitFun2::Evaluate,0.0,4.0,nPar,
		    "MyFitFun2","Evaluate");

  //TF1 *fg1 = new TF1("fg1","gaus+pol1",0.0,3.0);
  //TF1 *fg2a0 = new TF1("fg2a0","gaus+pol1",0.0,3.0);

  f1->SetParameter(0,20);
  f1->SetParameter(1,1.294);
  f1->SetParameter(2,0.055);
  f1->SetParameter(3,30);
  f1->SetParameter(4,1.4139);
  f1->SetParameter(5,0.0501);


  f2->SetParameter(0,30);
  //f2->SetParameter(1,1.2819);
  f2->SetParameter(1,1.276);
  f2->SetParameter(2,0.028);
  f2->SetParameter(3,12);
  f2->SetParameter(4,1.4263);
  f2->SetParameter(5,0.0545);
  f2->SetParameter(6,300);
  f2->SetParameter(7,1.36);
  f2->SetParameter(8,0.09);

  f2->SetParLimits(1,1.2819-0.02,1.2819+0.02);
  f2->SetParLimits(2,0.0227-0.011,0.0227+0.011);
  f2->SetParLimits(4,1.4263-0.09,1.4263+0.09);
  f2->SetParLimits(5,0.0545-0.026,0.0545+0.026);

  f2->SetParLimits(7,1.35,1.38);
  f2->SetParLimits(8,0.02,0.12);
  //ParLimit


  if (pdgLim == 1){
    f1->SetParLimits(1,1.294-0.004,1.294+0.004);
    f1->SetParLimits(2,0.055-0.005,0.055+0.005);
    f1->SetParLimits(4,1.4139-0.0017,1.4139+0.0017);
    f1->SetParLimits(5,0.048-004,0.048+004);
    
    f2->SetParLimits(1,1.2819-0.0005,1.2819+0.0005);
    f2->SetParLimits(2,0.0227-0.0011,0.0227+0.0011);
    f2->SetParLimits(4,1.4263-0.0009,1.4263+0.0009);
    f2->SetParLimits(5,0.0545-0.0026,0.0545+0.0026);
  }
  if (pdgLim == 2){
    f2->SetParLimits(1,1.2819-0.0005,1.2819+0.0005);
    f2->SetParLimits(2,0.0227-0.0011,0.0227+0.0011);
    f2->SetParLimits(4,1.4263-0.0009,1.4263+0.0009);
    f2->SetParLimits(5,0.0545-0.0026,0.0545+0.0026);

  }
  //fg1->SetParameter(3,0);
  //fg1->SetParameter(4,0);

  f1->SetLineColor(2);
  f2->SetLineColor(3);

  double fLow = 1.26;
  double fHigh = 1.45;

  hRealJ0->Fit("f1","RB","",fLow,fHigh);
  
  fLow = 1.25;
  fHigh = 1.45;
  //  return;
  hRealJ1->Fit("f2","RB","",fLow,fHigh);

  double ampJ0A = f1->GetParameter(0);f1->SetParameter(0,0);
  double ampJ0B = f1->GetParameter(3);f1->SetParameter(3,0);
  hJ0C->Eval(f1);

  double ampJ0C1 = f1->GetParameter(6);f1->SetParameter(6,0);
  double ampJ0C2 = f1->GetParameter(7);f1->SetParameter(7,0);
  f1->SetParameter(0,ampJ0A);
  hJ0A->Eval(f1);

  f1->SetParameter(0,0);
  f1->SetParameter(3,ampJ0B);
  hJ0B->Eval(f1);
  
  double ampJ1A = f2->GetParameter(0);f2->SetParameter(0,0);
  double ampJ1B = f2->GetParameter(3);f2->SetParameter(3,0);
  double ampJ1C = f2->GetParameter(6);f2->SetParameter(6,0);
  hJ1D->Eval(f2);

  double ampJ1D1 = f2->GetParameter(9);f2->SetParameter(9,0);
  double ampJ1D2 = f2->GetParameter(10);f2->SetParameter(10,0);
  f2->SetParameter(0,ampJ1A);
  hJ1A->Eval(f2);

  f2->SetParameter(0,0);
  f2->SetParameter(3,ampJ1B);
  hJ1B->Eval(f2);

  f2->SetParameter(3,0);
  f2->SetParameter(6,ampJ1C);
  hJ1C->Eval(f2);

  hReal3->SetMinimum(0);
  hReal3->GetXaxis()->SetRangeUser(1.25,1.45);
  //hJ0A->GetXaxis()->SetRangeUser(1.28,1.45);
  //hJ0B->GetXaxis()->SetRangeUser(1.28,1.45);
  //hJ0C->GetXaxis()->SetRangeUser(1.28,1.45);

  cout<<"----------"<<endl;
  double cJ0 = 1000*f1->GetParameter(1);
  double cJ0Err = 1000*f1->GetParError(1);
  cout<<"cJ0 = "<<cJ0<<" +/- "<<cJ0Err<<endl;

  double wJ0 = 1000*f1->GetParameter(2);
  double wJ0Err = 1000*f1->GetParError(2);
  cout<<"wJ0 = "<<wJ0<<" +/- "<<wJ0Err<<endl;

  //cout<<"----------"<<endl;
  double cJ1 = 1000*f2->GetParameter(1);
  double cJ1Err = 1000*f2->GetParError(1);
  cout<<"cJ1 = "<<cJ1<<" +/- "<<cJ1Err<<endl;

  double wJ1 = 1000*f2->GetParameter(2);
  double wJ1Err = 1000*f2->GetParError(2);
  cout<<"wJ1 = "<<wJ1<<" +/- "<<wJ1Err<<endl;

  hReal3->SetMarkerStyle(21);
  hReal3->GetXaxis()->SetTitle("Mass(K^{+}K^{-}#pi^{0})/GeV");
  hReal3->GetYaxis()->SetTitle("Counts");

  hReal1->SetMarkerStyle(20);
  hReal1->GetXaxis()->SetTitle("Mass(K^{+}K^{-}#pi^{0})/GeV");
  hReal1->GetYaxis()->SetTitle("Counts");

  hRealJ0->GetXaxis()->SetTitle("Mass(K^{+}K^{-}#pi^{0})/GeV");
  hRealJ0->GetYaxis()->SetTitle("Counts");

  //hReal3->GetXaxis()->SetRangeUser(1.24,1.31);
  hReal3->SetLineWidth(2);

  if (iShow == 0) {
    hReal3->Draw("e1");
  }


  if (iShow == 1) {
    hReal3->Draw("e1");
    hRealJ0->Draw("same");
    hRealJ1->Draw("same");

    hJ0A->Draw("chistsame");
    hJ0B->Draw("chistsame");
    hJ0C->Draw("chistsame");

    hJ1A->Draw("chistsame");
    hJ1B->Draw("chistsame");
    hJ1C->Draw("chistsame");
    hJ1D->Draw("chistsame");

    Double_t xl1=0.63, yl1 = 0.73, xl2=xl1+0.23, yl2=yl1+0.12;
    TLegend *legend = new TLegend(xl1,yl1,xl2,yl2);
    legend->SetFillColor(0);
    legend->SetLineColor(0);
    legend->SetTextSize(0.04);
    legend->AddEntry(hRealJ0,"J=0^{-}","f");
    legend->AddEntry(hRealJ1,"J=1^{+}","f");

    legend->Draw();

  }
  if (iShow == 11) {
    hReal3->Draw("e1");
    hRealJ0->Draw("same");
    //hRealJ1->Draw("same");

    Double_t xl1=0.13, yl1 = 0.73+0.06, xl2=xl1+0.23, yl2=yl1+0.06;
    TLegend *legend = new TLegend(xl1,yl1,xl2,yl2);
    legend->SetFillColor(0);
    legend->SetLineColor(0);
    legend->SetTextSize(0.04);
    legend->AddEntry(hRealJ0,"J=0^{-}","f");
    //legend->AddEntry(hRealJ1,"J=1","f");

    legend->Draw();

  }

  Double_t xl1=0.12, yl1 = 0.8, xl2=xl1+0.23, yl2=yl1+0.06;
  TLegend *legend = new TLegend(xl1,yl1,xl2,yl2);
  legend->SetFillColor(0);
  legend->SetLineColor(0);
  legend->SetTextSize(0.04);
 
  if (iShow == 2) {
    hRealJ0->Draw("e");
    legend->AddEntry(hRealJ0,"J=0","f");
    legend->Draw();
  }
  if (iShow == 3) {
    //hRealJ1L1S0->Draw("e");
    //legend->AddEntry(hRealJ1L1S0,"J=1,L=1,S=0","f");
    //legend->Draw();
  }
  if (iShow == 4) {
    //hRealJ1->Draw("e");
    //legend->AddEntry(hRealJ1,"J=1,L=1,S=0 a0","f");
    //legend->Draw();
  }
  if (iShow == 5) {
    //hRealJ1L0S1->Draw("e");
    //legend->AddEntry(hRealJ1L0S1,"J=1,L=0,S=1","f");
    //legend->Draw();
  }
  if (iShow == 6) {
    //hRealJ1L1S1->Draw("e");
    //legend->AddEntry(hRealJ1L1S1,"J=1,L=1,S=1","f");
    //legend->Draw();
  }
}