C++ TH1D::DrawCopy方法代码示例

void MCtoSTATerrors() {
  TLatex *tplus = labelLatex(0.1364943,0.8114407,"MC sample");

  TFile *file0 = TFile::Open("results/StandardCutflow/allmc.root");

  TCanvas *c1 = new TCanvas();
  //TCanvas *c2 = new TCanvas();
  TH1D * LPplus = (TH1D*)file0->Get("RECO_PolPlots_50toinf/RECO_ICVarPFPlus");
  TH1D * LPminus = (TH1D*)file0->Get("RECO_PolPlots_50toinf/RECO_ICVarPFMinus");

  LPplus->Rebin(20);
  LPminus->Rebin(20);

  LPplus->Scale(0.3);
  LPminus->Scale(0.3);

  TH1D * newLPplus = (TH1D*)LPplus->Clone();
  TH1D * newLPminus = (TH1D*)LPminus->Clone();

  for(unsigned int i=1; i < LPplus->GetXaxis()->GetNbins(); i++) {
    newLPplus->SetBinError(i, TMath::Sqrt(LPplus->GetBinContent(i)));
    newLPminus->SetBinError(i, TMath::Sqrt(LPminus->GetBinContent(i)));
  }

  c1->cd();
  newLPplus->GetYaxis()->SetTitle("Events / 300 nb^{-1}");
  //newLPplus->GetXaxis()->SetTitle("LP(#mu^{+})");
  newLPplus->GetXaxis()->SetTitle("LP(#mu)");
  newLPplus->GetXaxis()->SetRangeUser(-1.0,2.0);
  newLPplus->SetLineWidth(3);
  newLPplus->SetLineStyle(2);
  newLPplus->SetLineColor(kRed);
  newLPplus->DrawCopy();
  tplus->DrawClone("same");

//   c2->cd();
//   newLPminus->GetYaxis()->SetTitle("Events / 100 nb^{-1}");
//   newLPminus->GetXaxis()->SetTitle("LP(#mu^{-})");
  newLPminus->SetLineWidth(3);
  newLPminus->SetLineStyle(2);
  newLPminus->SetLineColor(kBlack);
  newLPminus->DrawCopy("same");
//   tplus->DrawClone("same");

  leg = new TLegend(0.762931,0.720339,0.8649425,0.8622881,NULL,"brNDC");
  leg->AddEntry(newLPplus,"#mu^{+}","l");
  leg->AddEntry(newLPminus,"#mu^{-}","l");

  leg->SetFillColor(kWhite);
  leg->SetBorderSize(0);
  leg->SetTextFont(62);
  leg->DrawClone();

  file0->Close();

  return;
}

