const Float_t kELMASS = 0.00051099907;
const Float_t kPMASS  = 0.93827231;
const Float_t kKMASS  = 0.49367;
const Float_t kPIMASS = 0.13957;

const Char_t *label_str[3]={"40 #font[12]{A}#upointGeV",
			    "80 #font[12]{A}#upointGeV",
			    "158 #font[12]{A}#upointGeV"};

void grey_palette()
{
  col=(TColor*) gROOT->GetListOfColors()->At(2);
  col->SetRGB(0.8,0.8,0.8);
  col=(TColor*) gROOT->GetListOfColors()->At(3);
  col->SetRGB(0.6,0.6,0.6);
  Int_t palette[20];
  for (Int_t i=0;i<20; i++) {
    col=(TColor*) gROOT->GetListOfColors()->At(50+i);
    col->SetRGB(1-0.05*i,1-0.05*i,1-0.05*i);
    palette[i]=50+i;
  }
  gStyle->SetPalette(15,palette);
}

void draw_numbers(TCanvas *c1,Int_t n_pad) {
  TLatex ltx;
  ltx.SetNDC();
  ltx.SetTextSize(0.07);
  for (Int_t i=0; i< n_pad; i++) {
    c1->cd(0);
    TVirtualPad::Pad()->cd(i+1);
    Char_t num[10];
    sprintf(num,"%c)",97+i);
    cout << TVirtualPad::Pad() << " " << num << endl;
    ltx.DrawLatex(0.22,0.8,num);
  }
}

Float_t get_txt_size(TVirtualPad *pad, const Float_t pix_size) {
  if (pad->UtoPixel(1) > pad->VtoPixel(0)) {
    // 'Horizontal pad'
    return pix_size/pad->VtoPixel(0);
  }
  else {
    // 'Vertical pad'
    return pix_size/pad->UtoPixel(1);
  }
}

void set_text_sizes(TVirtualPad *pad, TH1F* hist) {
  const Float_t label_size=9;
  const Float_t title_size=11;
  hist->SetLabelSize(get_txt_size(pad,label_size));
  hist->SetLabelSize(get_txt_size(pad,label_size),"Y");
  hist->SetTitleSize(get_txt_size(pad,title_size));
  hist->SetTitleSize(get_txt_size(pad,title_size),"Y");
}


TCanvas *make_canvas(const Char_t *name, const Char_t *title, Int_t n_x=1, Int_t n_y=1, Int_t share_axes=0, Int_t ww=0, Int_t wh=0){
  const Int_t width=350;
  const Int_t height=350;
  TCanvas *canvas;
  if (share_axes==0) {
    if (ww==0)
      ww=width*n_x;
    if (wh==0)
      wh=height*n_y;
    canvas=new TCanvas(name,title,ww,wh);
    canvas->Divide(n_x,n_y,0,0);
  }
  else {
    Float_t pix_width=(1-gStyle->GetPadLeftMargin()-gStyle->GetPadRightMargin())*width;
    Float_t pix_height=(1-gStyle->GetPadTopMargin()-gStyle->GetPadBottomMargin())*height;
    if (ww==0)
      ww=width+(n_x-1)*pix_width;
    if (wh==0)
      wh=height+(n_y-1)*pix_height;

    canvas=new TCanvas(name,title,ww,wh);
    Float_t tot_width;
    if (n_x>1) 
      tot_width=(n_x-2)+1./(1-gStyle->GetPadLeftMargin())
	+1./(1-gStyle->GetPadRightMargin());
    else
      tot_width=1./(1-gStyle->GetPadLeftMargin()-gStyle->GetPadRightMargin());
    Float_t tot_height;
    if (n_y>1) 
      tot_height=(n_y-2)+1./(1-gStyle->GetPadTopMargin())
	+1./(1-gStyle->GetPadBottomMargin());
    else
      tot_height=1./(1-gStyle->GetPadTopMargin()-gStyle->GetPadBottomMargin());

    //Int_t idx=n_x*n_y;
    
    for (Int_t j=0; j<n_y; j++) {
      for (Int_t i=0; i<n_x; i++) {
	//for (Int_t j=n_y-1; j>=0; j--) {
	//for (Int_t i=n_x-1; i>=0; i--) {
	Char_t tmp_str[256];
	Char_t p_title[256];
	Int_t idx=n_x*j+i+1;
	sprintf(tmp_str,"%s_%d",canvas->GetName(),idx);
	sprintf(p_title,"Pad %d",idx);
	Float_t x1=0,y1=0;
	Float_t x2=1,y2=1;
	if (n_x>1) {
	  if (i==0) 
	    x2=1./(1-gStyle->GetPadLeftMargin())/tot_width;
	  else {
	    x1=(1./(1-gStyle->GetPadLeftMargin())+i-1)/tot_width;
	    if (i<n_x-1)
	      x2=(1./(1-gStyle->GetPadLeftMargin())+i)/tot_width;
	  }
	}
	if (n_y>1) {
	  if (j==0) 
	    y1=1-1./(1-gStyle->GetPadTopMargin())/tot_height;
	  else {
	    y2=1-(1./(1-gStyle->GetPadTopMargin())+j-1)/tot_height;
	    if (j<n_y-1)
	      y1=1-(1./(1-gStyle->GetPadTopMargin())+j)/tot_height;
	  }
	}
	//cout << "x1 " << x1 << ", x2 " << x2 << endl;
	TPad *pad=new TPad(tmp_str,title,x1,y1,x2,y2);
	//pad->SetFillColor(idx+1);
	if (i>0)
	  pad->SetLeftMargin(0.004);
	if (i<n_x-1)
	  pad->SetRightMargin(0.004);
	if (j>0)
	  pad->SetTopMargin(0.004);
	if (j<n_y-1)
	  pad->SetBottomMargin(0.004);

	pad->SetNumber(idx);
	//pad->SetNumber(n_x*j+i+1);
	pad->Draw();
	//idx--;
	//idx++;
	//pad->SetP
      }
    }
  }

  return canvas;
}

