#include <TFile.h>
#include <TStyle.h>
#include <TSystem.h>
#include <TROOT.h>
#include <TCanvas.h>
#include <TGraphErrors.h>
#include <TMultiGraph.h>
#include <TH1.h>
#include <TH3.h>
#include <TString.h>
#include <TMath.h>
#include <TLatex.h>
#include <TPaveText.h>
#include <TLegend.h>
#include <TImage.h>

void SetGlobalStyle();
void SetStarStyle(TGraphErrors* grErr, int energy, int cent);
void SetTGraphStyle(TGraphErrors* grErr, int pid, int cent);
//void SetOldTGraphStyle(TGraphErrors* grErr, int pid, int cent);

//_________________
void PresentResults() {

  SetGlobalStyle();
  double mLamRange[2] = {0.1, 0.9};
  double mRadiiRange[2] = {0.6, 7.1};
  double mMtRange[2] = {0.15, 1.35};
  int mCanvSize[2] = {600,800};

  int nEnergies = 1;

  const int nPionKtBins = 7;
  int PionFitKtArray[7][2] = {{2,2}, {3,3}, {4,4}, {5,5}, {6,6}, {7,7}, {8,10}};
  const int nKaonKtBins = 5;
  int KaonFitKtArray[5][2] = {{0,3}, {4,5}, {6,7}, {7,8}, {9,10}};

  int nStarMtBins = 4;
  
  int nStar200CentBins = 5;
  //200 GeV
  double nStarRout200[5][4] = { {5.88623, 5.36104, 4.68839, 4.13888},
				{5.62342, 5.12386, 4.50129, 3.99388},
				{5.19246, 4.7173, 4.21355, 3.7131},
				{4.70398, 4.29394, 3.84329, 3.43377},
				{4.21394, 3.88334, 3.50813, 3.1131} };
  double nStarRoutErr200[5][4] = { {0.0235863, 0.0156105, 0.0166256, 0.019129},
				   {0.0225075, 0.0146684, 0.0168094, 0.0188995},
				   {0.0162752, 0.0108589, 0.0126066, 0.014325},
				   {0.0184947, 0.012548, 0.0142582, 0.0158074},
				   {0.0204706, 0.0137939, 0.0165987, 0.018784} };
  double nStarRside200[5][4] = { {5.20713, 4.9248, 4.532, 4.15774},
				 {4.93469, 4.66952, 4.34595, 3.97151},
				 {4.60779, 4.30443, 4.03165, 3.72077},
				 {4.12645, 3.90416, 3.68447, 3.41898},
				 {3.71263, 3.51214, 3.30994, 3.0682} };
  double nStarRsideErr200[5][4] = { {0.018666, 0.013728, 0.0152585, 0.0155588},
				    {0.0185243, 0.0124062, 0.0145618, 0.0165696},
				    {0.0125847, 0.00925802, 0.0110111, 0.0130999},
				    {0.0142614, 0.0111244, 0.0129176, 0.0143846},
				    {0.0154815, 0.011951, 0.0150637, 0.0168027} };
  double nStarRlong200[5][4] = { {6.64397, 5.80251, 4.96302, 4.34703},
				 {6.30332, 5.48477, 4.72166, 4.06245},
				 {5.84876, 5.04122, 4.36922, 3.74956},
				 {5.33634, 4.57213, 3.95305, 3.38169},
				 {4.82185, 4.1227, 3.58624, 3.05169} };
  double nStarRlongErr200[5][4] = { {0.0267233, 0.0160585, 0.0181593, 0.0178103},
				    {0.0261458, 0.015864, 0.0172435, 0.0183965},
				    {0.0177964, 0.0117965, 0.0129835, 0.0145198},
				    {0.0195327,0.0140871, 0.0150233, 0.0157312},
				    {0.0224088, 0.0151745, 0.0177424, 0.018657} };
  double nStarMt200[5][4] = { {0.262197, 0.328421, 0.417164, 0.524964},
			      {0.261902, 0.328167, 0.417082, 0.525045},
			      {0.261589, 0.327888, 0.41699, 0.525108},
			      {0.261225, 0.327544,0.416849,0.525124},
			      {0.260913, 0.327235, 0.416687, 0.52502} };
  double nStarMtErr200[5][4] = { {0.01, 0.01, 0.01, 0.01},
				 {0.01, 0.01, 0.01, 0.01},
				 {0.01, 0.01, 0.01, 0.01},
				 {0.01, 0.01, 0.01, 0.01},
				 {0.01, 0.01, 0.01, 0.01} };

  //19 GeV
  int nStar19CentBins = 2;
  double nStarRout19[2][4] = { {5.84087, 5.37106, 4.79477, 4.31619},
			       {4.56516, 4.26748, 3.87451, 3.45468} };
  double nStarRoutErr19[2][4] = { {0.0195966,0.0148756,0.0186787,0.0226372},
				  {0.0188045,0.0136939,0.0175341,0.0207569} };
  double nStarRside19[2][4] = { {4.83625,4.5331,4.23409,3.87759},
				{3.81836,3.59575,3.36289,3.13035} };
  double nStarRsideErr19[2][4] = { {0.0150974,0.0125819,0.0151165,0.018952},
				   {0.0135812,0.0116497,0.0145794,0.0176817} };
  double nStarRlong19[2][4] = { {5.7974,4.97864,4.23506,3.57462},
				{4.6873,4.02044,3.38523,2.85863} };
  double nStarRlongErr19[2][4] = { {0.0196311,0.0153343,0.0167063,0.0197918},
				   {0.0196429,0.0146235,0.0165197,0.018152} };
  double nStarMt19[2][4] = { {0.259659,0.326347,0.416095,0.523699},
			     {0.259161,0.325885,0.415924,0.523838} };
  double nStarMtErr19[2][4] = { {0.01, 0.01, 0.01, 0.01},
				{0.01, 0.01, 0.01, 0.01} };
 

  TGraphErrors *StarRoutMt200GrErr[5];
  TGraphErrors *StarRsideMt200GrErr[5];
  TGraphErrors *StarRlongMt200GrErr[5];



  //200 GeV
  for(int iCent=0; iCent<nStar200CentBins; iCent++) {

    float tmpRout[4],tmpRside[4],tmpRlong[4];
    float tmpRoutErr[4],tmpRsideErr[4],tmpRlongErr[4];
    float tmpMt[4], tmpMtErr[4];
    for(int iMt=0; iMt<4; iMt++) {
      tmpRout[iMt] = nStarRout200[iCent][iMt];
      tmpRoutErr[iMt] = nStarRoutErr200[iCent][iMt];
      tmpRside[iMt] = nStarRside200[iCent][iMt];
      tmpRsideErr[iMt] = nStarRsideErr200[iCent][iMt];
      tmpRlong[iMt] = nStarRlong200[iCent][iMt];
      tmpRlongErr[iMt] = nStarRlongErr200[iCent][iMt];
      tmpMt[iMt] = nStarMt200[iCent][iMt];
      tmpMtErr[iMt] = nStarMtErr200[iCent][iMt];
    }
    //Out
    StarRoutMt200GrErr[iCent] = new TGraphErrors(nStarMtBins,tmpMt,tmpRout,
						 tmpMtErr,tmpRoutErr);
    StarRoutMt200GrErr[iCent]->SetName(Form("StarRoutMt200GrErr_%d",iCent));
    SetStarStyle(StarRoutMt200GrErr[iCent],200,iCent);
    //Side
    StarRsideMt200GrErr[iCent] = new TGraphErrors(nStarMtBins,tmpMt,tmpRside,
						 tmpMtErr,tmpRsideErr);
    StarRsideMt200GrErr[iCent]->SetName(Form("StarRsideMt200GrErr_%d",iCent));
    SetStarStyle(StarRsideMt200GrErr[iCent],200,iCent);
    //Long
    StarRlongMt200GrErr[iCent] = new TGraphErrors(nStarMtBins,tmpMt,tmpRlong,
						 tmpMtErr,tmpRlongErr);
    StarRlongMt200GrErr[iCent]->SetName(Form("StarRlongMt200GrErr_%d",iCent));
    SetStarStyle(StarRlongMt200GrErr[iCent],200,iCent);

  }

  TGraphErrors *StarRoutMt19GrErr[2];
  TGraphErrors *StarRsideMt19GrErr[2];
  TGraphErrors *StarRlongMt19GrErr[2];
  //19 GeV
  for(int iCent=0; iCent<nStar19CentBins; iCent++) {
    
    float tmpRout[4],tmpRside[4],tmpRlong[4];
    float tmpRoutErr[4],tmpRsideErr[4],tmpRlongErr[4];
    float tmpMt[4], tmpMtErr[4];

    for(int iMt=0; iMt<4; iMt++) {
      tmpRout[iMt] = nStarRout19[iCent][iMt];
      tmpRoutErr[iMt] = nStarRoutErr19[iCent][iMt];
      tmpRside[iMt] = nStarRside19[iCent][iMt];
      tmpRsideErr[iMt] = nStarRsideErr19[iCent][iMt];
      tmpRlong[iMt] = nStarRlong19[iCent][iMt];
      tmpRlongErr[iMt] = nStarRlongErr19[iCent][iMt];
      tmpMt[iMt] = nStarMt19[iCent][iMt];
      tmpMtErr[iMt] = nStarMtErr19[iCent][iMt];
    }
    //Out
    StarRoutMt19GrErr[iCent] = new TGraphErrors(nStarMtBins,tmpMt,tmpRout,
						 tmpMtErr,tmpRoutErr);
    StarRoutMt19GrErr[iCent]->SetName(Form("StarRoutMt19GrErr_%d",iCent));
    SetStarStyle(StarRoutMt19GrErr[iCent],19,iCent);
    //Side
    StarRsideMt19GrErr[iCent] = new TGraphErrors(nStarMtBins,tmpMt,tmpRside,
						 tmpMtErr,tmpRsideErr);
    StarRsideMt19GrErr[iCent]->SetName(Form("StarRsideMt19GrErr_%d",iCent));
    SetStarStyle(StarRsideMt19GrErr[iCent],19,iCent);
    //Long
    StarRlongMt19GrErr[iCent] = new TGraphErrors(nStarMtBins,tmpMt,tmpRlong,
						 tmpMtErr,tmpRlongErr);
    StarRlongMt19GrErr[iCent]->SetName(Form("StarRlongMt19GrErr_%d",iCent));
    SetStarStyle(StarRlongMt19GrErr[iCent],19,iCent);
  }   

  
  double KtShift = 0.1;
  double KtStep = 0.1;
  double KtVal;
  double MtVal;


  TFile* kaonFile[2][5];
  TFile* pionFile[2][5];
  TString mKaonFileName[2][5];
  TString mPionFileName[2][5];
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      int energy;
      int centrality;
      switch(iEne) {
      case 0:
	energy = 200;
	break;
      case 1:
	energy = 19;
	break;
      default:
	energy = 200;
      };

      switch(iCent) {
      case 0:
	centrality = 8;
	break;
      case 1: 
	centrality = 7;
	break;
      case 2: 
	centrality = 6;
	break;
      case 3: 
	centrality = 5;
	break;
      case 4: 
	centrality = 4;
	break;
      default:
	centrality = 8;
      }

      //Kaons
      mKaonFileName[iEne][iCent] = "oAuAu";
      mKaonFileName[iEne][iCent] += energy;
      mKaonFileName[iEne][iCent] += "_kaons_3d_cent";
      mKaonFileName[iEne][iCent] += centrality;
      mKaonFileName[iEne][iCent] += ".root";

      //Pions
      mPionFileName[iEne][iCent] = "oAuAu";
      mPionFileName[iEne][iCent] += energy;
      mPionFileName[iEne][iCent] += "_pions_3d_cent";
      mPionFileName[iEne][iCent] += centrality;
      mPionFileName[iEne][iCent] += ".root";

      kaonFile[iEne][iCent] = TFile::Open(mKaonFileName[iEne][iCent].Data());
      pionFile[iEne][iCent] = TFile::Open(mPionFileName[iEne][iCent].Data());
    }
  }

  //Kaon projections
  TH1F* hKaonOutProj[2][5][nKaonKtBins];  //Energy, Centrality, kT
  TH1F* hKaonFitOutProj[2][5][nKaonKtBins];
  TH1F* hKaonSideProj[2][5][nKaonKtBins];
  TH1F* hKaonFitSideProj[2][5][nKaonKtBins];
  TH1F* hKaonLongProj[2][5][nKaonKtBins];
  TH1F* hKaonFitLongProj[2][5][nKaonKtBins];

  //Pion projections
  TH1F* hPionOutProj[2][5][nPionKtBins]; //Energy, Centrality, kT
  TH1F* hPionFitOutProj[2][5][nPionKtBins];
  TH1F* hPionSideProj[2][5][nPionKtBins];
  TH1F* hPionFitSideProj[2][5][nPionKtBins];
  TH1F* hPionLongProj[2][5][nPionKtBins];
  TH1F* hPionFitLongProj[2][5][nPionKtBins];

  //Kaon TGraphErrors
  TGraphErrors *mKaonLamKtGr[2][5][3]; //Energy, Centrality, kT
  TGraphErrors *mKaonRoutKtGr[2][5][3];
  TGraphErrors *mKaonRsideKtGr[2][5][3];
  TGraphErrors *mKaonRlongKtGr[2][5][3];

  TGraphErrors *mKaonLamMtGr[2][5][3]; //Energy, Centrality, kT
  TGraphErrors *mKaonRoutMtGr[2][5][3];
  TGraphErrors *mKaonRsideMtGr[2][5][3];
  TGraphErrors *mKaonRlongMtGr[2][5][3];

  //Pion TGraphErrors
  TGraphErrors *mPionLamKtGr[2][5][3]; //Energy, Centrality, kT
  TGraphErrors *mPionRoutKtGr[2][5][3];
  TGraphErrors *mPionRsideKtGr[2][5][3];
  TGraphErrors *mPionRlongKtGr[2][5][3];

  TGraphErrors *mPionLamMtGr[2][5][3]; //Energy, Centrality, kT
  TGraphErrors *mPionRoutMtGr[2][5][3];
  TGraphErrors *mPionRsideMtGr[2][5][3];
  TGraphErrors *mPionRlongMtGr[2][5][3];

  //Read TGraphErrors
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      for(int iCharge=0; iCharge<3; iCharge++) {

	//Kaons
	mKaonLamKtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("LamKtGrErr_%d",iCharge));
	mKaonLamKtGr[iEne][iCent][iCharge]->SetName(Form("mKaonLamKtGr_%d_%d_%d",iEne,iCent,iCharge));
	mKaonRoutKtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("RoutKtGrErr_%d",iCharge));
	mKaonRoutKtGr[iEne][iCent][iCharge]->SetName(Form("mKaonRoutKtGr_%d_%d_%d",iEne,iCent,iCharge));
	mKaonRsideKtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("RsideKtGrErr_%d",iCharge));
	mKaonRsideKtGr[iEne][iCent][iCharge]->SetName(Form("mKaonRsideKtGr_%d_%d_%d",iEne,iCent,iCharge));
	mKaonRlongKtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("RlongKtGrErr_%d",iCharge));
	mKaonRlongKtGr[iEne][iCent][iCharge]->SetName(Form("mKaonRlongKtGr_%d_%d_%d",iEne,iCent,iCharge));

	mKaonLamMtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("LamMtGrErr_%d",iCharge));
	mKaonLamMtGr[iEne][iCent][iCharge]->SetName(Form("mKaonLamMtGr_%d_%d_%d",iEne,iCent,iCharge));
	mKaonRoutMtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("RoutMtGrErr_%d",iCharge));
	mKaonRoutMtGr[iEne][iCent][iCharge]->SetName(Form("mKaonRoutMtGr_%d_%d_%d",iEne,iCent,iCharge));
	mKaonRsideMtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("RsideMtGrErr_%d",iCharge));
	mKaonRsideMtGr[iEne][iCent][iCharge]->SetName(Form("mKaonRsideMtGr_%d_%d_%d",iEne,iCent,iCharge));
	mKaonRlongMtGr[iEne][iCent][iCharge] = (TGraphErrors*)kaonFile[iEne][iCent]->Get(Form("RlongMtGrErr_%d",iCharge));
	mKaonRlongMtGr[iEne][iCent][iCharge]->SetName(Form("mKaonRlongMtGr_%d_%d_%d",iEne,iCent,iCharge));

	//Pions
	mPionLamKtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("LamKtGrErr_%d",iCharge));
	mPionLamKtGr[iEne][iCent][iCharge]->SetName(Form("mPionLamKtGr_%d_%d_%d",iEne,iCent,iCharge));
	mPionRoutKtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("RoutKtGrErr_%d",iCharge));
	mPionRoutKtGr[iEne][iCent][iCharge]->SetName(Form("mPionRoutKtGr_%d_%d_%d",iEne,iCent,iCharge));
	mPionRsideKtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("RsideKtGrErr_%d",iCharge));
	mPionRsideKtGr[iEne][iCent][iCharge]->SetName(Form("mPionRsideKtGr_%d_%d_%d",iEne,iCent,iCharge));
	mPionRlongKtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("RlongKtGrErr_%d",iCharge));
	mPionRlongKtGr[iEne][iCent][iCharge]->SetName(Form("mPionRlongKtGr_%d_%d_%d",iEne,iCent,iCharge));

	mPionLamMtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("LamMtGrErr_%d",iCharge));
	mPionLamMtGr[iEne][iCent][iCharge]->SetName(Form("mPionLamMtGr_%d_%d_%d",iEne,iCent,iCharge));
	mPionRoutMtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("RoutMtGrErr_%d",iCharge));
	mPionRoutMtGr[iEne][iCent][iCharge]->SetName(Form("mPionRoutMtGr_%d_%d_%d",iEne,iCent,iCharge));
	mPionRsideMtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("RsideMtGrErr_%d",iCharge));
	mPionRsideMtGr[iEne][iCent][iCharge]->SetName(Form("mPionRsideMtGr_%d_%d_%d",iEne,iCent,iCharge));
	mPionRlongMtGr[iEne][iCent][iCharge] = (TGraphErrors*)pionFile[iEne][iCent]->Get(Form("RlongMtGrErr_%d",iCharge));
	mPionRlongMtGr[iEne][iCent][iCharge]->SetName(Form("mPionRlongMtGr_%d_%d_%d",iEne,iCent,iCharge));
      }
    }
  }

  //Read kaon projections
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      for(int iKt=0; iKt<nKaonKtBins; iKt++) {
	//Out
	hKaonOutProj[iEne][iCent][iKt] = (TH1F*)kaonFile[iEne][iCent]->Get(Form("hOutProj_2_%d",iKt));
	hKaonOutProj[iEne][iCent][iKt]->SetName(Form("hKaonOutProj_%d_%d_%d",iEne,iCent,iKt));
	hKaonFitOutProj[iEne][iCent][iKt] = (TH1F*)kaonFile[iEne][iCent]->Get(Form("hFitOutProj_2_%d",iKt));
	hKaonFitOutProj[iEne][iCent][iKt]->SetName(Form("hKaonFitOutProj_%d_%d_%d",iEne,iCent,iKt));

	//Side
	hKaonSideProj[iEne][iCent][iKt] = (TH1F*)kaonFile[iEne][iCent]->Get(Form("hSideProj_2_%d",iKt));
	hKaonSideProj[iEne][iCent][iKt]->SetName(Form("hKaonSideProj_%d_%d_%d",iEne,iCent,iKt));
	hKaonFitSideProj[iEne][iCent][iKt] = (TH1F*)kaonFile[iEne][iCent]->Get(Form("hFitSideProj_2_%d",iKt));
	hKaonFitSideProj[iEne][iCent][iKt]->SetName(Form("hKaonFitSideProj_%d_%d_%d",iEne,iCent,iKt));

	//Long
	hKaonLongProj[iEne][iCent][iKt] = (TH1F*)kaonFile[iEne][iCent]->Get(Form("hLongProj_2_%d",iKt));
	hKaonLongProj[iEne][iCent][iKt]->SetName(Form("hKaonLongProj_%d_%d_%d",iEne,iCent,iKt));
	hKaonFitLongProj[iEne][iCent][iKt] = (TH1F*)kaonFile[iEne][iCent]->Get(Form("hFitLongProj_2_%d",iKt));
	hKaonFitLongProj[iEne][iCent][iKt]->SetName(Form("hKaonFitLongProj_%d_%d_%d",iEne,iCent,iKt));
      }
    }
  }

  //Read pion projections
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      for(int iKt=0; iKt<nPionKtBins; iKt++) {
	//Out
	hPionOutProj[iEne][iCent][iKt] = (TH1F*)pionFile[iEne][iCent]->Get(Form("hOutProj_2_%d",iKt));
	hPionOutProj[iEne][iCent][iKt]->SetName(Form("hPionOutProj_%d_%d_%d",iEne,iCent,iKt));
	hPionFitOutProj[iEne][iCent][iKt] = (TH1F*)pionFile[iEne][iCent]->Get(Form("hFitOutProj_2_%d",iKt));
	hPionFitOutProj[iEne][iCent][iKt]->SetName(Form("hPionFitOutProj_%d_%d_%d",iEne,iCent,iKt));

	//Side
	hPionSideProj[iEne][iCent][iKt] = (TH1F*)pionFile[iEne][iCent]->Get(Form("hSideProj_2_%d",iKt));
	hPionSideProj[iEne][iCent][iKt]->SetName(Form("hPionSideProj_%d_%d_%d",iEne,iCent,iKt));
	hPionFitSideProj[iEne][iCent][iKt] = (TH1F*)pionFile[iEne][iCent]->Get(Form("hFitSideProj_2_%d",iKt));
	hPionFitSideProj[iEne][iCent][iKt]->SetName(Form("hPionFitSideProj_%d_%d_%d",iEne,iCent,iKt));

	//Long
	hPionLongProj[iEne][iCent][iKt] = (TH1F*)pionFile[iEne][iCent]->Get(Form("hLongProj_2_%d",iKt));
	hPionLongProj[iEne][iCent][iKt]->SetName(Form("hPionLongProj_%d_%d_%d",iEne,iCent,iKt));
	hPionFitLongProj[iEne][iCent][iKt] = (TH1F*)pionFile[iEne][iCent]->Get(Form("hFitLongProj_2_%d",iKt));
	hPionFitLongProj[iEne][iCent][iKt]->SetName(Form("hPionFitLongProj_%d_%d_%d",iEne,iCent,iKt));
      }
    }
  }

  double margin = 0.;
  //Plot pion projections and fits
  TCanvas *PionProjCanvas[2][5];
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      PionProjCanvas[iEne][iCent] = new TCanvas(Form("PionProjCanvas_%d_%d",iEne,iCent),
						Form("PionProjCanvas_%d_%d",iEne,iCent),
						1200,800);
      PionProjCanvas[iEne][iCent]->Divide(nPionKtBins, 3, margin, margin);
    }
  }

  TLatex latex;
  latex.SetTextSize(0.05);

  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      for(int iProj=0; iProj<3; iProj++) {
	for(int iKt=0; iKt<nPionKtBins; iKt++) {

	  int iPad = iProj*nPionKtBins + iKt + 1;
	  PionProjCanvas[iEne][iCent]->cd(iPad);
	  double KtLo=KtShift + PionFitKtArray[iKt][0]*KtStep;
	  double KtHi=KtShift + PionFitKtArray[iKt][1]*KtStep+KtStep;

	  if(iProj==0) { //Out
	    hPionOutProj[iEne][iCent][iKt]->Draw();
	    hPionFitOutProj[iEne][iCent][iKt]->Draw("sameC");
	    latex.DrawLatex(0.02, 0.94, "STAR preliminary");
	    latex.DrawLatex(0.02, 1.4, Form("%3.2f#leqk_{T}#leq%3.2f",KtLo,KtHi) );
	    latex.DrawLatex(0.1, 1.2, Form("Out Projection") );
	    
	  }
	  else if(iProj==1) { //Side
	    hPionSideProj[iEne][iCent][iKt]->Draw();
	    hPionFitSideProj[iEne][iCent][iKt]->Draw("sameC");
	    latex.DrawLatex(0.02,0.94,"STAR preliminary");
	    latex.DrawLatex(0.1, 1.2, Form("Side Projection") );
	  }
	  else if(iProj==2) { //Long
	    hPionLongProj[iEne][iCent][iKt]->Draw();
	    hPionFitLongProj[iEne][iCent][iKt]->Draw("sameC");
	    latex.DrawLatex(0.02,0.94,"STAR preliminary");
	    latex.DrawLatex(0.1, 1.2, Form("Long Projection") );
	  }      
	}
      }
    }
  }

  //Plot kaon projections and fits
  TCanvas *KaonProjCanvas[2][5];
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      KaonProjCanvas[iEne][iCent] = new TCanvas(Form("KaonProjCanvas_%d_%d",iEne,iCent),
						Form("KaonProjCanvas_%d_%d",iEne,iCent),
						1200,800);
      KaonProjCanvas[iEne][iCent]->Divide(nKaonKtBins, 3, margin, margin);
    }
  }

  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      for(int iProj=0; iProj<3; iProj++) {
	for(int iKt=0; iKt<nKaonKtBins; iKt++) {

	  int iPad = iProj*nKaonKtBins + iKt + 1;
	  KaonProjCanvas[iEne][iCent]->cd(iPad);
	  double KtLo=KtShift + KaonFitKtArray[iKt][0]*KtStep;
	  double KtHi=KtShift + KaonFitKtArray[iKt][1]*KtStep+KtStep;

	  if(iProj==0) { //Out
	    hKaonOutProj[iEne][iCent][iKt]->Draw();
	    hKaonFitOutProj[iEne][iCent][iKt]->Draw("sameC");
	    latex.DrawLatex(0.02, 0.94, "STAR preliminary");
	    latex.DrawLatex(0.02, 1.4, Form("%3.2f#leqk_{T}#leq%3.2f",KtLo,KtHi) );
	    latex.DrawLatex(0.1, 1.3, Form("Out Projection") );
	  }
	  else if(iProj==1) { //Side
	    hKaonSideProj[iEne][iCent][iKt]->Draw();
	    hKaonFitSideProj[iEne][iCent][iKt]->Draw("sameC");
	    latex.DrawLatex(0.02,0.94,"STAR preliminary");
	    latex.DrawLatex(0.1, 1.3, Form("Side Projection") );
	  }
	  else if(iProj==2) { //Long
	    hKaonLongProj[iEne][iCent][iKt]->Draw();
	    hKaonFitLongProj[iEne][iCent][iKt]->Draw("sameC");
	    latex.DrawLatex(0.02,0.94,"STAR preliminary");
	    latex.DrawLatex(0.1, 1.3, Form("Long Projection") );
	  }      
	}
      }
    }
  }

  /////////////////////////////////////////////
  //               kT dependence             //
  /////////////////////////////////////////////
  TMultiGraph* mRoutKtMGr[2][5];   //iEnergy,Centrality
  TMultiGraph* mRsideKtMGr[2][5];  //iEnergy,Centrality
  TMultiGraph* mRlongKtMGr[2][5];  //iEnergy,Centrality
  TMultiGraph* mLamKtMGr[2][5];    //iEnergy,Centrality

  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {

      mRoutKtMGr[iEne][iCent] = new TMultiGraph();
      mRsideKtMGr[iEne][iCent] = new TMultiGraph();
      mRlongKtMGr[iEne][iCent] = new TMultiGraph();
      mLamKtMGr[iEne][iCent] = new TMultiGraph();

      //for(int iCent=0; iCent<nStar200CentBins; iCent++) {
      //Out
      SetTGraphStyle(mKaonRoutKtGr[iEne][iCent][2],4,iCent);
      mRoutKtMGr[iEne][iCent]->Add(mKaonRoutKtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionRoutKtGr[iEne][iCent][2],3,iCent);
      mRoutKtMGr[iEne][iCent]->Add(mPionRoutKtGr[iEne][iCent][2]);
      //Side
      SetTGraphStyle(mKaonRsideKtGr[iEne][iCent][2],4,iCent);
      mRsideKtMGr[iEne][iCent]->Add(mKaonRsideKtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionRsideKtGr[iEne][iCent][2],3,iCent);
      mRsideKtMGr[iEne][iCent]->Add(mPionRsideKtGr[iEne][iCent][2]);
      //Long
      SetTGraphStyle(mKaonRlongKtGr[iEne][iCent][2],4,iCent);
      mRlongKtMGr[iEne][iCent]->Add(mKaonRlongKtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionRlongKtGr[iEne][iCent][2],3,iCent);
      mRlongKtMGr[iEne][iCent]->Add(mPionRlongKtGr[iEne][iCent][2]);
      //Lambda
      SetTGraphStyle(mKaonLamKtGr[iEne][iCent][2],4,iCent);
      mLamKtMGr[iEne][iCent]->Add(mKaonLamKtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionLamKtGr[iEne][iCent][2],3,iCent);
      mLamKtMGr[iEne][iCent]->Add(mPionLamKtGr[iEne][iCent][2]);
      //}
    }
  }

  TCanvas *KtDepCanv[2];
  for(int iEne=0; iEne<nEnergies; iEne++) {
    KtDepCanv[iEne]= new TCanvas(Form("KtDepCanv_%d",iEne),
				 Form("KtDepCanv_%d",iEne),
				 mCanvSize[0],mCanvSize[1]);
    KtDepCanv[iEne]->Divide(4, nStar200CentBins, 0.0, 0.0);
  }

  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {

      KtDepCanv[iEne]->cd(iCent*4+1);
      gPad->SetBorderSize(2);
      mLamKtMGr[iEne][iCent]->Draw("AP");
      mLamKtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mLamRange[0],mLamRange[1]);
      mLamKtMGr[iEne][iCent]->GetYaxis()->SetTitle("#lambda");
      mLamKtMGr[iEne][iCent]->GetXaxis()->SetTitle("k_{T} (GeV/c)");
      mLamKtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mLamKtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      //latex.DrawLatex(0.4,0.3, "Red symbols : 0-5%   centrality");
      //latex.DrawLatex(0.4,0.1, "Blue symbols: 20-30% centrality");

      KtDepCanv[iEne]->cd(iCent*4+2);
      gPad->SetBorderSize(2);
      mRoutKtMGr[iEne][iCent]->Draw("AP");
      mRoutKtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mRadiiRange[0],mRadiiRange[1]);
      mRoutKtMGr[iEne][iCent]->GetYaxis()->SetTitle("R_{out} (fm)");
      mRoutKtMGr[iEne][iCent]->GetXaxis()->SetTitle("k_{T} (GeV/c)");
      mRoutKtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mRoutKtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      latex.DrawLatex(0.4,1.5,"STAR preliminary");
      TLegend *leg1 = new TLegend(0.65,0.65, 0.85, 0.85);
      leg1->AddEntry(mKaonRoutKtGr[iEne][0][2],"KK","p");
      leg1->AddEntry(mPionRoutKtGr[iEne][0][2],Form("#pi#pi"),"p");
      leg1->SetLineColor(0);
      leg1->Draw("same");

      KtDepCanv[iEne]->cd(iCent*4+3);
      gPad->SetBorderSize(2);
      mRsideKtMGr[iEne][iCent]->Draw("AP");
      mRsideKtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mRadiiRange[0],mRadiiRange[1]);
      mRsideKtMGr[iEne][iCent]->GetYaxis()->SetTitle("R_{side} (fm)");
      mRsideKtMGr[iEne][iCent]->GetXaxis()->SetTitle("k_{T} (GeV/c)");
      mRsideKtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mRsideKtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      latex.DrawLatex(0.4, 1.5,"STAR preliminary");
      if(iEne==0) {
	latex.DrawLatex(0.4, 6., "Au+Au #sqrt{S_{NN}}=200 GeV");
      }
      else if(iEne==1) {
	latex.DrawLatex(0.4, 6., "Au+Au #sqrt{S_{NN}}=19.6 GeV");
      }

      KtDepCanv[iEne]->cd(iCent*4+4);
      gPad->SetBorderSize(2);
      mRlongKtMGr[iEne][iCent]->Draw("AP");
      mRlongKtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mRadiiRange[0],mRadiiRange[1]);
      mRlongKtMGr[iEne][iCent]->GetYaxis()->SetTitle("R_{long} (fm)");
      mRlongKtMGr[iEne][iCent]->GetXaxis()->SetTitle("k_{T} (GeV/c)");
      mRlongKtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mRlongKtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      latex.DrawLatex(0.4,1.5,"STAR preliminary");
    }
  }

  /////////////////////////////////////////////
  //               mT dependence             //
  /////////////////////////////////////////////

  TMultiGraph* mRoutMtMGr[2][5];   //iEnergy,Centrality
  TMultiGraph* mRsideMtMGr[2][5];  //iEnergy,Centrality
  TMultiGraph* mRlongMtMGr[2][5];  //iEnergy,Centrality
  TMultiGraph* mLamMtMGr[2][5];    //iEnergy,Centrality

  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {

      mRoutMtMGr[iEne][iCent] = new TMultiGraph();
      mRsideMtMGr[iEne][iCent] = new TMultiGraph();
      mRlongMtMGr[iEne][iCent] = new TMultiGraph();
      mLamMtMGr[iEne][iCent] = new TMultiGraph();
      
      //for(int iCent=0; iCent<nStar200CentBins; iCent++) {

      //Out
      SetTGraphStyle(mKaonRoutMtGr[iEne][iCent][2],4,iCent);
      mRoutMtMGr[iEne][iCent]->Add(mKaonRoutMtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionRoutMtGr[iEne][iCent][2],3,iCent);
      mRoutMtMGr[iEne][iCent]->Add(mPionRoutMtGr[iEne][iCent][2]);
      if(iEne == 0) {
	mRoutMtMGr[iEne][iCent]->Add(StarRoutMt200GrErr[iCent]);
      }
      else {
	mRoutMtMGr[iEne][iCent]->Add(StarRoutMt19GrErr[iCent]);
      }
      //Side
      SetTGraphStyle(mKaonRsideMtGr[iEne][iCent][2],4,iCent);
      mRsideMtMGr[iEne][iCent]->Add(mKaonRsideMtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionRsideMtGr[iEne][iCent][2],3,iCent);
      mRsideMtMGr[iEne][iCent]->Add(mPionRsideMtGr[iEne][iCent][2]);
      if(iEne == 0) {
	mRsideMtMGr[iEne][iCent]->Add(StarRsideMt200GrErr[iCent]);
      }
      else {
	mRsideMtMGr[iEne][iCent]->Add(StarRsideMt19GrErr[iCent]);
      }
      //Long
      SetTGraphStyle(mKaonRlongMtGr[iEne][iCent][2],4,iCent);
      mRlongMtMGr[iEne][iCent]->Add(mKaonRlongMtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionRlongMtGr[iEne][iCent][2],3,iCent);
      mRlongMtMGr[iEne][iCent]->Add(mPionRlongMtGr[iEne][iCent][2]);
      if(iEne == 0) {
	mRlongMtMGr[iEne][iCent]->Add(StarRlongMt200GrErr[iCent]);
      }
      else {
	mRlongMtMGr[iEne][iCent]->Add(StarRlongMt19GrErr[iCent]);
      }
      //Lambda
      SetTGraphStyle(mKaonLamMtGr[iEne][iCent][2],4,iCent);
      mLamMtMGr[iEne][iCent]->Add(mKaonLamMtGr[iEne][iCent][2]);
      SetTGraphStyle(mPionLamMtGr[iEne][iCent][2],3,iCent);
      mLamMtMGr[iEne][iCent]->Add(mPionLamMtGr[iEne][iCent][2]);
    }
  }

  TCanvas *MtDepCanv[2];
  for(int iEne=0; iEne<nEnergies; iEne++) {
    MtDepCanv[iEne]= new TCanvas(Form("MtDepCanv_%d",iEne),
				 Form("MtDepCanv_%d",iEne),
				 mCanvSize[0],mCanvSize[1]);
    MtDepCanv[iEne]->Divide(4,nStar200CentBins,0.,0.);

  }

  MtDepCanv[0]->Update();
  TString str;
  latex.SetTextSize(0.1);
  char mIndex = 'a';
  for(int iEne=0; iEne<nEnergies; iEne++) {
    for(int iCent=0; iCent<nStar200CentBins; iCent++) {

      MtDepCanv[iEne]->cd(iCent*4+1);
      gPad->SetTicks(1,1);
      mLamMtMGr[iEne][iCent]->Draw("AP");
      mLamMtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mLamRange[0],mLamRange[1]);
      mLamMtMGr[iEne][iCent]->GetXaxis()->SetLimits(mMtRange[0],mMtRange[1]);
      mLamMtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      mLamMtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      
      switch(iCent) {
      case 0:
	str = "0-5%";
	break;
      case 1:
	str = "5-10%";
	break;
      case 2:
	str = "10-20%";
	break;
      case 3:
	str = "20-30%";
	break;
      case 4:
	str = "30-40%";
	break;
      default:
	str = "unknown";
	break;
      }

      latex.SetTextSize(0.1);
      latex.DrawLatex(0.4,0.8, str.Data());
      latex.SetTextSize(0.1);

      switch(iCent) {
      case 0:
	str = "(a)";
	break;
      case 1:
	str = "(e)";
	break;
      case 2:
	str = "(i)";
	break;
      case 3:
	str = "(m)";
	break;
      case 4:
	str = "(q)";
	break;
      define:
	str = "(?)";
      };
      latex.DrawLatex(1.1,0.2,str.Data());

      if(iCent == 2) {
	TLegend *leg2 = new TLegend(0.3,0.3,0.73, 0.4);
	leg2->SetTextSize(0.08);
	leg2->AddEntry(StarRsideMt200GrErr[iCent],Form("#pi#pi (STAR)"),"p");
	leg2->SetLineColor(0);
	leg2->Draw("same");
      }

      if(iCent == 0) {
	TLegend *leg1 = new TLegend(0.35,0.08,0.7, 0.4);
	leg1->SetTextSize(0.15);
	leg1->AddEntry(mKaonRoutMtGr[iEne][0][2],"KK","p");
	leg1->AddEntry(mPionRoutMtGr[iEne][0][2],Form("#pi#pi"),"p");
	//leg1->AddEntry(StarRoutMt200GrErr[0],Form("#pi#pi: STAR PRC92(2015)014904"),"p");
	leg1->SetLineColor(0);
	leg1->Draw("same");
      }
      else if(iCent == 2) {
	latex.SetTextAngle(90);
	latex.SetTextSize(0.2);
	latex.DrawLatex(-0.1,0.3,"#lambda");
	latex.SetTextSize(0.15);
	latex.SetTextAngle(0);
      }

      MtDepCanv[iEne]->cd(iCent*4+2);
      gPad->SetTicks(1,1);
      //gPad->SetTicky(1,1);
      mRoutMtMGr[iEne][iCent]->Draw("AP");
      mRoutMtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mRadiiRange[0],mRadiiRange[1]);
      mRoutMtMGr[iEne][iCent]->GetXaxis()->SetLimits(mMtRange[0],mMtRange[1]);
      mRoutMtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mRoutMtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      mRoutMtMGr[iEne][iCent]->GetYaxis()->SetLabelOffset(999);
      if(iCent==0) {
	latex.DrawLatex(0.9,5.8,"#mu=out");
      }

      if(iCent == 1) {
	if(iEne==0) {
	  latex.SetTextSize(0.09);
	  latex.DrawLatex(0.4, 6.0, "Au+Au #sqrt{s_{NN}}=200 GeV");
	  latex.SetTextSize(0.1);
	}
	else if(iEne==1) {
	  latex.DrawLatex(0.3, 5.5, "Au+Au #sqrt{S_{NN}}=19.6 GeV");
	}
      }

      switch(iCent) {
      case 0:
	str = "(b)";
	break;
      case 1:
	str = "(f)";
	break;
      case 2:
	str = "(j)";
	break;
      case 3:
	str = "(n)";
	break;
      case 4:
	str = "(r)";
	break;
      default:
	str = "(?)";
      };
      latex.DrawLatex(1.1,0.8,str.Data());

      MtDepCanv[iEne]->cd(iCent*4+3);
      gPad->SetTicks(1,1);
      //gPad->SetTicky(2);
      mRsideMtMGr[iEne][iCent]->Draw("AP");
      mRsideMtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mRadiiRange[0],mRadiiRange[1]);
      mRsideMtMGr[iEne][iCent]->GetXaxis()->SetLimits(mMtRange[0],mMtRange[1]);
      mRsideMtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mRsideMtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      mRsideMtMGr[iEne][iCent]->GetYaxis()->SetLabelOffset(999);
      if(iCent==0) {
	latex.DrawLatex(0.9,5.8,"#mu=side");
      }
      if(iCent == 1) {
	latex.DrawLatex(0.3, 1.5,"STAR preliminary");
      }

      switch(iCent) {
      case 0:
	str = "(c)";
	break;
      case 1:
	str = "(g)";
	break;
      case 2:
	str = "(k)";
	break;
      case 3:
	str = "(o)";
	break;
      case 4:
	str = "(s)";
	break;
      default:
	str = "(?)";
      };
      latex.DrawLatex(1.1,0.8,str.Data());

      MtDepCanv[iEne]->cd(iCent*4+4);
      //gPad->SetTicky(2);
      mRlongMtMGr[iEne][iCent]->Draw("APY+");
      gPad->SetTicks(1,1);
      mRlongMtMGr[iEne][iCent]->GetYaxis()->SetRangeUser(mRadiiRange[0],mRadiiRange[1]);
      mRlongMtMGr[iEne][iCent]->GetXaxis()->SetLimits(mMtRange[0],mMtRange[1]);
      mRlongMtMGr[iEne][iCent]->GetXaxis()->SetNdivisions(505);
      mRlongMtMGr[iEne][iCent]->GetYaxis()->SetNdivisions(505);
      if(iCent==0) {
	latex.DrawLatex(0.9,5.8,"#mu=long");
      }

      switch(iCent) {
      case 0:
	str = "(d)";
	break;
      case 1:
	str = "(h)";
	break;
      case 2:
	str = "(l)";
	break;
      case 3:
	str = "(p)";
	break;
      case 4:
	str = "(t)";
	break;
      default:
	str = "(?)";
      };
      latex.DrawLatex(1.1,0.8,str.Data());

    }
  }

  //MtDepCanv[0]->Modified();

  MtDepCanv[0]->Print("Try.png");
  
  //TFile* oFile = new TFile(Form("oResults_auau200_cent%d",cent),"recreate");
  //PionProjCanvas->Write();
  //KaonProjCanvas->Write();
}