void zeeAnaDiff(TH1D* hsee, TCut addCut,TH1D* hseeMc) {
  
  TFile *f1=new TFile(fname1.Data());
  TTree *zp = (TTree*)f1->Get("tz");
  zp->AddFriend("yEvt=yongsunHiEvt"       ,fname1.Data());
  zp->AddFriend("ySkim=yongsunSkimTree"   ,fname1.Data());
  zp->AddFriend("yHlt=yongsunHltTree"     ,fname1.Data());
  zp->AddFriend("tgj",                     fname1.Data());
  
  TH1D* hmass = new TH1D("hmass",";inv. mass (GeV);event fraction",34,30,200);
  hmass->Sumw2();

  zp->Draw("invm>>hmass",addCut && "isEle1==1 && isEle2==1");
  TH1D* htemp = (TH1D*)hsee->Clone("htemp");   htemp->Reset();
  zp->Draw("see1>>htemp",addCut && "isEle1==1 && isEle2==1" && "invm>80 && invm<100" );
  hsee->Reset();
  hsee->Add(htemp);
  handsomeTH1(hmass,2);
  scaleInt(hmass);
  
  TH1D* htempMC;
  TFile *f1mc;
  TTree *zpmc ;
  f1mc=new TFile(fname1MC.Data());
  zpmc = (TTree*)f1mc->Get("tz");
  zpmc->AddFriend("yEvt=yongsunHiEvt"       ,fname1MC.Data());
  zpmc->AddFriend("ySkim=yongsunSkimTree"   ,fname1MC.Data());
  zpmc->AddFriend("yHlt=yongsunHltTree"     ,fname1MC.Data());
  zpmc->AddFriend("yPho=yongsunPhotonTree"     ,fname1MC.Data());
  zpmc->AddFriend("tgj",                     fname1MC.Data());
  hseeMc->Sumw2();
  htempMC = (TH1D*)hseeMc->Clone("htempMc");   htempMC->Reset();
  TH1D* hmassMC = new TH1D("hmassMC11",";inv. mass (GeV);event fraction",34,30,200);
  hmassMC->Sumw2();
  zpmc->Draw("invm>>hmassMC11",Form("(%s) *yPho.ncoll",addCut.GetTitle()));
  zpmc->Draw("see1>>htempMc",Form("(%s) *yPho.ncoll", (addCut && "invm>80 && invm<100").GetTitle()) );
  hseeMc->Reset();
  hseeMc->Add(htempMC);
  handsomeTH1(hmassMC,1);
  scaleInt(hmassMC);
  
  
  hmassMC->SetAxisRange(0,0.6,"Y");
  hmassMC->DrawCopy();
  hmass->DrawCopy("same");                                                                                             
  TLegend* leg1 =  new TLegend(0.5580645,0.5381356,0.8923387,0.7648305,NULL,"brNDC");
  easyLeg(leg1,"Z->ee mass");
  leg1->AddEntry(hmass,"DATA","pl");
  leg1->AddEntry(hmassMC,"MC","pl");
  leg1->Draw();
  
}

void mc2ibd_R2(TChain *tMC, TFile *fBgnd, TCanvas *cv)
{
	char str[1024];
	TLatex *txt = new TLatex();

	TH1D *hExp = (TH1D *) fBgnd->Get("hR2A-diff");
	if (!hExp) {
		printf("Histogram hR2A-diff not found in %s\n", fBgnd->GetName());
		return;
	}
	hExp->SetTitle("Distance between positron and neutron, 3D case;cm;Events/4cm");
	hExp->SetLineColor(kBlack);
	hExp->SetLineWidth(3);

	gROOT->cd();
	TH1D *hMC = new TH1D("hR2MC", "Distance between positron and neutron, 3D case (MC);cm;Events/4cm", 40, 0, 160);
	hMC->SetLineColor(kBlue);
	tMC->Project(hMC->GetName(), "Distance", cX && cY && cZ && cRXY && c20 && cGamma && cGammaMax && cPe && cN);
	hMC->Sumw2();
	hMC->Scale(hExp->Integral() / hMC->Integral());
	
	cv->Clear();
	hExp->DrawCopy();
	hMC->Draw("hist,same");
	TLegend *lg = new TLegend(0.65, 0.8, 0.89, 0.89);
	lg->AddEntry(hExp, "IBD", "LE");
	lg->AddEntry(hMC, "MC", "L");
	lg->Draw();
}

void PlotWriteErrors(TFile* input, TFile* output, TString output_folder, TString canvas, TString Hist, TString Opt="")
{

  // Setup the canvas
  TCanvas *c1= new TCanvas(canvas,canvas,800,700);

  // Get the histograms from the files
  TH1D *Data = (TH1D*)input->Get(Hist);

  //check to make sure there are some events before setting the log scale
  if(Data->Integral() == 0 ) c1->SetLogy(0);
  else  c1->SetLogy(1);

  // Fill for histogram
  Data->SetFillColor(kBlue);
  Data->GetXaxis()->SetTitleSize(0.06);
  Data->GetXaxis()->SetTitleOffset(0.75);
  Data->GetYaxis()->SetTitleSize(0.05);
  Data->GetYaxis()->SetTitleOffset(1.00);

  // plot them
  Data->DrawCopy("hist");

  gPad->RedrawAxis();

  //write canvas as png
  c1->Print(TString(output_folder+canvas+".png"));  

  //write canvas to output file
  output->cd();
  c1->Write();
  return;
}