TPolyLine3D *draw_func_3d(TF1 *f1) {
  Double_t xmin,xmax;
  Bool_t logx=TVirtualPad::Pad()->GetLogx();
  Bool_t logy=TVirtualPad::Pad()->GetLogy();
  Float_t vmin[3],vmax[3];
  TView *view=TVirtualPad::Pad()->GetView();
  if (view==0){
    cout << "ERROR in draw_func_3d : no view found" << endl; 
    return 0;
  }
  view->GetRange(vmin,vmax);
  Int_t npmax=f1->GetNpx();
  Float_t *xx=new Float_t[npmax];
  Float_t *yy=new Float_t[npmax];
  Float_t *zz=new Float_t[npmax];
  f1->GetRange(xmin,xmax);
  if (logx) {
    if (xmin<=0)
      xmin=vmin[0];
    else
      xmin=log10(xmin);
    xmax=log10(xmax);
  }  
  Float_t step=(xmax-xmin)/(npmax-1);
  Int_t np=0;
  for (Int_t i=0; i < npmax; i++) {
    Float_t x=xmin+step*i;
    Float_t y;
    if (logx)
      y=f1->Eval(pow(10,x));
    else
      y=f1->Eval(x);
    if (logy) {
      if (y<=0)
	y=vmin[1];
      else
	y=log10(y);
    }
    if (y>=vmin[1] && y<=vmax[1] &&
	x>=vmin[0] && x<=vmax[0]) {
      xx[np]=x;
      yy[np]=y;
      zz[np]=vmin[2];
      np++;
    }
  }
  TPolyLine3D *poly=new TPolyLine3D(np,xx,yy,zz);
  poly->SetLineStyle(f1->GetLineStyle());
  poly->SetLineColor(f1->GetLineColor());
  poly->Draw();
  delete [] xx;
  delete [] yy;
  delete [] zz;
  return poly;
}

TH2F *rebin_2d(TH2F *hist,Int_t nx=2,Int_t ny=2,const Char_t *name) { 
  
  Float_t xmin=hist->GetXaxis()->GetXmin();
  Float_t xmax=hist->GetXaxis()->GetXmax();
  Float_t ymin=hist->GetYaxis()->GetXmin();
  Float_t ymax=hist->GetYaxis()->GetXmax();
  
  Int_t nb_x=hist->GetXaxis()->GetNbins();
  nb_x=nb_x%nx?nb_x/nx+1:nb_x/nx;
  xmax=xmin+nb_x*nx*hist->GetXaxis()->GetBinWidth(1);
  Int_t nb_y=hist->GetYaxis()->GetNbins();
  nb_y=nb_x%ny?nb_y/ny+1:nb_y/ny;
  ymax=ymin+nb_y*ny*hist->GetYaxis()->GetBinWidth(1);
  
  TH2F *n_hist=new TH2F(name,hist->GetTitle(),nb_x,xmin,xmax,nb_y,ymin,ymax);

  for (Int_t bx=1; bx<=nb_x; bx++)
    for (Int_t by=1; by<=nb_y; by++) 
      n_hist->Fill(hist->GetXaxis()->GetBinCenter(bx),hist->GetYaxis()->GetBinCenter(by),hist->GetCellContent(bx,by));
   
  return n_hist;
}

