void plot()
{

  gROOT->ProcessLine(".x $VMCWORKDIR/gconfig/rootlogon.C");
  //  TFile *inputfile = new TFile(Form("$RECOOUT/GenE/run_%d.reco.out.root",run));

  TCut GoodMom = "mom>0.05";
  //TCut ValidFit = "fitflag==0"; 
  TCut ValidFit = "fitflag==0&&cdcflag==0";
  TCut ValidChi = "chi<5.&&chi>0.6"; 
  TCut GoodVert = "abs(cvz)>0.0001";
  TCut GoodVel = "rpcv>0.1";
  TCut PionP = "mass>0.03&&mass<0.7&&charge>0";
  //TCut PionP = "bmass>0.03&&bmass<0.7&&charge>0";
  //TCut Proton = "mass>0.8&&mass<1.25";
  TCut Proton ="bmass>0.7&&bmass<1.2";
  TCut PionN = "mass>0.03&&mass<0.8&&charge<0";
  //TCut PionN = "bmass>0.03&&bmass<0.8&&charge<0";
  TCut KaonP  = "mass>0.42&&mass<0.7&&charge>0";
  //TCut KaonP  = "bmass>0.32&&bmass<0.75&&charge>0";
  
  //TCut Target = "abs(z0)<25.";
  TCut Target = "abs(z0)<25.&&abs(d0)<2.";
  
  const Float_t ProtMass=0.938;
  const Float_t PionMass=0.140;
  const Float_t KaonMass=0.494;
  const Float_t c=29.98;
  const Int_t nPi = 8;
  const Int_t nP = 10;


  TH1F* histodiff = new TH1F("histodiff","CDC+TPC - CDC Momentum Diff. ", 500,-1.,1.);
  TH2F* histo2d = new TH2F("histo2d","CDC+TPC vs CDC Momentum",200,0.,5.,200.,0.,5.);
  TH2F* histdiff2d = new TH2F("histdiff2d","CDC+TPC - CDC Momentum Diff",
			      100,0.,2., 50,-0.5,0.5);
  TH2F* hProtonDiff = new TH2F("hProtonDiff","Combined - CDC Momentum Diff, Protons",
			      200,0.,2., 100,-0.5,0.5);
  TH2F* hPionDiff = new TH2F("hPionDiff","Combined - CDC Momentum Diff, Pions",
			      200,0.,2., 100,-0.5,0.5);
  TH2F* hdiffperc = new TH2F("hdiffperc","CDC+TPC - CDC Relative Diff", 
			     100,0.,2.,100,-1.,1.);
  TH2F* hProtonRelDiff = new TH2F("hProtonRelDiff","Combined - CDC Rel. Momentum Diff, Protons",
				  100,0.,2., 100,-0.5,0.5);
  TH2F* hPionRelDiff = new TH2F("hPionRelDiff","Combined - CDC Rel. Momentum Diff, Pions",
			     100,0.,2., 200,-0.5,0.5);
  //
  //
  TCanvas * CanvasCombinedFit = new TCanvas("CanvasCombinedFit", "Momentum Distribution",
				      100,10,800,900);

  //Double par[6];
  TF1 * g1 = new TF1("g1","gaus",-0.1,0.1);
  //TF2 * g2 = new TF1("g2","gaus",-0.05,0.05);
  TF1 * gx1 = new TF1("gx1","gaus",-1.,1.);
  TF1 * gx2 = new TF1("gx2","gaus",-1.,1.);
  TF1 * gx3 = new TF1("gx3","gaus",-0.5,0.5);
  TF1 * gy1 = new TF1("gy1","gaus",-1.,1.);
  TF1 * gy2 = new TF1("gy2","gaus",-1.,1.);
  TF1 * gy3 = new TF1("gy3","gaus",-0.5,0.5);
  TF1 * gz1 = new TF1("gz1","gaus",-5.,5.);
  TF1 * gz2 = new TF1("gz2","gaus",-5.,25.);
  TF1 * gz3 = new TF1("gz3","gaus",-5.,5.);
  //TF1 * gtot = new TF1("gtot",gaus(0) + gaus(3),-0.2,0.2);

  //gStyle->SetPalette(8);
  CanvasCombinedFit->Divide(1,3);
  CanvasCombinedFit->cd(1);
  dispTrT->Draw("mom:momb>>histo2d");
  histo2d->GetXaxis()->SetTitle("CDC Momentum [GeV]");
  histo2d->GetYaxis()->SetTitle("TPC+CDC Momentum [GeV]");
  CanvasCombinedFit->cd(2);
  dispTrT->Draw("mom-momb>>histodiff",ValidFit);
  histodiff->Draw();
  g1->SetLineColor(2);
  histodiff->Fit(g1,"R");
  //histodiff->Fit(g2,"R+");
  //g1->GetParameters(&par[0]);
  //g2->GetParameters(&par[3]);  //histodiff->Fit(gtot,"R");
  CanvasCombinedFit->cd(3);
  dispTrT->Draw("mom-momb:momb>>histdiff2d",ValidFit);
  histdiff2d->Draw();
  histdiff2d->GetXaxis()->SetTitle("CDC Momentum [GeV]");
  histdiff2d->GetYaxis()->SetTitle("Momentum Diff. [Gev]");

  histodiff->GetXaxis()->SetTitle("Momentum Diff. [Gev]");
  histodiff->GetYaxis()->SetTitle("Counts");

  TH1F* VertexX = new TH1F("VertexX","X Vertex ", 200,-5.,5.);
  TH1F* VertexCdcX = new TH1F("VertexCdcX","X Vertex ", 200,-5.,5.);
  TH1F* TpcVertexX = new TH1F("TpcVertexX","X Vertex ", 200,-5.,5.);
  TH1F* VertexY = new TH1F("VertexY","Y Vertex ", 200,-5.,5.);
  TH1F* VertexCdcY = new TH1F("VertexCdcY","Y Vertex ", 200,-5.,5.);
  TH1F* TpcVertexY = new TH1F("TpcVertexY","Y Vertex ", 200,-5.,5.);
  TH1F* VertexZ = new TH1F("VertexZ","Z Vertex ", 200,-50.,50.);
  TH1F* VertexCdcZ = new TH1F("VertexCdcZ","Z Vertex ", 200,-50.,50.);
  TH1F* TpcVertexZ = new TH1F("TpcVertexZ","Z Vertex ", 200,-50.,50.);
  //
  TF1 * gammabetaP = new TF1("gammabetaP","(0.938/29.98)*x/sqrt(1-pow(x/29.98,2))",
			     1.,29.);
  TF1 * gammabetaPi = new TF1("gammabetaPi","(0.14/29.98)*x/sqrt(1-pow(x/29.98,2))",
			      1.,29.);
   TF1 * gammabetaPiN = new TF1("gammabetaPiN","-(0.14/29.98)*x/sqrt(1-pow(x/29.98,2))",
			      1.,29.);
  TF1 * gammabetaK = new TF1("gammabetaK","(0.494/29.98)*x/sqrt(1-pow(x/29.98,2))",
			     1.,29.);
  TF1 * gammabetaKN = new TF1("gammabetaKN","-(0.494/29.98)*x/sqrt(1-pow(x/29.98,2))",
			     1.,29.);
  TF1 * gammabetaD = new TF1("gammabetaD","(1.876/29.98)*x/sqrt(1-pow(x/29.98,2))",
			     1.,29.);

  /*
  TCanvas * CanvasVertex = new TCanvas("CanvasVertex", "Vertex Distribution",
				       100,10,600,800);
 
  TF1* gaussCdcX = new TF1("gaussCdcX","gaus",-1.,1.);
  gaussCdcX->SetParameter(0,0.055);
  gaussCdcX->SetParameter(1,0.016);
  gaussCdcX->SetParameter(2,0.314);
  //
  TF1* gaussCdcY = new TF1("gaussCdcY","gaus",-1.,1.);
  gaussCdcY->SetParameter(0,0.062);
  gaussCdcY->SetParameter(1,0.190);
  gaussCdcY->SetParameter(2,0.278);
  //
  TF1* gaussCdcZ = new TF1("gaussCdcZ","gaus",-10.,10.);
  gaussCdcZ->SetParameter(0,0.0185);
  gaussCdcZ->SetParameter(1,1.75);
  gaussCdcZ->SetParameter(2,5.07);
  //
  TF1* gaussVertZ = new TF1("gaussVertZ","gaus",-10.,10.);
  gaussVertZ->SetParameter(0,0.046); //0.053 for "first-step" align
  gaussVertZ->SetParameter(1,2.367); //0.106 for "first-step" align
  gaussVertZ->SetParameter(2,2.05); //1.68 for "first-step" align
  //
  CanvasVertex->Divide(1,3);
  CanvasVertex->cd(1);
  dispEvT->Draw("cvx>>VertexCdcX",GoodVert);
  Int_t evts = VertexCdcX->GetSum();
  dispTrT->Draw("tvx>>TpcVertexX",ValidFit);
  dispTrT->Draw("vx>>VertexX",ValidFit);
  Int_t counts = VertexX->GetSum();
  VertexX->Scale(1./counts);
  VertexX->SetLineColor(4);
  VertexX->SetAxisRange(-2.5,2.5,"X");
  VertexX->Draw("");
  VertexX->Fit("gx3","R");
  TpcVertexX->Scale(1./counts);
  //TpcVertexX->Fit(gx2,"R");
  TpcVertexX->Draw("same");
  //TpcVertexX->Fit(gx2,"R+");
  VertexCdcX->Scale(1./evts);
  VertexCdcX->SetLineColor(2);
  VertexCdcX->Draw("same");
  //gaussCdcX->Draw("same");
  //VertexCdcX->Fit(gx1,"R+");
  VertexX->GetXaxis()->SetTitle("x [cm]");
  VertexX->GetYaxis()->SetTitle("Norm. Counts");
  //dispEvT->Draw("cvx");
  CanvasVertex->cd(2);
  dispTrT->Draw("tvy>>TpcVertexY",ValidFit);
  dispTrT->Draw("vy>>VertexY",ValidFit);
  dispEvT->Draw("cvy>>VertexCdcY",GoodVert);
  VertexY->Scale(1./counts);
  VertexY->SetLineColor(4);
  VertexY->SetAxisRange(-2.5,2.5,"X");
  VertexY->Draw("");
  VertexY->Fit("gy3","R");
  TpcVertexY->Scale(1./counts);
  TpcVertexY->Draw("same");
  //TpcVertexY->Fit(gy2,"R+");
  VertexCdcY->Scale(1./evts);
  VertexCdcY->SetLineColor(2);
  VertexCdcY->Draw("same");
  //gaussCdcY->Draw("same");
  //VertexCdcY->Fit(gy1,"R+");
  VertexY->GetXaxis()->SetTitle("y [cm]");
  VertexY->GetYaxis()->SetTitle("Norm. Counts");
  CanvasVertex->cd(3);
  dispTrT->Draw("vz>>VertexZ",ValidFit);
  dispTrT->Draw("tvz>>TpcVertexZ",ValidFit);
  dispEvT->Draw("cvz>>VertexCdcZ",GoodVert);
  TpcVertexZ->Scale(1./counts);
  TpcVertexZ->SetAxisRange(-10.,40.,"X");
  TpcVertexZ->Draw("");
  TpcVertexZ->Fit(gz2,"R");
  VertexZ->Scale(1./counts);
  VertexZ->SetLineColor(4);
  VertexZ->Draw("same");
  //gaussVertZ->Draw("same");
  VertexCdcZ->Scale(1./evts);
  VertexCdcZ->SetLineColor(2);
  VertexCdcZ->Draw("same");
  //gaussCdcZ->Draw("same");
  //VertexCdcZ->Fit(gz1,"R");
  //VertexZ->Fit(gz3,"R");
  TpcVertexZ->GetXaxis()->SetTitle("z [cm]");
  TpcVertexZ->GetYaxis()->SetTitle("Norm. Counts");
  */

  TProfile * hProf = new TProfile("hProf","TPC+CDC - CDC Momentum vs CDC Momentum",
				  10,0.,1., -0.1,0.1); 
  TProfile * ProfRPC = new TProfile("ProfRPC","Valid RPC",
				  10,0.,1., -0.15,0.25); 
  TProfile* ProfProtCdc = new TProfile("ProfProtCdc","CDC, Proton",
				       10,0.250,1.250,-0.25,0.25);
  TProfile* ProfProtGf = new TProfile("ProfProtGf","GF, Proton",
				       10,0.250,1.250,-0.25,0.25);
  TProfile* ProfProtComb = new TProfile("ProfProtComb","Combined, Proton",
				       10,0.250,1.250,-0.25,0.25);
  //

  TH2F* MomVsVel = new TH2F("MomVsVel","Combined Momentum vs RPC Velocity",175,0.,35.,75,0.05,1.55);
  TH2F* MomBvsVel = new TH2F("MomBvsVel","CDC Momentum vs RPC Velocity",175,0.,35.,75,0.05,1.55);

  TH2F* ChMomVsVel = new TH2F("ChMomVsVel","Combined Momentum vs RPC Velocity",175,0.,35.,200,-2.,2.);
  TH2F* ChMomVsVelK = new TH2F("ChMomVsVelK","Combined Momentum vs RPC Velocity, kaon mass",175,0.,35.,200,-2.,2.);
  TH2F* ChMomBvsVel = new TH2F("ChMomBvsVel","CDC Momentum vs RPC Velocity",175,0.,35.,200,-2.,2.);
  TH2F* ChMomBvsVelK = new TH2F("ChMomBvsVelK","CDC Momentum vs RPC Velocity, kaon mass",175,0.,35.,200,-2.,2.);
  TH2F* ChMomgfVsVelK = new TH2F("ChMomgfVsVelK","GF Momentum vs RPC vel, kaon",175,0.,35.,200,-2.,2.);
  TH2F* GfCorr = new TH2F("GfCorr","GF Momentum vs RPC Velocity",175,0.,35.,200,-2.,2);
  //
  //
  /*
    TH2F * histoPComb = new TH2F("histoPComb","Proton, Combined",
                                120,0.05,1.25,100,-1.,1.);
    TH2F * histoPCdc = new TH2F("histoPCdc","Proton, CDC ",
                                120,0.05,1.25,100,-1.,1.);
  */
  TH2F * histoPComb = new TH2F("histoPComb","Proton, Combined",
   			     24,0.025,1.225,800,-1.,1.);
  TH2F * histoPCdc = new TH2F("histoPCdc","Proton, CDC ",
  			      24,0.025,1.225,800,-1.,1.);
  TH2F * histoPGf = new TH2F("histoPGf","Proton, GENFIT",
  			      24,0.025,1.225,800,-1.,1.);
  //
  //
  TH2F * histoPiPComb = new TH2F("histoPiPComb","Pion+, Combined",
			      20,0.025,1.025,400,-1.,1.);
  TH2F * histoPiPCdc = new TH2F("histoPiPCdc","Pion+, CDC ",
			       20,0.025,1.025,400,-1.,1.);
  TH2F * histoPiPGf = new TH2F("histoPiPGf","Pion+, GENFIT",
			       20,0.025,1.025,400,-1.,1.);
  //
  TH2F * histoPiNComb = new TH2F("histoPiNComb","Pion-, Combined",
			       20,0.025,1.025,400,-1.,1.);
  TH2F * histoPiNCdc = new TH2F("histoPiNCdc","Pion-, CDC ",
				20,0.025,1.025,400,-1.,1.);
  TH2F * histoPiNGf = new TH2F("histoPiNGf","Pion-, GENFIT",
				20,0.025,1.025,400,-1.,1.);
  //
  TH2F * histoKComb = new TH2F("histoKComb","Kaon, Combined",
			     24,0.025,1.225,50,-1.,1.);
  TH2F * histoKCdc = new TH2F("histoKCdc","Kaon, CDC",
			      24,0.025,1.225,50,-1.,1.);
  TH2F * histoKGf = new TH2F("histoKGf","Kaon, GENFIT",
			      24,0.025,1.225,50,-1.,1.);
  //
  //
  //  TH2F * histoPiPCombT = new TH2F("histoPiPCombT","Pion, Combined Test",
  //		     50,0.05,1.05,200,-2.,2.);
  //TH2F * histoPiPCdcT = new TH2F("histoPiPCdcT","Pion, CDC Test",
  //			     50,0.05,1.05,200,-2.,2.);

  //  TH2F * histoPCombT = new TH2F("histoPCombT","Proton, Combined Test",
  //			     100,0.05,1.05,100,-1.,1.);
  //TH2F * histoPiPCdcT = new TH2F("histoPiCdcT","Pion, CDC Test",
  //			     100,0.05,1.05,100,-1.,1.);
  

  TH1F* hMass = new TH1F("hMass","CDC Mass",400,0.,4.);

  TH1F* hChargeM = new TH1F("hChargeM","CDC Mass",800,-4.,4.);
  TH1F* BarrMass = new TH1F("BarrMass","Barrel Mass", 800,-4.,4.);
 
  TH1F* ChiSquared = new TH1F("ChiSquared","Chi Squared", 1000,0.,100.);
  TH1F* ChiSquaredP = new TH1F("ChiSquaredP","Chi Squared Protons", 1000,0.,100.);
  dispTrT->Draw("chi>>ChiSquared",ValidFit&&Target);
  dispTrT->Draw("chi>>ChiSquaredP",ValidFit&&Target&&GoodVel&&Proton);

  TCanvas * ProfCanvas = new TCanvas("ProfCanvas", "Momentum Dependence",
				     200,50,700,900);
  ProfCanvas->Divide(1,3);
  ProfCanvas->cd(1);
  dispTrT->Draw("mass>>hMass");
  hMass->GetXaxis()->SetTitle("CDC Mass [GeV]");
  hMass->GetYaxis()->SetTitle("Counts");
  hChargeM->GetXaxis()->SetTitle("CDC Mass*Charge [GeV]");
  hChargeM->GetYaxis()->SetTitle("Counts");
  dispTrT->Draw("mass*charge>>hChargeM");
  BarrMass->GetXaxis()->SetTitle("Barrel Mass*Charge [GeV]");
  BarrMass->GetYaxis()->SetTitle("Counts");
  dispTrT->Draw("bmass*charge>>BarrMass","!bmass==0");
  dispTrT->Draw("(mom-momb)/momb:momb>>ProfRPC",ValidFit&&GoodVel&&Target&&ValidChi);
  ProfRPC->GetXaxis()->SetTitle("CDC Momentum [GeV]");
  ProfRPC->GetYaxis()->SetTitle("Momentum Relavite Diff.");
  //
  dispTrT->Draw("(mom-momb)/momb:momb>>hProf",ValidFit&&Target&&ValidChi);
  hProf->GetXaxis()->SetTitle("CDC Momentum [GeV]");
  hProf->GetYaxis()->SetTitle("Momentum Relavite Diff.");
  dispTrT->Draw("momb:rpcv>>MomBvsVel",GoodVel&&ValidFit&&Target&&ValidChi);
  MomBvsVel->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  MomBvsVel->GetYaxis()->SetTitle("CDC Momentum [GeV]");
  MomBvsVel->SetContour(20);
  /*  
  MomBvsVel->SetContourLevel(0,0);
  MomBvsVel->SetContourLevel(1,50);
  MomBvsVel->SetContourLevel(2,100);
  MomBvsVel->SetContourLevel(3,120);
  MomBvsVel->SetContourLevel(4,140);
  MomBvsVel->SetContourLevel(5,160);
  MomBvsVel->SetContourLevel(6,180);
  MomBvsVel->SetContourLevel(7,200);
  MomBvsVel->SetContourLevel(8,220);
  MomBvsVel->SetContourLevel(9,240);
  MomBvsVel->SetContourLevel(10,260);
  MomBvsVel->SetContourLevel(11,280);
  MomBvsVel->SetContourLevel(12,300);
  MomBvsVel->SetContourLevel(13,320);
  MomBvsVel->SetContourLevel(14,340);
  MomBvsVel->SetContourLevel(15,360);
  MomBvsVel->SetContourLevel(16,390);
  MomBvsVel->SetContourLevel(17,420);
  MomBvsVel->SetContourLevel(18,450);
  MomBvsVel->SetContourLevel(19,500);
  */
  hProf->Draw();
  //ProfRPC->Draw();
  //MomBvsVel->Draw("colz");
  //
  ProfCanvas->cd(2);
  dispTrT->Draw("(mom-momb)/momb:momb>>hdiffperc",ValidFit&&Target&&ValidChi);
  hdiffperc->GetXaxis()->SetTitle("CDC Momentum [GeV]");
  hdiffperc->GetYaxis()->SetTitle("Momentum Relative Diff.");
  //
  dispTrT->Draw("mom:rpcv>>MomVsVel",GoodVel&&ValidFit&&Target&&ValidChi);
  MomVsVel->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  MomVsVel->GetYaxis()->SetTitle("Combined Momentum [GeV]");
  dispTrT->Draw("momb*charge:rpcv>>ChMomBvsVel",GoodVel&&ValidFit&&Target&&ValidChi);
  dispTrT->Draw("mom*charge:rpcv>>ChMomVsVel",GoodVel&&ValidFit&&Target&&ValidChi);
  dispTrT->Draw("mom*charge:rpcv>>ChMomVsVel",GoodVel&&ValidFit&&Target&&ValidChi);
  dispTrT->Draw("momb*charge:rpcv>>ChMomBvsVelK",GoodVel&&ValidFit&&Target&&ValidChi&&KaonP);
  dispTrT->Draw("mom*charge:rpcv>>ChMomVsVelK",GoodVel&&ValidFit&&Target&&ValidChi&&KaonP);
  dispTrT->Draw("GfCdcM*charge:rpcv>>ChMomgfVsVelK",GoodVel&&ValidFit&&Target&&ValidChi&&KaonP);
  dispTrT->Draw("GfCdcM*charge:rpcv>>GfCorr",GoodVel&&ValidFit&&Target&&ValidChi);
  ChMomBvsVelK->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  ChMomVsVelK->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  ChMomBvsVelK->GetYaxis()->SetTitle("CDC Momentum [GeV/c]");
  GfCorr->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  GfCorr->GetYaxis()->SetTitle("GENFIT Momentum [GeV/c]");
  ChMomgfVsVelK->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  ChMomgfVsVelK->GetYaxis()->SetTitle("GENFIT Momentum [GeV/c]");
  //
  MomVsVel->SetContour(20);
  /*
  MomVsVel->SetContourLevel(0,0);
  MomVsVel->SetContourLevel(1,50);
  MomVsVel->SetContourLevel(2,100);
  MomVsVel->SetContourLevel(3,120);
  MomVsVel->SetContourLevel(4,140);
  MomVsVel->SetContourLevel(5,160);
  MomVsVel->SetContourLevel(6,180);
  MomVsVel->SetContourLevel(7,200);
  MomVsVel->SetContourLevel(8,220);
  MomVsVel->SetContourLevel(9,240);
  MomVsVel->SetContourLevel(10,260);
  MomVsVel->SetContourLevel(11,280);
  MomVsVel->SetContourLevel(12,300);
  MomVsVel->SetContourLevel(13,320);
  MomVsVel->SetContourLevel(14,340);
  MomVsVel->SetContourLevel(15,360);
  MomVsVel->SetContourLevel(16,390);
  MomVsVel->SetContourLevel(17,420);
  MomVsVel->SetContourLevel(18,450);
  MomVsVel->SetContourLevel(19,500);
  */
  hdiffperc->Draw("colz");
  //MomVsVel->Draw("colz");


  TCanvas * CorrelCanvas = new TCanvas("CorrelCanvas","Correlation Plots ",
				     200,20,800,900);

  CorrelCanvas->Divide(1,3);
  CorrelCanvas->cd(1);
  Int_t colors[] = {0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19};
  float levels[20] = {0,25,50,75,100,125,150,175,200,225,250,275,300,325,350,375,400,
		      425,450,480};
  //
  ChMomBvsVel->SetContour(20);
  /*
  ChMomBvsVel->SetContourLevel(0,0);
  ChMomBvsVel->SetContourLevel(1,50);
  ChMomBvsVel->SetContourLevel(2,100);
  ChMomBvsVel->SetContourLevel(3,120);
  ChMomBvsVel->SetContourLevel(4,140);
  ChMomBvsVel->SetContourLevel(5,160);
  ChMomBvsVel->SetContourLevel(6,180);
  ChMomBvsVel->SetContourLevel(7,200);
  ChMomBvsVel->SetContourLevel(8,220);
  ChMomBvsVel->SetContourLevel(9,240);
  ChMomBvsVel->SetContourLevel(10,260);
  ChMomBvsVel->SetContourLevel(11,280);
  ChMomBvsVel->SetContourLevel(12,300);
  ChMomBvsVel->SetContourLevel(13,320);
  ChMomBvsVel->SetContourLevel(14,340);
  ChMomBvsVel->SetContourLevel(15,360);
  ChMomBvsVel->SetContourLevel(16,390);
  ChMomBvsVel->SetContourLevel(17,420);
  ChMomBvsVel->SetContourLevel(18,450);
  ChMomBvsVel->SetContourLevel(19,500);
  */
  ChMomBvsVel->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  //ChMomBvsVel->GetXaxis()->SetTitleSize(0.55);
  //ChMomBvsVel->GetXaxis()->SetTitleOffset(1.0);
  ChMomBvsVel->GetYaxis()->SetTitle("CDC Momentum [GeV/c]");
  ChMomBvsVel->Draw("colz");
  CorrelCanvas->cd(3);
  ChMomVsVel->SetContour(20);
  /*
  ChMomVsVel->SetContourLevel(0,0);
  ChMomVsVel->SetContourLevel(1,50);
  ChMomVsVel->SetContourLevel(2,100);
  ChMomVsVel->SetContourLevel(3,120);
  ChMomVsVel->SetContourLevel(4,140);
  ChMomVsVel->SetContourLevel(5,160);
  ChMomVsVel->SetContourLevel(6,180);
  ChMomVsVel->SetContourLevel(7,200);
  ChMomVsVel->SetContourLevel(8,220);
  ChMomVsVel->SetContourLevel(9,240);
  ChMomVsVel->SetContourLevel(10,260);
  ChMomVsVel->SetContourLevel(11,280);
  ChMomVsVel->SetContourLevel(12,300);
  ChMomVsVel->SetContourLevel(13,320);
  ChMomVsVel->SetContourLevel(14,340);
  ChMomVsVel->SetContourLevel(15,360);
  ChMomVsVel->SetContourLevel(16,390);
  ChMomVsVel->SetContourLevel(17,420);
  ChMomVsVel->SetContourLevel(18,450);
  ChMomVsVel->SetContourLevel(19,500);
  */
  ChMomVsVel->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  //ChMomVsVel->GetXaxis()->SetTitleSize(0.55);
  //ChMomVsVel->GetXaxis()->SetTitleOffset(1.0);
  ChMomVsVel->GetYaxis()->SetTitle("Combined Momentum [GeV/c]");
  ChMomVsVel->Draw("colz");
  //
  gammabetaP->SetLineWidth(2);
  gammabetaP->SetLineStyle(2);
  gammabetaP->SetRange(5,26);
  //gammabetaP->DrawCopy("same");
  gammabetaPi->SetLineWidth(2);
  gammabetaPi->SetLineStyle(2);
  gammabetaPi->SetRange(5,30);
  gammabetaPiN->SetLineWidth(2);
  gammabetaPiN->SetLineStyle(2);
  gammabetaPiN->SetRange(5,30);
  gammabetaK->SetLineWidth(2);
  gammabetaK->SetLineStyle(2);
  gammabetaK->SetRange(5,27.5);
  gammabetaD->SetLineWidth(2);
  gammabetaD->SetLineStyle(2);
  gammabetaD->SetRange(5,20);
  gammabetaKN->SetLineWidth(2);
  gammabetaKN->SetLineStyle(2);
  gammabetaKN->SetRange(5,27.5);

  CorrelCanvas->cd(2);
  GfCorr->Draw("colz");

  TCanvas * ProjCanvas = new TCanvas("ProjCanvas", "Momentum Projection ",
				     200,10,800,900);
  ProjCanvas->Divide(1,3);
  ProjCanvas->cd(1);
  //
  dispTrT->Draw("(momb-(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPiPCdc",PionP&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(mom-(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPiPComb",PionP&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(GfCdcM-(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPiPGf",PionP&&ValidFit&&Target&&ValidChi&&GoodVel);
  //
  dispTrT->Draw("(momb-(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPiNCdc",PionN&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(mom-(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPiNComb",PionN&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(GfCdcM-(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.140/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPiNGf",PionN&&ValidFit&&Target&&ValidChi&&GoodVel);
  //
  dispTrT->Draw("(momb-(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPCdc",Proton&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(momb-(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>ProfProtCdc",Proton&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(mom-(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPComb",Proton&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(mom-(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>ProfProtComb",Proton&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(GfCdcM-(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoPGf",Proton&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(GfCdcM-(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.938/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>ProfProtGf",Proton&&ValidFit&&Target&&ValidChi&&GoodVel);
  // 
  //
   dispTrT->Draw("(momb-(0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoKCdc",KaonP&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(mom-(0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoKComb",KaonP&&ValidFit&&Target&&ValidChi&&GoodVel);
  dispTrT->Draw("(GfCdcM-(0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2)))/((0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))):(0.494/29.98)*rpcv/sqrt(1-pow(rpcv/29.98,2))>>histoKGf",KaonP&&ValidFit&&Target&&ValidChi&&GoodVel);
  //
  ProfProtCdc->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  ProfProtComb->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  ProfProtGf->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  ProfProtCdc->GetYaxis()->SetTitle("CDC-RPC Relative Momentum");
  ProfProtGf->GetYaxis()->SetTitle("GF-RPC Relative Momentum");
  ProfProtComb->GetYaxis()->SetTitle("combined-RPC Relative Momentum");
  ProfCanvas->cd(1);
  ProfProtCdc->Draw();
  ProfCanvas->cd(2);
  ProfProtGf->Draw();
  ProfCanvas->cd(3);
  ProfProtComb->Draw();
  //
  ProjCanvas->cd(1);
  //
  histoPiPComb->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPiPComb->GetYaxis()->SetTitle("Combined-RPC Relative Momentum");
  histoPiPGf->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPiPGf->GetYaxis()->SetTitle("GENFIT-RPC Relative Momentum");
  histoPiPCdc->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPiPCdc->GetYaxis()->SetTitle("CDC-RPC Relative Momentum");
  //
  histoPComb->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPComb->GetYaxis()->SetTitle("Combined-RPC Relative Momentum");
  histoPGf->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPGf->GetYaxis()->SetTitle("GENFIT-RPC Relative Momentum");
  histoPCdc->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPCdc->GetYaxis()->SetTitle("CDC-RPC Relative Momentum");
  //
  histoPiNComb->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPiNComb->GetYaxis()->SetTitle("Combined-RPC Relative Momentum");
  histoPiNGf->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPiNGf->GetYaxis()->SetTitle("GENFIT-RPC Relative Momentum");
  histoPiNCdc->GetXaxis()->SetTitle("RPC Momentum [GeV/c]");
  histoPiNCdc->GetYaxis()->SetTitle("CDC-RPC Relative Momentum");
  //
  /****************************************************************
  THE FIRST BIN OF THE PROJECTION IS BIN 1, NOT BIN 0 EVEN IF 0 
  DOES NOT RETURN ANY ERROR!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  ****************************************************************/
  
  // POSITIVE PIONS
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_150",3,3);
  histoPiPGf->ProjectionY("ProjPiPGf_150",3,3);
  histoPiPComb->ProjectionY("ProjPiPComb_150",3,3);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_200",4,4);
  histoPiPGf->ProjectionY("ProjPiPGf_200",4,4);
  histoPiPComb->ProjectionY("ProjPiPComb_200",4,4);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_250",5,5);
  histoPiPGf->ProjectionY("ProjPiPGf_250",5,5);
  histoPiPComb->ProjectionY("ProjPiPComb_250",5,5);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_300",6,6);
  histoPiPGf->ProjectionY("ProjPiPGf_300",6,6);
  histoPiPComb->ProjectionY("ProjPiPComb_300",6,6);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_350",7,7);
  histoPiPGf->ProjectionY("ProjPiPGf_350",7,7);
  histoPiPComb->ProjectionY("ProjPiPComb_350",7,7);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_400",8,8);
  histoPiPGf->ProjectionY("ProjPiPGf_400",8,8);
  histoPiPComb->ProjectionY("ProjPiPComb_400",8,8);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_450",9,9);
  histoPiPGf->ProjectionY("ProjPiPGf_450",9,9);
  histoPiPComb->ProjectionY("ProjPiPComb_450",9,9);
  //
  histoPiPCdc->ProjectionY("ProjPiPCdc_500",10,10);
  histoPiPGf->ProjectionY("ProjPiPGf_500",10,10);
  histoPiPComb->ProjectionY("ProjPiPComb_500",10,10);
  //
  // NEGATIVE PIONS
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_150",3,3);
  histoPiNGf->ProjectionY("ProjPiNGf_150",3,3);
  histoPiNComb->ProjectionY("ProjPiNComb_150",3,3);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_200",4,4);
  histoPiNGf->ProjectionY("ProjPiNGf_200",4,4);
  histoPiNComb->ProjectionY("ProjPiNComb_200",4,4);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_250",5,5);
  histoPiNGf->ProjectionY("ProjPiNGf_250",5,5);
  histoPiNComb->ProjectionY("ProjPiNComb_250",5,5);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_300",6,6);
  histoPiNGf->ProjectionY("ProjPiNGf_300",6,6);
  histoPiNComb->ProjectionY("ProjPiNComb_300",6,6);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_350",7,7);
  histoPiNGf->ProjectionY("ProjPiNGf_350",7,7);
  histoPiNComb->ProjectionY("ProjPiNComb_350",7,7);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_400",8,8);
  histoPiNGf->ProjectionY("ProjPiNGf_400",8,8);
  histoPiNComb->ProjectionY("ProjPiNComb_400",8,8);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_450",9,9);
  histoPiNGf->ProjectionY("ProjPiNGf_450",9,9);
  histoPiNComb->ProjectionY("ProjPiNComb_450",9,9);
  //
  histoPiNCdc->ProjectionY("ProjPiNCdc_500",10,10);
  histoPiNGf->ProjectionY("ProjPiNGf_500",10,10);
  histoPiNComb->ProjectionY("ProjPiNComb_500",10,10);
  //
  /* 
  // This is valid for 10MeV bin, starting at 50MeV
  histoPCdc->ProjectionY("ProjPCdc_300",26,26);
  histoPComb->ProjectionY("ProjPComb_300",26,26);
  histoPCdc->ProjectionY("ProjPCdc_350",31,31);
  histoPComb->ProjectionY("ProjPComb_350",31,31);
  histoPCdc->ProjectionY("ProjPCdc_400",36,36);
  histoPComb->ProjectionY("ProjPComb_400",36,36);
  histoPComb->ProjectionY("ProjPComb_440",40,40);
  histoPCdc->ProjectionY("ProjPCdc_450",41,41);
  histoPComb->ProjectionY("ProjPComb_450",41,41);
  histoPComb->ProjectionY("ProjPComb_460",42,42);
  histoPComb->ProjectionY("ProjPComb_470",43,43);
  histoPComb->ProjectionY("ProjPComb_480",44,44);
  histoPComb->ProjectionY("ProjPComb_490",45,45);
  histoPCdc->ProjectionY("ProjPCdc_500",46,46);
  histoPComb->ProjectionY("ProjPComb_500",46,46);
  histoPCdc->ProjectionY("ProjPCdc_550",51,51);
  histoPComb->ProjectionY("ProjPComb_550",51,51);
  histoPCdc->ProjectionY("ProjPCdc_600",56,56);
  histoPComb->ProjectionY("ProjPComb_600",56,56);
  histoPCdc->ProjectionY("ProjPCdc_650",61,61);
  histoPComb->ProjectionY("ProjPComb_650",61,61);
  histoPCdc->ProjectionY("ProjPCdc_700",66,66);
  histoPComb->ProjectionY("ProjPComb_700",66,66);
  histoPCdc->ProjectionY("ProjPCdc_800",86,86);
  histoPComb->ProjectionY("ProjPComb_800",86,86);
  histoPCdc->ProjectionY("ProjPCdc_900",96,96);
  histoPComb->ProjectionY("ProjPComb_900",96,96);
  histoPCdc->ProjectionY("ProjPCdc_1000",106,106);
  histoPComb->ProjectionY("ProjPComb_1000",106,106);
  */
  // PROTONS 
  //
  // bin at 50MeV, starting at 25MeV, so first bin is 50MeV
  histoPCdc->ProjectionY("ProjPCdc_50",1,1);
  histoPGf->ProjectionY("ProjPGf_50",1,1);
  histoPComb->ProjectionY("ProjPComb_50",1,1);
  //
  histoPCdc->ProjectionY("ProjPCdc_300",6,6);
  histoPGf->ProjectionY("ProjPGf_300",6,6);
  histoPComb->ProjectionY("ProjPComb_300",6,6);
  //
  histoPCdc->ProjectionY("ProjPCdc_350",7,7);
  histoPGf->ProjectionY("ProjPGf_350",7,7);
  histoPComb->ProjectionY("ProjPComb_350",7,7);
  //
  histoPCdc->ProjectionY("ProjPCdc_400",8,8);
  histoPGf->ProjectionY("ProjPGf_400",8,8);
  histoPComb->ProjectionY("ProjPComb_400",8,8);
  //
  histoPCdc->ProjectionY("ProjPCdc_450",9,9);
  histoPGf->ProjectionY("ProjPGf_450",9,9);
  histoPComb->ProjectionY("ProjPComb_450",9,9);
  //
  histoPCdc->ProjectionY("ProjPCdc_500",10,10);
  histoPGf->ProjectionY("ProjPGf_500",10,10);
  histoPComb->ProjectionY("ProjPComb_500",10,10);
  //
  histoPCdc->ProjectionY("ProjPCdc_550",11,11);
  histoPGf->ProjectionY("ProjPGf_550",11,11);
  histoPComb->ProjectionY("ProjPComb_550",11,11);
  //
  histoPCdc->ProjectionY("ProjPCdc_600",12,12);
  histoPGf->ProjectionY("ProjPGf_600",12,12);
  histoPComb->ProjectionY("ProjPComb_600",12,12);
  //
  histoPCdc->ProjectionY("ProjPCdc_700",14,14);
  histoPGf->ProjectionY("ProjPGf_700",14,14);
  histoPComb->ProjectionY("ProjPComb_700",14,14);
  //
  histoPCdc->ProjectionY("ProjPCdc_800",16,16);
  histoPGf->ProjectionY("ProjPGf_800",16,16);
  histoPComb->ProjectionY("ProjPComb_800",16,16);
  //
  histoPCdc->ProjectionY("ProjPCdc_1000",20,20);
  histoPGf->ProjectionY("ProjPGf_1000",20,20);
  histoPComb->ProjectionY("ProjPComb_1000",20,20);
  //
  histoPCdc->ProjectionY("ProjPCdc_1200",24,24);
  histoPGf->ProjectionY("ProjPGf_1200",24,24);
  histoPComb->ProjectionY("ProjPComb_1200",24,24);
  //
  // KAONS
  // 
  histoKCdc->ProjectionY("ProjKCdc_200",4,4);
  histoKCdc->ProjectionY("ProjKCdc_250",5,5);
  histoKCdc->ProjectionY("ProjKCdc_300",6,6);
  histoKCdc->ProjectionY("ProjKCdc_350",7,7);
  histoKCdc->ProjectionY("ProjKCdc_400",8,8);
  histoKCdc->ProjectionY("ProjKCdc_one",5,7);
  //histoKCdc->ProjectionY("ProjKCdc_450",9,9);
  //histoKCdc->ProjectionY("ProjKCdc_b",6,9);
  //
  histoKGf->ProjectionY("ProjKGf_200",4,4);
  histoKGf->ProjectionY("ProjKGf_250",5,5);
  histoKGf->ProjectionY("ProjKGf_300",6,6);
  histoKGf->ProjectionY("ProjKGf_350",7,7);
  histoKGf->ProjectionY("ProjKGf_400",8,8);
  histoKGf->ProjectionY("ProjKGf_one",5,7);
  //histoKGf->ProjectionY("ProjKGf_450",9,9);
  //
  histoKComb->ProjectionY("ProjKComb_200",4,4);
  histoKComb->ProjectionY("ProjKComb_250",5,5);
  histoKComb->ProjectionY("ProjKComb_300",6,6);
  histoKComb->ProjectionY("ProjKComb_350",7,7);
  histoKComb->ProjectionY("ProjKComb_400",8,8);
  histoKComb->ProjectionY("ProjKComb_one",5,7);

  //histoKComb->ProjectionY("ProjKComb_450",9,9);
  //
  /*
  ChMomVsVel->ProjectionY("projComb_10ns",50,50);
  ChMomBvsVel->ProjectionY("projCdc_10ns",50,50);
  ChMomVsVel->ProjectionY("projComb_11ns",55,55);
  ChMomBvsVel->ProjectionY("projCdc_11ns",55,55);
  ChMomVsVel->ProjectionY("projComb_12ns",60,60);
  ChMomBvsVel->ProjectionY("projCdc_12ns",60,60);
  ChMomVsVel->ProjectionY("projComb_13ns",65,65);
  ChMomBvsVel->ProjectionY("projCdc_13ns",65,65);
  ChMomVsVel->ProjectionY("projComb_14ns",70,70);
  ChMomBvsVel->ProjectionY("projCdc_14ns",70,70);
  ChMomVsVel->ProjectionY("projComb_15ns",75,75);
  ChMomBvsVel->ProjectionY("projCdc_15ns",75,75);
  ChMomVsVel->ProjectionY("projComb_16ns",80,80);
  ChMomBvsVel->ProjectionY("projCdc_16ns",80,80);
  ChMomVsVel->ProjectionY("projComb_17ns",85,85);
  ChMomBvsVel->ProjectionY("projCdc_17ns",85,85);
  ChMomVsVel->ProjectionY("projComb_18ns",90,90);
  ChMomBvsVel->ProjectionY("projCdc_18ns",90,90);
  ChMomVsVel->ProjectionY("projComb_19ns",95,95);
  ChMomBvsVel->ProjectionY("projCdc_19ns",95,95);
  ChMomVsVel->ProjectionY("projComb_20ns",100,100);
  ChMomBvsVel->ProjectionY("projCdc_20ns",100,100);
  ChMomVsVel->ProjectionY("projComb_21ns",105,105);
  ChMomBvsVel->ProjectionY("projCdc_21ns",105,105);
  ChMomVsVel->ProjectionY("projComb_22ns",110,110);
  ChMomBvsVel->ProjectionY("projCdc_22ns",110,110);
  ChMomVsVel->ProjectionY("projComb_23ns",115,115);
  ChMomBvsVel->ProjectionY("projCdc_23ns",115,115);
  ChMomVsVel->ProjectionY("projComb_24ns",120,120);
  ChMomBvsVel->ProjectionY("projCdc_24ns",120,120);
  ChMomVsVel->ProjectionY("projComb_25ns",125,125);
  ChMomBvsVel->ProjectionY("projCdc_25ns",125,125);
  ChMomVsVel->ProjectionY("projComb_26ns",130,130);
  ChMomBvsVel->ProjectionY("projCdc_26ns",130,130);
  ChMomVsVel->ProjectionY("projComb_27ns",135,135);
  ChMomBvsVel->ProjectionY("projCdc_27ns",135,135);
  ChMomVsVel->ProjectionY("projComb_274ns",137,137);
  ChMomBvsVel->ProjectionY("projCdc_274ns",137,137);
  ChMomVsVel->ProjectionY("projComb_28ns",140,140);
  ChMomBvsVel->ProjectionY("projCdc_28ns",140,140);
  ChMomVsVel->ProjectionY("projComb_29ns",145,145);
  ChMomBvsVel->ProjectionY("projCdc_29ns",145,145);
  ChMomVsVelK->ProjectionY("projKcomb_a",60,69);
  ChMomVsVelK->ProjectionY("projKcomb_b",70,79);
  ChMomVsVelK->ProjectionY("projKcomb_c",80,89);
  ChMomVsVelK->ProjectionY("projKcomb_d",90,99);
  */
  //
  
  ProjCanvas->cd(1);
  ProjPiNCdc_200->Draw();
  ProjCanvas->cd(2);
  ProjPiNGf_200->Draw();
  ProjCanvas->cd(3);
  ProjPiNComb_200->Draw();

  TCanvas * FreeCanvas = new TCanvas("FreeCanvas","Resolution Plots ",
				     200,10,600,700);
  FreeCanvas->Divide(1,2);
  FreeCanvas->cd(1);
  
  Float_t yCdcP[nP] = {20./305., -14./345., -14./395., -12./445., -8./495., -4./545., 
		       10./595., 13./645., 22./695., 29./745.};
  Float_t xP[nP] = {300, 350, 400, 450, 500, 550, 600, 650, 700, 800};
  Float_t eyCdcP[nP] = {25./305., 28./345., 32./395., 36./445., 38./495., 48./545., 
			49./595., 58./645., 68./695., 77./745.};
  Float_t exP[nP] = {5., 5., 5., 5., 5., 5., 5., 5., 5., 5.};
  Float_t yCdcPi[nPi] = {30./105., 31./115., 36./125., 38./135., 46./145., 50./155., 
			 53./165., 77./175.};
  Float_t xPi[nPi] = {105, 115, 125., 135, 145., 155., 165., 175.};

  Float_t eyCdcPi[nPi] = {9./105., 11./115., 16./125., 21./135., 29./145., 30./155., 
			  34./165., 71./175.};
  Float_t exPi[nPi] = {5., 5., 5., 5., 5., 5., 5., 5.};
  Float_t yGfP[nP] = {-16./305., -20./345., -21./395., -17./445., -14./495, -5./545., 
0./595., 6./645., 25./695., 42./745.};
  Float_t eyGfP[nP] = {16./305., 18./345., 21./395., 25./445., 32./495., 43./545., 
		       48./595., 57./645., 73./695.., 102./745.};
  Float_t yGfPi[nPi] = {23./105., 23./115., 31./125., 37./135., 41./145., 45./155., 
			55./165., 75./175.};
  Float_t eyGfPi[nPi] = {12./105., 13./115., 19./125., 23./135., 29./145., 37./155, 
			 33./165., 73./175.};

  TGraphErrors *grCdcP = new TGraphErrors(nP,xP, yCdcP,exP,eyCdcP);
  TGraphErrors *grCdcPi = new TGraphErrors(nPi,xPi,yCdcPi,exPi,eyCdcPi);
  TGraphErrors *grGfP = new TGraphErrors(nP,xP,yGfP,exP,eyGfP);
  TGraphErrors *grGfPi = new TGraphErrors(nPi,xPi,yGfPi,xPi,eyGfPi);
  grCdcP->SetTitle("Momentum Difference and Width");
  TH2F* hFrame = new TH2F("hFrame","Momentum Resolution",1000,50.,1050.,400,-2.,2);
  hFrame->GetXaxis()->SetTitle("RPC Momentum [MeV/c]");
  hFrame->GetYaxis()->SetTitle("Relative Momentum Difference");
  hFrame->Draw();
  grCdcP->Draw("same");
  grCdcPi->Draw("same");
  grGfP->Draw("same");
  grGfPi->Draw("same");
  
  TH2F * dEdxH = new TH2F("dEdxH","dE/dx in CDC",200,0.,2.,1500,0.,15000.);
  dispTrT->Draw("dEdx:momb>>dEdxH",Target&&ValidFit);
  dEdxH->GetXaxis()->SetTitle("CDC Momentum [GeV/c]");
  dEdxH->GetYaxis()->SetTitle("CDC dE/dx [a.u.]");
  //
  TH2F * dEdxTpc = new TH2F("dEdxTpc","dE/dx in TPC",200,0.,2.,2500,0.,10000.);
  dispTrT->Draw("dEdxTpc:momb>>dEdxTpc",Target&&ValidFit);
  dEdxTpc->GetXaxis()->SetTitle("CDC Momentum [GeV/c]");
  dEdxTpc->GetYaxis()->SetTitle("TPC dE/dx [a.u.]");
  //
  TH2F* dEdxKaonsTpc = new TH2F("dEdxKaonsTpc","dE/dx TPC, Kaons", 200,0.,2.,2000,0.,4000.);
  dispTrT->Draw("dEdxTpc:momb>>dEdxKaonsTpc",Target&&ValidFit&&KaonP);
  dEdxKaonsTpc->GetXaxis()->SetTitle("CDC Momentum [GeV/c]");
  dEdxKaonsTpc->GetYaxis()->SetTitle("TPC dE/dx [a.u.]");
  //
  TH2F* dEdxTpcSelection = new TH2F("dEdxTpcSelection","dE/dx TPC, good tracks", 200,0.,2.,2500,0.,10000.);
  //dispTrT->Draw("dEdxTpc:momb>>dEdxTpcSelection",Target&&ValidFit&&"ntpchit>13");
  //dispTrT->Draw("dEdxTpc:momb>>dEdxTpcSelection",Target&&ValidFit&&"chi>0&&chi<10&&ntpchit>12");
  dispTrT->Draw("dEdxTpc:momb>>dEdxTpcSelection",Target&&ValidFit&&"chi>0&&chi<5&&theta<1.");
  dEdxTpcSelection->GetXaxis()->SetTitle("CDC Momentum [GeV/c]");
  dEdxTpcSelection->GetYaxis()->SetTitle("TPC dE/dx [a.u.]");
  //
  histoPComb->Draw("colz");  
  //
  dEdxTpc->Draw("colz");
  FreeCanvas->cd(2);
  dEdxTpcSelection->Draw("colz");
  //
}

void GraphPlot()
{
  /*
  const Int_t  nP = 10; 
  const Int_t nPi = 5;
  Float_t xP[nP] = {300., 350., 400., 450., 500., 550., 600, 650, 700, 800};
  Float_t eX[nP] = {5., 5., 5., 5., 5., 5., 5., 5., 5.};
  Float_t yCdcP[nP] = {0.0753, 0.084, 0.087, 0.086, 0.094, 0.091, 0.099, 0.101, 0.101, 0.128};
  Float_t eyCdcP[nP] = {0.003, 0.002, 0.002, 0.002, 0.002, 0.002, 0.002, 0.002, 0.003, 0.004};
  //
  Float_t yGfP[nP] = {0.056, 0.0587, 0.0605, 0.062, 0.064, 0.067, 0.07, 0.072, 0.078, 0.099};
  Float_t eyGfP[nP] = {0.002, 0.0015, 0.0012, 0.002, 0.001, 0.001, 0.001, 0.001, 0.002, 0.003};
  //
  Float_t xPi[nP] = {150., 200., 250., 300., 350.};
  Float_t exPi[nP] = {5., 5., 5., 5., 5.};
  //
  Float_t yCdcPi[nP] = {0.140, 0.174, 0.211, 0.262, 0.268};
  Float_t eyCdcPi[nP] = {0.003, 0.003, 0.004, 0.005, 0.005};
  //
  Float_t yGfPi[nP] = {0.105, 0.157, 0.205, 0.261, 0.259};
  Float_t eyGfPi[nP] = {0.003, 0.004, 0.004, 0.005, 0.005};
  */
  const Int_t  nP = 6; 
  const Int_t nPiP = 8;
  const Int_t nPiN = 5;
  const Int_t nK = 4;
  Float_t xP[nP] = {330., 370., 410., 540., 700., 840.};
  Float_t eX[nP] = {5., 5., 5., 5., 5., 5.};
  Float_t yCdcP[nP] = {0.089, 0.077, 0.088, 0.081, 0.100, 0.113};
  Float_t eyCdcP[nP] = {0.003, 0.002, 0.002, 0.002, 0.002, 0.002};
  Float_t AveCdcP[nP] = {0.036, 0.014, 0.012, 0.0, 0.004, -0.015};
  //
  Float_t yGfP[nP] = {0.066, 0.061, 0.063, 0.068, 0.077, 0.087};
  Float_t eyGfP[nP] = {0.002, 0.0015, 0.0012, 0.002, 0.001, 0.001};
  Float_t AveGfP[nP] = {-0.037, -0.041, -0.044, -0.033, -0.031, -0.049};
  //
  Float_t xPiP[nPiP] = {150., 170., 190., 210., 245., 290., 365., 535. };
  Float_t exPiP[nPiP] = {5., 5., 5., 5., 5., 5., 5., 5.};
  //Float_t xPiP[nPiP] = {290., 365., 535. };
  //Float_t exPiP[nPiP] = {5., 5., 5.};
  Float_t xPiN[nPiN] = {150., 170., 190., 210., 245.};
  Float_t exPiN[nPiN] = {5., 5., 5., 5., 5.};
  Float_t xK[nK] = {238., 285., 340., 405.};
  Float_t exK[nK] = {10., 10., 10., 10.};
  Float_t exK[nK] = {0.01, 0.01, 0.01, 0.01};
  //
  Float_t yCdcPiP[nPiP] = {0.148, 0.131, 0.164, 0.162, 0.205, 0.137, 0.156, 0.124};
  Float_t eyCdcPiP[nPiP] = {0.003, 0.003, 0.004, 0.005, 0.005, 0.005, 0.005, 0.005};
  Float_t AveCdcPiP[nPiP] =  {0.140, 0.068, 0.095, 0.076, 0.094, 0.014, -0.016, -0.217};
  
  //
  Float_t yGfPiP[nPiP] = {0.119, 0.099, 0.130, 0.137, 0.156, 0.125, 0.147, 0.121};
  Float_t eyGfPiP[nPiP] = {0.003, 0.004, 0.004, 0.005, 0.005, 0.005, 0.005, 0.005};
  Float_t AveGfPiP[nPiP] =  {0.008, -0.039, -0.014, 0.004, 0.012, -0.028, -0.055, -0.243};
  //
  Float_t yCdcPiN[nPiN] = {0.126, 0.138, 0.156, 0.145, 0.168};
  Float_t eyCdcPiN[nPiN] = {0.003, 0.003, 0.004, 0.005, 0.005};
  Float_t AveCdcPiN[nPiN] =  {0.125, 0.085, 0.092, 0.093, 0.069};
  //
  Float_t yGfPiN[nPiN] = {0.098, 0.114, 0.133, 0.152, 0.158};
  Float_t eyGfPiN[nPiN] = {0.003, 0.004, 0.004, 0.005, 0.005};
  Float_t AveGfPiN[nPiN] =  {0.053, -0.001, 0.012, 0.044, 0.043};
  //
  Float_t yCdcK[nK] = {0.132, 0.140, 0.128, 0.152};
  Float_t eyCdcK[nK] = {0.003, 0.003, 0.004, 0.005};
  Float_t AveCdcK[nK] =  {0.077, 0.035, 0.026, -0.032};
  //
  Float_t yGfK[nK] = {0.093, 0.108, 0.103, 0.137};
  Float_t eyGfK[nK] = {0.003, 0.004, 0.004, 0.005};
  Float_t AveGfK[nK] =  {0.006, -0.000, -0.015, -0.066};

  TGraphErrors * gAveCdcP = new TGraphErrors(nP,xP,AveCdcP,eX,eyCdcP);
  TGraphErrors * gAveGfP = new TGraphErrors(nP,xP,AveGfP,eX,eyGfP);
  TGraphErrors * gYCdcP = new TGraphErrors(nP,xP,yCdcP,eX,eyCdcP);
  TGraphErrors * gYGfP = new TGraphErrors(nP,xP,yGfP,eX,eyCdcP);
  TGraphErrors * gAveCdcPiP = new TGraphErrors(nPiN,xPiN,AveCdcPiP,exPiN,eyCdcPiP);
  TGraphErrors * gAveGfPiP = new TGraphErrors(nPiN,xPiN,AveGfPiP,exPiN,eyGfPiP);
  TGraphErrors * gYCdcPiP = new TGraphErrors(nPiP,xPiP,yCdcPiP,exPiN,eyCdcPiP);
  TGraphErrors * gYGfPiP = new TGraphErrors(nPiP,xPiP,yGfPiP,exPiN,eyGfPiP);
  TGraphErrors * gAveCdcPiN = new TGraphErrors(nPiN,xPiN,AveCdcPiN,exPiN,eyCdcPiN);
  TGraphErrors * gAveGfPiN = new TGraphErrors(nPiN,xPiN,AveGfPiN,exPiN,eyGfPiN);
  TGraphErrors * gYCdcPiN = new TGraphErrors(nPiN,xPiN,yCdcPiN,exPiN,eyCdcPiN);
  TGraphErrors * gYGfPiN = new TGraphErrors(nPiN,xPiN,yGfPiN,exPiN,eyGfPiN);
  TGraphErrors * gAveCdcK = new TGraphErrors(nK,xK,AveCdcK,exK,eyCdcK);
  TGraphErrors * gAveGfK = new TGraphErrors(nK,xK,AveGfK,exK,eyGfK);
  TGraphErrors * gYCdcK = new TGraphErrors(nK,xK,yCdcK,exK,eyGfK);
  TGraphErrors * gYGfK = new TGraphErrors(nK,xK,yGfK,exK,eyGfK);
  /*
  TGraphErrors *gCdcProtW = new TGraphErrors(nP,xP,yCdcP,eX,eyCdcP);
  TGraphErrors *gGfProtW = new TGraphErrors(nP,xP,yGfP,eX,eyCdcP);
  TGraphErrors *gCdcPiW = new TGraphErrors(nPi,xPi,yCdcPi,exPi,eyCdcPi);
  TGraphErrors *gGfPiW = new TGraphErrors(nPi,xPi,yGfPi,exPi,eyGfPi);
  */
  TH2F* hResFramePr = new TH2F("hResFramePr","Momentum Resolution, Proton",
			       1000,50.,1050.,100,-0.5,0.5);
  TH2F* hResFramePi = new TH2F("hResFramePi","Momentum Resolution, Pion+",
			       350,50.,550.,100,-0.5,0.5);
  TH2F* hResFramePiN = new TH2F("hResFramePiN","Momentum Resolution, Pion-",
			       300,50.,500.,100,-0.5,0.5);
  TH2F* hResFrameK = new TH2F("hResFrameK","Momentum Resolution, K+",
			       300,50.,500.,100,-0.5,0.5);
  hResFramePr->GetXaxis()->SetTitle("Momentum [MeV/c]");
  //hResFramePr->GetYaxis()->SetTitle("Relative Momentum Resolution");
  hResFramePr->GetYaxis()->SetTitle("Relative Momentum Offset");
  hResFramePi->GetXaxis()->SetTitle("Momentum [MeV/c]");
  //hResFramePi->GetYaxis()->SetTitle("Relative Momentum Resolution");
  hResFramePi->GetYaxis()->SetTitle("Relative Momentum Offset");
  hResFrameK->GetXaxis()->SetTitle("Momentum [MeV/c]");
  //hResFrameK->GetYaxis()->SetTitle("Relative Momentum Resolution");
  hResFrameK->GetYaxis()->SetTitle("Relative Momentum Offset");
  hResFramePiN->GetXaxis()->SetTitle("Momentum [MeV/c]");
  hResFramePiN->GetYaxis()->SetTitle("Relative Momentum Offset");
  //hResFramePiN->GetYaxis()->SetTitle("Relative Momentum Resolution");

  TCanvas * MomResolCanvas = new TCanvas("MomResolCanvas", "Momentum Resolution ",
				     200,10,600,800);
  MomResolCanvas->Divide(1,2);
  MomResolCanvas->cd(1);
  hResFramePr->Draw();
  //gYGfP->Draw("");
  //gYCdcP->Draw("same");
  gAveGfP->Draw("");
  gAveCdcP->Draw("same");
  //gGfProtW->Draw();
  //gCdcProtW->Draw("same");
  MomResolCanvas->cd(2);
  hResFramePi->Draw();
  //gYGfPiP->Draw("");
  //gYCdcPiP->Draw("same");
  gAveGfPiP->Draw("");
  gAveCdcPiP->Draw("same");
  //gGfPiW->Draw();
  //gCdcPiW->Draw("same");
  
  TCanvas * MomResolCanvasK = new TCanvas("MomResolCanvasK", "Momentum Resolution ",
				     200,10,600,800);
  MomResolCanvasK->Divide(1,2);
  MomResolCanvasK->cd(1);
  hResFrameK->Draw();
  //gYGfK->Draw("");
  //gYCdcK->Draw("same");
  gAveGfK->Draw("");
  gAveCdcK->Draw("same");
  //gGfProtW->Draw();
  //gCdcProtW->Draw("same");
  MomResolCanvasK->cd(2);
  hResFramePiN->Draw();
  //gYGfPiN->Draw("");
  //gYCdcPiN->Draw("same");
  gAveGfPiN->Draw("");
  gAveCdcPiN->Draw("same");
  //gGfPiW->Draw();
  //gCdcPiW->Draw("same");

}

void MartinPlot()
{
 const Int_t  nP = 5; 
 Float_t xP[nP] = {300., 400., 500., 600., 700.};
 Float_t exP[nP] = {25., 25., 25., 25., 25.};
 Float_t yCdcP[nP] = {0.074, 0.08, 0.092, 0.093, 0.098};
 Float_t eyCdcP[nP] = {0.003, 0.003, 0.002, 0.002, 0.003};
 Float_t yGfP[nP] = {0.054, 0.060, 0.066, 0.070, 0.074};
 Float_t eyGfP[nP] = {0.002, 0.002, 0.002, 0.002, 0.002};

 TGraphErrors * gYCdcP = new TGraphErrors(nP,xP,yCdcP,exP,eyCdcP);
 TGraphErrors * gYGfP = new TGraphErrors(nP,xP,yGfP,exP,eyGfP);

  TH2F* hResFramePr = new TH2F("hResFramePr","Momentum Resolution, Proton",
			       1000,50.,1050.,100,-0.5,0.5);

  hResFramePr->GetXaxis()->SetTitle("Momentum [MeV/c]");
  hResFramePr->GetYaxis()->SetTitle("Relative Momentum Resolution");

  TCanvas * MomResolCanvas = new TCanvas("MomResolCanvas", "Momentum Resolution ",
				     200,10,600,800);

  hResFramePr->Draw();
  gYCdcP->Draw("P");
  gYGfP->Draw("Psame");  

}

void TriplePlots()
{
  const Int_t  nP = 7; 
  Float_t xP[nP] = {300., 400., 450., 500., 550., 600., 700.};
  Float_t exP[nP] = {25., 25., 25., 25., 25., 25., 25.};
  //
  Float_t yCdcP[nP] = {0.067, 0.071, 0.070, 0.071, 0.073, 0.076, 0.077}; 
  Float_t eyCdcP[nP] = {0.005, 0.003, 0.002, 0.002, 0.003, 0.003, 0.003};
  Float_t AveCdcP[nP] = {0.028, 0.017, 0.015, 0.009, 0.004, 0.001, -0.001}; 
  Float_t eAveCdcP[nP] = {0.002, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001};
  //
  Float_t yGfP[nP] = {0.057, 0.063, 0.063, 0.066, 0.067, 0.068, 0.073};
  Float_t eyGfP[nP] = {0.004, 0.002, 0.002, 0.002, 0.002, 0.002, 0.002};
  Float_t AveGfP[nP] = {-0.023, -0.026, -0.024, -0.026, -0.029, -0.029, -0.031}; 
  Float_t eAveGfP[nP] = {0.002, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001};
  //
  Float_t ycombP[nP] = {0.053, 0.055, 0.055, 0.056, 0.058, 0.060, 0.064};
  Float_t eycombP[nP] = {0.005, 0.003, 0.002, 0.002, 0.003, 0.003, 0.003};
  Float_t AvecombP[nP] = {-0.045, -0.043, -0.039, -0.039, -0.038, -0.037, -0.038}; 
  Float_t eAvecombP[nP] = {0.002, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001};
  //
  
  TGraphErrors * gYCdcP = new TGraphErrors(nP,xP,yCdcP,exP,eyCdcP);
  TGraphErrors * gYGfP = new TGraphErrors(nP,xP,yGfP,exP,eyGfP);
  TGraphErrors * gYcombP = new TGraphErrors(nP,xP,ycombP,exP,eycombP);
  //
  TGraphErrors * gAveCdcP = new TGraphErrors(nP,xP,AveCdcP,exP,eAveCdcP);
  TGraphErrors * gAveGfP = new TGraphErrors(nP,xP,AveGfP,exP,eAveGfP);
  TGraphErrors * gAvecombP = new TGraphErrors(nP,xP,AvecombP,exP,eAvecombP);
  
  
  TH2F* hResFramePr = new TH2F("hResFramePr","Momentum Resolution, Protons",
			       1000,50.,1050.,100,-0.5,0.5);
  TH2F* hOffFramePr = new TH2F("hOffFramePr","Momentum Offset, Protons",
			       1000,50.,1050.,100,-0.5,0.5);
  
  hResFramePr->GetXaxis()->SetTitle("Momentum [MeV/c]");
  hResFramePr->GetYaxis()->SetTitle("Relative Momentum Resolution");
  hOffFramePr->GetXaxis()->SetTitle("Momentum [MeV/c]");
  hOffFramePr->GetYaxis()->SetTitle("Relative Momentum Offset");

  TCanvas * MomResolCanvas = new TCanvas("MomResolCanvas", "Momentum Resolution",
				     200,10,600,800);
  
  hResFramePr->Draw();
  gYCdcP->Draw("P");
  gYGfP->Draw("Psame");
  gYcombP->Draw("Psame");

  TCanvas * MomOffCanvas = new TCanvas("MomOffCanvas", "Momentum Offset",
				       200,10,600,800);
  
  hOffFramePr->Draw();
  gAveCdcP->Draw("P");
  gAveGfP->Draw("Psame");
  gAvecombP->Draw("Psame");

  const Int_t kaonP = 3;
  Float_t xKaon[nP] = {225., 300., 350.};
  Float_t exK[nP] = {50., 25., 25.};
  //
  Float_t yCdcK[nP] = {0.089, 0.086, 0.089}; 
  Float_t eyCdcK[nP] = {0.011, 0.009, 0.009};
  Float_t AveCdcK[nP] = {0.044, 0.027, 0.013}; 
  Float_t eAveCdcK[nP] = {0.012, 0.011, 0.010};
  //
  Float_t yGfK[nP] = {0.10, 0.077, 0.081};
  Float_t eyGfK[nP] = {0.012, 0.011, 0.010};
  Float_t AveGfK[nP] = {0.008, 0.0, -0.033}; 
  Float_t eAveGfK[nP] = {0.012, 0.011, 0.010};
  //
  Float_t ycombK[nP] = {0.090, 0.063, 0.075};
  Float_t eycombK[nP] = {0.012, 0.011, 0.010};
  Float_t AvecombK[nP] = {-0.016, -0.012, -0.038}; 
  Float_t eAvecombK[nP] = {0.012, 0.011, 0.010};
  //
  
  TGraphErrors * gYCdcK = new TGraphErrors(kaonP,xKaon,yCdcK,exK,eyCdcK);
  TGraphErrors * gYGfK = new TGraphErrors(kaonP,xKaon,yGfK,exK,eyGfK);
  TGraphErrors * gYcombK = new TGraphErrors(kaonP,xKaon,ycombK,exK,eycombK);
  //
  TGraphErrors * gAveCdcK = new TGraphErrors(kaonP,xKaon,AveCdcK,exK,eAveCdcK);
  TGraphErrors * gAveGfK = new TGraphErrors(kaonP,xKaon,AveGfK,exK,eAveGfK);
  TGraphErrors * gAvecombK = new TGraphErrors(kaonP,xKaon,AvecombK,exK,eAvecombK);
  
  
  TH2F* hResFrameK = new TH2F("hResFrameK","Momentum Resolution, Kaons",
			       1000,50.,1050.,100,-0.5,0.5);
  TH2F* hOffFrameK = new TH2F("hOffFrameK","Momentum Offset, Kaons",
			       1000,50.,1050.,100,-0.5,0.5);
  
  hResFrameK->GetXaxis()->SetTitle("Momentum [MeV/c]");
  hResFrameK->GetYaxis()->SetTitle("Relative Momentum Resolution");
  hOffFrameK->GetXaxis()->SetTitle("Momentum [MeV/c]");
  hOffFrameK->GetYaxis()->SetTitle("Relative Momentum Offset");

  //TCanvas * MomResolCanvas = new TCanvas("MomResolCanvas", "Momentum Resolution",
  //				     200,10,600,800);
  MomResolCanvas->cd();

  hResFrameK->Draw();
  gYCdcK->Draw("P");
  gYGfK->Draw("Psame");
  gYcombK->Draw("Psame");

  //  TCanvas * MomOffCanvas = new TCanvas("MomOffCanvas", "Momentum Offset",
  //				       200,10,600,800);
  MomOffCanvas->cd();
  
  hOffFrameK->Draw();
  gAveCdcK->Draw("P");
  gAveGfK->Draw("Psame");
  gAvecombK->Draw("Psame");

}

void CdcPlot()
{

  TCut GoodVel = "RpcTrackBr.fVel>0.1";
  TCut GoodMom = "CdcTrackBr.fMom>0.05";

  TF1 * gammabetaP = new TF1("gammabetaP","(0.938/29.98)*x/sqrt(1-pow(x/29.98,2))",
			      1.,29.9);
 TF1 * gammabetaPi = new TF1("gammabetaPi","(0.14/29.98)*x/sqrt(1-pow(x/29.98,2))",
                              1.,29.9);
 TF1 * gammabetaPiN = new TF1("gammabetaPiN","-(0.14/29.98)*x/sqrt(1-pow(x/29.98,2))",
			      1.,29.9);
 
 TF1 * gammabetaK = new TF1("gammabetaK","(0.494/29.98)*x/sqrt(1-pow(x/29.98,2))",
			    1.,29.9);
 TF1 * gammabetaKN = new TF1("gammabetaKN","-(0.494/29.98)*x/sqrt(1-pow(x/29.98,2))",
			     1.,29.9);
 TF1 * gammabetaD = new TF1("gammabetaD","(1.876/29.98)*x/sqrt(1-pow(x/29.98,2))",
			    1.,29.);                                

  TH2F* MomVsVel = new TH2F("MomVsVel","CDC Momentum vs RPC Velocity",175,0.,35.,75.,0.05,1.55);
  TH2F* ChMomVsVel = new TH2F("ChMomVsVel","CDC Momentum vs RPC Velocity",175,0.,35.,200.,-2.,2.);
  FOPITree->Draw("CdcTrackBr.fMom*CdcTrackBr.fCharge:RpcTrackBr.fVel>>ChMomVsVel",GoodVel&&GoodMom);
  ChMomVsVel->GetXaxis()->SetTitle("RPC Velocity [cm/ns]");
  ChMomVsVel->GetYaxis()->SetTitle("CDC Momentum [GeV]");

}