void mc2ibd_NE(TChain *tMC, TFile *fBgnd, TCanvas *cv)
{
	char str[1024];
	TLatex *txt = new TLatex();

	TH1D *hExp = (TH1D *) fBgnd->Get("hNEA-diff");
	if (!hExp) {
		printf("Histogram hNEA-diff not found in %s\n", fBgnd->GetName());
		return;
	}
	hExp->SetTitle("Delayed event energy;MeV;Events/200 keV");
	hExp->SetLineColor(kBlack);
	hExp->SetLineWidth(3);

	gROOT->cd();
	TH1D *hMC = new TH1D("hNEMC", "Delayed event energy (MC);MeV;Events/200 keV", 45, 3, 12);
	hMC->SetLineColor(kBlue);
	tMC->Project(hMC->GetName(), "NeutronEnergy", cX && cY && cZ && cR && c20 && cGamma && cGammaMax && cPe);
	hMC->Sumw2();
	hMC->Scale(hExp->Integral(15, 45) / hMC->Integral(15,45));
	
	cv->Clear();
	hMC->Draw("hist");
	hExp->DrawCopy("same");
	TLegend *lg = new TLegend(0.65, 0.8, 0.89, 0.89);
	lg->AddEntry(hExp, "IBD", "LE");
	lg->AddEntry(hMC, "MC", "L");
	lg->Draw();
}

void DrawHisto(TH1* h_orig, TString name,  Option_t *drawopt){
	TH1D* h = (TH1D*)h_orig->Clone(name);
	TString canvname = "canv_"+name;
	TCanvas *col = new TCanvas(canvname, "", 0, 0, 500, 500);
	col->SetFillStyle(0);
	col->SetFrameFillStyle(0);
	gPad->SetFillStyle(0);
	gStyle->SetPalette(1);
	gPad->SetRightMargin(0.15);
	h->DrawCopy(drawopt);
}

TCanvas* pMenu(std::string what)
{
  TCanvas* c = new TCanvas(("cMenu"+what).c_str(),("cMenu"+what).c_str(),1500,400); 
  c->SetBottomMargin(0.6);
  c->SetLeftMargin(0.04);
  c->SetRightMargin(0.01);
  TH1D* h = (TH1D*)gROOT->FindObject(("hMenuAlgos"+what).c_str());
  h->SetLineColor(4);
  h->SetFillColor(4);
  h->SetFillStyle(3003);
  h->LabelsOption("v");
  h->GetYaxis()->SetLabelOffset(0.005);
//  h->GetXaxis()->SetRange(550,680);
  h->DrawCopy();
  return c;
}

void PlotErrors(TFile* data, TString Hist, TString Opt="")
{

    // Get the histograms from the files
    TH1D *Data = (TH1D*)data->Get(Hist);

    //check to make sure there are some events for log scale
    if(Data->Integral() == 0 ) c1->SetLogy(0);
    else  c1->SetLogy(1);

    // Fill for histogram
    Data->SetFillColor(kBlue);

    // plot them
    Data->DrawCopy("hist");

    gPad->RedrawAxis();

}