TPad *overlay_2d(const Char_t *name="ovl_p",Char_t *title="Overlay pad"){
  TPad *ovl_p=new TPad(name,title,0,0,1,1);
  ovl_p->SetFillStyle(4000);
  ovl_p->RangeAxis(0,0,2,2);
  Float_t dxr=2./(1-TVirtualPad::Pad()->GetLeftMargin()-TVirtualPad::Pad()->GetRightMargin());
  Float_t dyr=2./(1-TVirtualPad::Pad()->GetTopMargin()-TVirtualPad::Pad()->GetBottomMargin());
  ovl_p->Range(-TVirtualPad::Pad()->GetLeftMargin()*dxr,-TVirtualPad::Pad()->GetBottomMargin()*dyr,
	       2+TVirtualPad::Pad()->GetRightMargin()*dxr,2+TVirtualPad::Pad()->GetTopMargin()*dyr);
  ovl_p->RangeAxis(0,0,2,2);
  ovl_p->Draw();
  ovl_p->cd();
  return ovl_p;
}

TGraph *draw_graph_error(TGraphErrors *gr) {
  Int_t np=gr->GetN();
  TGraph *tmp_gr=new TGraph(2*np);
  for (Int_t i=0; i<np; i++) {
    Double_t x,y;
    gr->GetPoint(i,x,y);
    tmp_gr->SetPoint(i,x,y+gr->GetErrorY(i));
    tmp_gr->SetPoint(2*np-i-1,x,y-gr->GetErrorY(i));
  }
  tmp_gr->Draw("f");
  tmp_gr->Draw("l");
  gr->Draw("pX");
  return tmp_gr;
}

void energy_label(Int_t i, Float_t x=0,Float_t y=0, TLatex *ltx=0) {
  if (ltx==0) {
    ltx=new TLatex();
    ltx->SetTextSize(18);
    ltx->SetTextAlign(21);
    ltx->SetNDC();
  }
  if (x==0) 
    x=gPad->GetLeftMargin()+0.5*(1-gPad->GetRightMargin()-gPad->GetLeftMargin());
  if (y==0) 
    y=1-gPad->GetTopMargin()-0.10;

  ltx->DrawLatex(x,y,label_str[i]);
}

Float_t ndc_x(Float_t x) {
  return (1-TVirtualPad::Pad()->GetLeftMargin()-TVirtualPad::Pad()->GetRightMargin())*x+TVirtualPad::Pad()->GetLeftMargin();
}

Float_t ndc_y(Float_t y) {
  return (1-TVirtualPad::Pad()->GetTopMargin()-TVirtualPad::Pad()->GetBottomMargin())*y+TVirtualPad::Pad()->GetBottomMargin();
}

void draw_fill_graph(TGraph *gr, Char_t *opt="p") {
  TGraph *gr_fill=(TGraph*) gr->Clone();
  gr_fill->SetMarkerStyle(gr->GetMarkerStyle()-4);
  gr_fill->SetMarkerColor(10);
  Char_t optx[25];
  sprintf(optx,"%sx",opt);
  gr_fill->Draw(optx);
  gr->Draw(opt);
}

void draw_legend_m(Float_t x, Float_t y, TAttMarker *att_m, const Char_t *label,Int_t i_ent, Int_t colored=1) {
  y-=0.07*i_ent;
  Float_t x_real=TVirtualPad::Pad()->GetUxmin()+(x+0.035)*(TVirtualPad::Pad()->GetUxmax()-TVirtualPad::Pad()->GetUxmin());
  Float_t y_real=TVirtualPad::Pad()->GetUymin()+y*(TVirtualPad::Pad()->GetUymax()-TVirtualPad::Pad()->GetUymin());
  Int_t m_style=att_m->GetMarkerStyle();
  TMarker *mrk=new TMarker(x_real,y_real,m_style);
  Float_t m_siz=att_m->GetMarkerSize();
  mrk->SetMarkerSize(m_siz);
  if (colored)
    mrk->SetMarkerColor(att_m->GetMarkerColor());
  mrk->Draw();
  TLatex *ltx=new TLatex();
  ltx->SetTextFont(133);
  ltx->SetTextSize(16);
  ltx->SetTextAlign(12);
  if (colored)
    ltx->SetTextColor(att_m->GetMarkerColor());
  //ltx->SetNDC();
  ltx->DrawLatex(x+0.08,y,label);
}