//_________________
void SetStarStyle(TGraphErrors* grErr, int energy, int cent) {

  int mLineColor = 1;
  int mMarkerColor = 3;
  double mMarkerSize = 1.6;
  int mLineWidth = 2;
  int mMarkerStyle = 29;
  
  mLineColor = mMarkerColor;
  grErr->SetLineWidth(mLineWidth);
  grErr->SetLineColor(mLineColor);
  grErr->SetMarkerStyle(mMarkerStyle);
  grErr->SetMarkerSize(mMarkerSize);
  grErr->SetMarkerColor(mMarkerColor);
}

//_________________
void SetTGraphStyle(TGraphErrors* grErr, int pid, int cent) {
  
  int mLineColor = 1;
  int mMarkerColor;
  double mMarkerSize = 1.2;
  int mLineWidth = 2.;
  int mMarkerStyle;

  if(pid==3) {
    mMarkerColor = 2;
    mMarkerStyle = 20;
  }
  else if(pid==4) {
    mMarkerColor = 4;
    mMarkerStyle = 22;
  }
  mLineColor = mMarkerColor;
  grErr->SetLineWidth(mLineWidth);
  grErr->SetLineColor(mLineColor);
  grErr->SetMarkerSize(mMarkerSize);
  grErr->SetMarkerColor(mMarkerColor);
  grErr->SetMarkerStyle(mMarkerStyle);
}


//_________________
void SetGlobalStyle() {

  //gROOT->SetStyle("Plain");
  gStyle->SetOptStat(0);
  gStyle->SetOptTitle(0);
  gStyle->SetLineWidth(2);
  gStyle->SetPadBorderSize(2);
  gStyle->SetPadTopMargin(0.3);
  gStyle->SetPadLeftMargin(0.3);
  gStyle->SetPadRightMargin(0.3);
  gStyle->SetPadBottomMargin(0.3);
  gStyle->SetLabelSize(0.08,"X");
  gStyle->SetLabelSize(0.08,"Y");
  gStyle->SetLabelOffset(0.05,"X");
  gStyle->SetLabelOffset(0.05,"Y");
  //gStyle->SetPadTickX(1);
  //gStyle->SetPadTickY(1);
  gStyle->SetFrameLineWidth(2);
}