void photonEnergyValidate   (sampleType collision = kHIMC) {
  TH1::SetDefaultSumw2(); 
  
  TFile* f = new TFile("photonEnergyScaleTable_lowPt_v6.root");
  TF1* curve[2][4];  // [ conversion 0 or 1 ] [ centrality bin 1,2,3] 
  for ( int icent = 1 ; icent<=3 ; icent++) {
    curve[0][icent] = (TF1*)f->Get(Form("fit_hscale_r9gt94_%d",icent));   // Not converted
    curve[1][icent] = (TF1*)f->Get(Form("fit_hscale_r9lt94_%d",icent));  // Converted
    curve[1][icent]->SetLineStyle(7);
  }
  TH1D* htemp = new TH1D("htemp",";p_{T}^{#gamma} (GeV); Residual correction",200,35,200);
  
  TCanvas* c1 = new TCanvas("c1","",500,500);

  handsomeTH1(htemp,0);
  htemp->SetAxisRange(0.97,1.12,"Y");
  htemp->DrawCopy();
  
  for (int icent = 1 ;icent<=3 ; icent++) {
    curve[0][icent]->Draw("same");
    curve[1][icent]->Draw("same");
    //    onSun(35,1,200,1);
  }
  TLegend *l0 =  new TLegend(0.4910859,0.5132576,1,0.7102273,NULL,"brNDC");
  easyLeg(l0,"r9 > 0.94 : Not Conv");
  l0->AddEntry(curve[0][1],"0-10%","l");
  l0->AddEntry(curve[0][2],"10-30%","l");
  l0->AddEntry(curve[0][3],"30-100%","l");
  l0->Draw();
  
  TLegend *l1 = new TLegend(0.4910859,0.7443182,1,0.9318182,NULL,"brNDC");
   easyLeg(l1,"r9 < 0.94 : Conv");
  l1->AddEntry(curve[1][1],"0-10%","l");
  l1->AddEntry(curve[1][2],"10-30%","l");
  l1->AddEntry(curve[1][3],"30-100%","l");
  l1->Draw();
  

       

}

void drawJacob(string fnm="../../output/thinRadiator/o_msc_L001um_03MeV_1e5ev.root"){
  TFile *fin=TFile::Open(fnm.c_str(),"READ");
  TTree *t=(TTree*)fin->Get("t");
  int material, pType, trackID, parentID;
  double pX, pY, pZ;
  
  t->SetBranchAddress("pType", &pType);
  t->SetBranchAddress("material", &material);
  t->SetBranchAddress("trackID", &trackID);
  t->SetBranchAddress("parentID", &parentID);
  t->SetBranchAddress("postMomX", &pX);
  t->SetBranchAddress("postMomY", &pY);
  t->SetBranchAddress("postMomZ", &pZ);

  int nev = t->GetEntries();
  TH1D *scatAng = new TH1D("scatAng","Rate/sin(theta); scattering angle [deg]", 400,0,100);
  double pi = acos(-1);
  cout<<"pi "<<pi<<endl;
  for(int i=0;i<nev;i++){
    t->GetEntry(i);

    if(pType!=11) continue;
    if(material != 1) continue;
    if(trackID != 1 || parentID!=0) continue;

    double r = sqrt( pX*pX + pY*pY + pZ*pZ);
    if(r==0) continue;

    double thetaRad = acos(pZ/r);
    double theta = thetaRad * 180/pi;
    if(theta<0.01) continue;
    
    scatAng->Fill(theta, 1./sin(thetaRad));
    
  }
  gStyle->SetOptStat("eMRou");
  scatAng->DrawCopy();
  fin->Close();
}