void draw_legend_l(Float_t x, Float_t y, TAttLine *att_l, const Char_t *label,Int_t i_ent, Size_t txt_size=16) {
  //y-=0.07*i_ent;
  Float_t pix_y_siz=TVirtualPad::Pad()->AbsPixeltoY(1)-TVirtualPad::Pad()->AbsPixeltoY(2);
  y-=pix_y_siz*txt_size*1.05*i_ent/(TVirtualPad::Pad()->GetUymax()-TVirtualPad::Pad()->GetUymin());

  Float_t x_min_real=TVirtualPad::Pad()->GetUxmin()+(x+0.01)*(TVirtualPad::Pad()->GetUxmax()-TVirtualPad::Pad()->GetUxmin());
  Float_t x_max_real=TVirtualPad::Pad()->GetUxmin()+(x+0.08)*(TVirtualPad::Pad()->GetUxmax()-TVirtualPad::Pad()->GetUxmin());
  Float_t y_real=TVirtualPad::Pad()->GetUymin()+y*(TVirtualPad::Pad()->GetUymax()-TVirtualPad::Pad()->GetUymin());
  Int_t l_style=att_l->GetLineStyle();
  TLine *line=new TLine(x_min_real,y_real,x_max_real,y_real);
  Float_t l_siz=att_l->GetLineWidth();
  line->SetLineWidth(l_siz);
  line->SetLineStyle(l_style);
  line->SetLineColor(att_l->GetLineColor());
  line->Draw();
  TLatex *ltx=new TLatex();
  ltx->SetTextFont(133);
  ltx->SetTextSize(txt_size);
  ltx->SetTextAlign(12);
  ltx->SetTextColor(att_l->GetLineColor());
  ltx->SetNDC();
  ltx->DrawLatex(ndc_x(x+0.1),ndc_y(y),label);
}

void scale_graph(TGraphErrors *gr,Float_t scale) {
  for (Int_t i=0; i<gr->GetN(); i++) {
    Double_t x,y;
    gr->GetPoint(i,x,y);
    gr->SetPoint(i,x,scale*y);
    gr->SetPointError(i,gr->GetErrorX(i),scale*gr->GetErrorY(i));
  }
}

void scale_graph(TGraph *gr,Float_t scale) {
  for (Int_t i=0; i<gr->GetN(); i++) {
    Double_t x,y;
    gr->GetPoint(i,x,y);
    gr->SetPoint(i,x,scale*y);
  }
}

void set_ratio(TGraphErrors *gr, Int_t pnt, Float_t x_val, Float_t num, Float_t den, Float_t num_err=0, Float_t den_err=0) {
  Float_t rat=num/den;
  gr->SetPoint(i,x_val,rat);
  if (rat!=0) {
    Float_t rat_err=rat*sqrt(num_err*num_err/num/num+den_err*den_err/den/den);
    gr->SetPointError(i,0,rat_err);
  }
}

void set_chi_levels(TH2 *h1, const Int_t n_cont, Double_t *c_levels, Double_t *x_min=0, Double_t *y_min=0) {

  Float_t min_chi=1e10;
  Int_t min_xbin=0, min_ybin=0;
  
  for (Int_t i=1; i <= h1->GetNbinsX(); i++) {
    for (Int_t j=1; j <= h1->GetNbinsY(); j++) {
      if (h1->GetCellContent(i,j) < min_chi) {
	min_xbin=i;
	min_ybin=j;
	min_chi=h1->GetCellContent(i,j);
      }
    }
  }
  
  cout << "min_xbin " << min_xbin << ", " << min_ybin << ", min_chi " << min_chi << endl;

  Double_t *c_levels_tmp=new Double_t[n_cont];
  for (Int_t i=0; i<n_cont;i++) {
    c_levels_tmp[i]=c_levels[i]+min_chi;
    cout << "c_levels[" <<i << "] " << c_levels[i] << ", tmp " << c_levels_tmp[i] << endl;
  }
  h1->SetContour(n_cont,c_levels_tmp);
  //delete c_levels_tmp;
  
  if (x_min) 
    *x_min=h1->GetXaxis()->GetBinCenter(min_xbin);
  if (y_min) 
    *y_min=h1->GetYaxis()->GetBinCenter(min_ybin);

}


void make_gr_file(const Char_t *name,TGraphErrors *gr,const Char_t mode='w') {
  //ios::open_mode m_flag=ios::out|ios::trunc;
  int m_flag=ios_base::out|ios_base::trunc;
  if (mode=='a')
    m_flag=ios_base::out|ios_base::app;
  //  m_flag=ios::out|ios::app;
  ofstream fout(name,m_flag);
  for (Int_t i=0; i< gr->GetN(); i++) {
    Double_t x,y;
    gr->GetPoint(i,x,y);
    fout << x << " " << y;
    fout << " " << gr->GetErrorY(i)<< endl;
  }
  fout.close();
}