void draw_graphs(void)
{
	const char title[4][20] = {"Na_mean;um", "Na_sigma;um", "Co_mean;um", "Co_sigma;um"};
	const double tgt[4] = {1.787, 0.394, 2.082, 0.473};
	const double mc[4][4] = {{2.008, 1.985, 1.974, 1.933}, {0.319, 0.326, 0.330, 0.336}, 
		{2.233, 2.211, 2.185, 2.149}, {0.377, 0.385, 0.394, 0.405}};
	const double dl[4] = {0, 50, 100, 200};
	const double mm[4][2] = {{1.7, 2.1}, {0.3, 0.4}, {2.0, 2.3}, {0.35, 0.5}};
	TGraph *gr;
	TLine *ln;
	TH1D *h;
	TLegend *lg;
	int i;
	
	TCanvas *cv = new TCanvas;
	cv->Divide(2,2);
	gr = new TGraph;
	ln = new TLine;
	gr->SetMarkerStyle(kFullCircle);
	gr->SetMarkerSize(1.5);
	gr->SetMarkerColor(kBlue);
	ln->SetLineWidth(4);
	ln->SetLineColor(kRed);
	lg = new TLegend(0.65, 0.65, 0.95, 0.85);
	lg->AddEntry(ln, "Exp.", "L");
	lg->AddEntry(gr, "MC", "P");
	for (i=0; i<4; i++) {
		cv->cd(i+1);
		
		h = new TH1D(title[i], title[i], 100, 0, 200);
		h->SetMinimum(mm[i][0]);
		h->SetMaximum(mm[i][1]);
		h->DrawCopy();
		gr->DrawGraph(4, dl, mc[i], "P");
		ln->DrawLine(0, tgt[i], 200, tgt[i]);
		lg->Draw();
		delete h;
	}
}

//*************************************************************************************************
//Main part of the macro
//*************************************************************************************************
void AnalyzeYYPhotonEvent(const string InputFilename) {

  
  TTree* YYTree = getTreeFromFile(InputFilename.c_str());
  assert(YYTree);    
  YYPhotonEvent photon(YYTree);
      
  TH1D *plot = new TH1D("hPlot", ";XAxis; Number of Events ", 6, -3,3);

  for (int n=0;n<YYTree->GetEntries();n++) { 
    photon.GetEntry(n);

    Double_t weight = photon.photon_branch_weight;    
    //cout << weight << " " << photon.photon_branch_NPhotons << endl;
    if (photon.photon_branch_NPhotons > 0 ) {
      plot->Fill(photon.photon_branch_Photon1Eta);
    }

  }

  TCanvas *c = new TCanvas("canvas" , "canvas" , 800, 600);
  plot->DrawCopy();
  
}

//.........这里部分代码省略.........
  }
  
  c3->cd(1);
  tempJ->Draw();
  handsomeTH1(meanJetPt[2],1 );
//  meanJetPt[2]->Draw();
  drawSys(meanJetPt[2], meanJetPtSys[2], kGreen);
  drawSys(meanJetPt[4], meanJetPtSys[4]);
  meanJetPt[2]->Draw("same");
  handsomeTH1(meanJetPt[4],2 );
  meanJetPt[4]->SetMarkerStyle(24);
  meanJetPt[4]->Draw("same");
  if ( 1==1 ){
    TLegend *ly = new  TLegend(0.4937249,0.262193,0.8493139,0.4740574,NULL,"brNDC");
    easyLeg(ly,"2.76TeV");
    ly->AddEntry(meanJetPt[4],"PbPb 30-100%","p");
    ly->AddEntry(meanJetPt[2],"pp (smeared)","p");
    ly->Draw();
  }
  c3->SaveAs("pT_dependence_jetPt_pp_pbpb.pdf");
  
  TH1D* hTempPt = new TH1D("hTemp",";p_{T}^{#gamma} (GeV);",200,10,300);
  TCanvas* c5 = new TCanvas("c5","",1200,700);
  makeMultiPanelCanvas(c5,nPtBin,2,0.0,0.0,0.2,0.15,0.02);
  for ( int ipt = 1 ; ipt<=nPtBin  ; ipt++) {
    c5->cd(ipt+nPtBin);
    // draw pp                                                                                                                                           
    hTempPt->SetXTitle("p_{T}^{Jet} (GeV)");
    hTempPt->SetYTitle("#frac{dN}{dp_{T}} #frac{1}{N}");
    hTempPt->SetAxisRange(10,150,"X");
    //        hTempPt->SetAxisRange(0,0.025,"Y");                                                                                                        
    hTempPt->SetAxisRange(0,0.07,"Y");
    handsomeTH1(hTempPt,0);
    hTempPt->DrawCopy();

    drawSys(dNdJetPt[1][ipt], dNdJetPtSys[1][ipt], kGreen);
    drawSys(dNdJetPt[3][ipt], dNdJetPtSys[3][ipt], kYellow);
    handsomeTH1(dNdJetPt[3][ipt],2);
    dNdJetPt[1][ipt]->Draw("same");
    dNdJetPt[3][ipt]->Draw("same");
    if ( ipt == 1 ) {
      TLegend *ly = new TLegend(0.4484643,0.6988445,0.9140673,0.9102941,NULL,"brNDC");
      easyLeg(ly,"2.76TeV");
      ly->AddEntry(meanJetPt[1],"PbPb 0-30%","p");
      ly->AddEntry(meanJetPt[3],"pp (smeared)","p");
      ly->Draw();
    }
 

    c5->cd(ipt);
    hTempPt->DrawCopy();

    drawSys(dNdJetPt[2][ipt], dNdJetPtSys[2][ipt], kGreen);
    drawSys(dNdJetPt[4][ipt], dNdJetPtSys[4][ipt], kYellow);
    handsomeTH1(dNdJetPt[4][ipt],2);
    dNdJetPt[4][ipt]->SetMarkerStyle(24);
    dNdJetPt[2][ipt]->Draw("same");
    dNdJetPt[4][ipt]->Draw("same");
    if ( ipt==1 ){
      TLegend *ly = new TLegend(0.4484643,0.5088445,0.9140673,0.75102941,NULL,"brNDC");
      easyLeg(ly,"2.76TeV");
      ly->AddEntry(meanJetPt[4],"PbPb 30-100%","p");
      ly->AddEntry(meanJetPt[2],"pp (smeared)","p");
      ly->Draw();
    }

void DrawTHStack(){
        static const int nColor = 9;
        static const int colorCode[nColor] = {
        1,2,3,6,4,8,1,9,5};
	gStyle->SetOptStat(kFALSE);
	int sth=0; int Gth=0;
	TFile *f = TFile::Open(outG);
	TObjArray *hs[3];
        std::vector<TH1D*> h[3];
	TCanvas *c1 = new TCanvas("c1","stacked hists",10,10,600,600);
	c1->SetLogy();
	TString dirname;
	TString name;	
        TH1D* hFrame = new TH1D("","",4000,0,4000);
	hFrame->GetXaxis()->SetRangeUser(0,500);
	hFrame->GetYaxis()->SetRangeUser(1,1e5);
	hFrame->GetXaxis()->SetTitle("N_{part}");
	hFrame->GetYaxis()->SetTitle("# of Events");
        hFrame->DrawCopy();
       	if(Gth==0)
        	name = "G0";
        else if(Gth<nGlau && Gth>0)
        	name = Form("Glau_%d",Gth);
        else
        	name = Form("bin_%d",Gth-nGlau+1);
        TVectorD *centbin = (TVectorD*)f->Get(Form("std/G0/centbin"));
        TVectorD *kpoint = (TVectorD*)f->Get(Form("std/G0/kpoint"));

	for(int sth=0; sth<1; sth++){
        	if(sth==0){dirname = "std";}
        	else if(sth==1){dirname ="Gri055";}
        	else {dirname ="Gri101";}
		hs[sth] = (TObjArray*)f->Get(Form("%s/%s/Npartdis",dirname.Data(),name.Data()));
                h[sth].resize(N);
                for(int i=0;i<N;i++){
//		h[sth][i] = (TH1D*)hs[sth]->GetStack()->At(N-1);
		h[sth][i] = (TH1D*)hs[sth]->At(i);
		h[sth][i]->SetLineColor(colorCode[i]);
		h[sth][i]->SetMarkerSize(1.0);
		h[sth][i]->SetMarkerStyle(20);
		h[sth][i]->SetMarkerColor(colorCode[i]);
		h[sth][i]->Draw("Psame");
//                hs[sth]->Draw("Psame");
                }
	}
    std::vector<TString> label(N);
	TLegend *leg0 = new TLegend(0.60,0.60,0.90,0.88);
    leg0->SetFillColor(10);
    leg0->SetBorderSize(0);
    leg0->SetTextFont(42);
    leg0->SetTextSize(0.030);
    for(int i=0;i<N;i++){
	if(method==0)label[i] = Form("%.2f-%.2f%%",(*centbin)[i]*100,(*centbin)[i+1]*100);
	else label[i] = Form("%.2f-%.2f",(*kpoint)[i],(*kpoint)[i+1]);
//    leg0->AddEntry(h[0][i],"standard","p");
    leg0->AddEntry(h[0][i],label[i],"p");
    }
 //   leg0->AddEntry(h[1],"Gribov #Omega=0.55","p");
 //   leg0->AddEntry(h[2],"Gribov #Omega=1.01","p");
        leg0->Draw();
c1->Print("Npartdis.pdf");
c1->Print("Npartdis.png");
}

//.........这里部分代码省略.........
                sum += NCharUps->GetBinContent(30000-i)*(NCharUps->GetBinCenter(30000-i));
                bins += NCharUps->GetBinContent(30000-i); }
        TrackAverage[c] = sum/bins;
	TotalUps += bins;
	cout << "Point " << c << ": " << bins << " +- " << sqrt(bins) << "   Upsilons" << "\n";
        min = i; }
	cout << NCharUps->GetBinCenter(100) << "\n";
       	// cout << 10*c << "-" << 10*(c+1) << "%" << "     " << TrackAverage[c] << "\n";  }
  cout << "Total number of Upsilons: " << TotalUps << "\n";
  cout << "\n" << "Corrected X/Y Track values" << "\n";
  double TrackYpoints[20];
  for(i = 0; i < 12; i++) {
	cout << (Ranges[i+1] - Ranges[i])*MinBiasEvents << "\n"; }

	//Take Average Values and Divide by ND Average to Get X points
 // for(c = 0; c < 4; c++) {
//	TrackAverage[c] = TrackAverage[c]/Tracks->GetMean();
//	TrackYpoints[c] = (50./NCharUps->Integral(0,25000))*NCharUps->Integral(TrackCent[c+1], TrackCent[c]-1 );
//	cout << "(" <<TrackAverage[c] << "," << TrackYpoints[c] <<  ")" <<  "\n"; }
  //for(c = 4; c < 7; c++) {
    //    TrackAverage[c] = TrackAverage[c]/Tracks->GetMean();
  //      TrackYpoints[c] = (25./NCharUps->Integral(0,25000))*NCharUps->Integral(TrackCent[c+1], TrackCent[c]-1 );
    //    cout << "(" <<TrackAverage[c] << "," << TrackYpoints[c] <<  ")" <<  "\n"; }
//  for(c = 7; c < 10; c++) {
 //       TrackAverage[c] = TrackAverage[c]/Tracks->GetMean();
   //     TrackYpoints[c] = (20./NCharUps->Integral(0,25000))*NCharUps->Integral(TrackCent[c+1], TrackCent[c]-1 );
//	cout << "(" <<TrackAverage[c] << "," << TrackYpoints[c] <<  ")" <<  "\n"; }


TrackAverage[0] = TrackAverage[0]/Tracks->GetMean();
        TrackYpoints[0] = (50./NCharUps->Integral(0,30000))*NCharUps->Integral(100*TrackCent[1], 100*(TrackCent[0]-1) );
	cout << "(" <<TrackAverage[0] << "," << TrackYpoints[0] <<  ")" <<  "\n";

TrackAverage[1] = TrackAverage[1]/Tracks->GetMean();
        TrackYpoints[1] = (16.67/NCharUps->Integral(0,30000))*NCharUps->Integral(100*TrackCent[2], 100*TrackCent[1]-100 );
	cout << "(" <<TrackAverage[1] << "," << TrackYpoints[1] <<  ")" <<  "\n";

TrackAverage[2] = TrackAverage[2]/Tracks->GetMean();
        TrackYpoints[2] = (8.33/NCharUps->Integral(0,30000))*NCharUps->Integral(100*TrackCent[3], 100*TrackCent[2]-100 );
	cout << "(" <<TrackAverage[2] << "," << TrackYpoints[2] <<  ")" <<  "\n";

 for(c = 3; c < 14; c++) {
        TrackAverage[c] = TrackAverage[c]/Tracks->GetMean();
        TrackYpoints[c] = (10./NCharUps->Integral(0,30000))*NCharUps->Integral(100*TrackCent[c+1], 100*TrackCent[c]-100 );
        cout << "(" <<TrackAverage[c] << "," << TrackYpoints[c] <<  ")" <<  "\n"; } 
  TF1* line = new TF1("line", "x", 0, 4.5);
  TF1* line2 = new TF1("line2", "1.8*x", 0, 4.5);
  line2->SetLineColor(kWhite);
  sum = 0;
  double scaledSumE = 0, scaleSumT = 0;
  scaleSumT = NCharUps->Integral(0,30000);
  TCanvas* can2 = new TCanvas("can2", "Centrality");
  can2->cd();
  line->SetLineColor(1);
  line2->GetYaxis()->SetTitle("#varUpsilon(1S)/<#varUpsilon(1S)>");
  line2->GetXaxis()->SetTitle("#Sigma E_{T}^{|#eta| < 2.4}/<#Sigma E_{E}^{|#eta| < 2.4}>");
  line2->SetTitle("#varUpsilon(1S) Cross Section vs Track Multiplicity");
  line2->Draw();
  line->SetLineStyle(3);
  line->Draw("same");
   for(c = 0; c < 18; c++) { cout << TrackCent[c] << "           Bin Average: " << 10.32*TrackAverage[c] <<  "\n";}
  double Test, Testaverage, Testsum, Testbin, count = 0, TestY = 0;
  Test = NCharUps->Integral(38,100);
  for(i = 38; i < 100; i++) {
	Testsum += NCharUps->GetBinContent(i)*NCharUps->GetBinCenter(i);
	Testbin += NCharUps->GetBinContent(i); }
  Testaverage = Testsum/Testbin;
  TestY = (10./NCharUps->Integral(0,30000))*NCharUps->Integral(38, 100);
  cout << "\n" << "\n" << TestY << "\n" << "\n";

  TH1D* modelT = new TH1D("modelT", "Pythia model", 5000, 0, 5);
  //Fill Histogram with points from Model
  for(i = 0; i < 20; i++) {
	//Track Values
	modelT->SetBinContent(1000*TrackAverage[i], TrackYpoints[i]);
        modelT->SetBinError(1000*TrackAverage[i], .01); }
  modelT->SetMarkerStyle(20);
  modelT->SetMarkerColor(kBlue);
  modelT->SetLineColor(kBlue);
  modelT->Draw("sameE");
  //Histogram for STAR Data Points
  TH1D* data = new TH1D("data", "CMS Data", 4000, 0, 4);
  data->SetMarkerStyle(20);
  data->SetMarkerColor(kRed);
  data->SetLineColor(kRed);
  data->SetBinContent(3262, 6.67);
  data->SetBinError(3262, .26);
  data->SetBinContent(2010, 3.12);
  data->SetBinError(2010, .15);
  data->SetBinContent(1240, 1.4);
  data->SetBinError(1240, .06);
  data->SetBinContent(630, .240);
  data->SetBinError(630, .01);
  data->Draw("sameE");
  TLegend* leg = new TLegend(.1,.74,.30,.9);
  leg->SetHeader("For |#eta| < 2.4");
  leg->AddEntry(modelT->DrawCopy("same"), "Pythia 8 Model", "l");
  leg->AddEntry(data->DrawCopy("same"), "CMS data", "l");
  leg->Draw("same");
  return; }