void chisq_integral_compare(void)
{
	// Get Histograms from files
	TFile *f = new TFile("chisq_integral_option1.root");
	TH1D *option1 = (TH1D*)gROOT->FindObject("chisq_corrected_integral");
	
	f = new TFile("chisq_integral_option2.root");
	TH1D *option2 = (TH1D*)gROOT->FindObject("chisq_corrected_integral");

	f = new TFile("chisq_integral_close_packed.root");
	TH1D *close_packed = (TH1D*)gROOT->FindObject("chisq_corrected_integral");

	f = new TFile("chisq_integral_phase_shifted.root");
	TH1D *phase_shifted = (TH1D*)gROOT->FindObject("chisq_corrected_integral");

	f = new TFile("chisq_integral_close_packed_less_stereo.root");
	TH1D *close_packed_less_stereo = (TH1D*)gROOT->FindObject("chisq_corrected_integral");
	
	option1->SetFillStyle(3000);
	option2->SetFillStyle(3000);
	close_packed->SetFillStyle(3000);
	phase_shifted->SetFillStyle(3000);
	close_packed_less_stereo->SetFillStyle(3000);
	
	option1->SetFillColor(kGreen);
	option2->SetFillColor(kCyan);
	close_packed->SetFillColor(kBlue);
	phase_shifted->SetFillColor(kMagenta);
	close_packed_less_stereo->SetFillColor(kRed);

	// Create canvas and axes to draw on
	TCanvas *c1 = new TCanvas("c1");
	c1->SetTickx();
	c1->SetTicky();
	c1->SetGridy();
	c1->SetGridx();
	c1->SetLogx();
	
	TH2D *axes = new TH2D("axes","Efficiency vs. cut on #chi^{2}/dof (scaled)", 100, 0.7, 1000.0, 100, 0.5, 1.05);
	axes->SetStats(0);
	axes->SetXTitle("#chi^{2}/dof cut using uncertainties scaled by Gaussian fit to pull distributions");
	axes->SetYTitle("Efficiency (fraction)");
	axes->Draw();

	close_packed->Draw("C same");
	close_packed_less_stereo->Draw("C same");
	option2->Draw("C same");
	phase_shifted->Draw("C same");
	option1->Draw("C same");

	TLegend *leg = new TLegend(0.50, 0.20, 0.85, 0.45);
	leg->AddEntry(close_packed	, "CDC close packed  (N_{stereo}=11)");
	leg->AddEntry(close_packed_less_stereo	, "CDC close packed  (N_{stereo}=8)");
	leg->AddEntry(option2		, "CDC option-2_{shifted}   (N_{stereo}=9)");
	leg->AddEntry(phase_shifted, "CDC option-1_{shifted}    (N_{stereo}=8)");
	leg->AddEntry(option1		, "CDC option-1          (N_{stereo}=8)");
	leg->Draw();
	
	AddStandardLabels(axes);

	c1->SaveAs("chisq_integral_compare.gif");
	c1->SaveAs("chisq_integral_compare.pdf");
}

//----------------
// AddStandardLabels
//----------------
void AddStandardLabels(TH2D *axes=NULL)
{
	// This will draw a label or two on the
	// current plot using the NDC coordinates.
	// It is put here to make sure all plots have
	// a consistent labeling.
	
	// Date, Author
	TLatex *lab = new TLatex(0.7, 0.675, "August 21, 2008 DL");
	ConvertFromNDC(lab, axes);
	lab->SetTextSize(0.03);
	lab->SetTextAlign(33);
	lab->Draw();

	// SVN Revision
	lab = new TLatex(0.7, 0.645, "svn revision: 4096");
	ConvertFromNDC(lab, axes);
	lab->SetTextSize(0.02);
	lab->SetTextAlign(31);
	lab->Draw();

	// Candidate type
	lab = new TLatex(0.7, 0.635, "Candidates from finder");
	ConvertFromNDC(lab, axes);
	lab->SetTextSize(0.02);
	lab->SetTextAlign(32);
	lab->Draw();
	
	// Event type
	lab = new TLatex(0.7, 0.610, "Single #pi^{+} with MULS, LOSS, DRAY on");
	ConvertFromNDC(lab, axes);
	lab->SetTextSize(0.03);
	lab->SetTextAlign(31);
	lab->Draw();
}

//----------------
// ConvertFromNDC
//----------------
void ConvertFromNDC(TLatex *obj, TH2D *h=NULL)
{
	// Bugs in ROOT make it hard to plot labels consistently.
	// For 1D plots, the histogram axes define the coordinate
	// system. For 2D plots, we seem to be forced to use the
	// NDC. There does not seem to be an obvious way to tell
	// which we're using so we pass the information in in the
	// form of the "axes" histogram. If it is not NULL, then
	// we use it to define the limits. Otherwise, we do nothing.

	if(h==NULL)return;

	TAxis *xaxis = h->GetXaxis();
	int Nbinsx = xaxis->GetNbins();
	double xmin = xaxis->GetBinLowEdge(1);
	double xmax = xmin + xaxis->GetBinLowEdge(Nbinsx);

	TAxis *yaxis = h->GetYaxis();
	int Nbinsy = yaxis->GetNbins();
	double ymin = yaxis->GetBinLowEdge(1);
	double ymax = ymin + yaxis->GetBinLowEdge(Nbinsy);
	
	
	double x = obj->GetX();
	double y = obj->GetY();
	cout<<" in:  x="<<x<<"  y="<<y<<endl;
	
	x = xmin + (xmax-xmin)*(0.5+x/1.15);
	y = ymin + (ymax-ymin)*(0.0+y/1.15);
	
	obj->SetX(x);
	obj->SetY(y);

	cout<<"out:  x="<<x<<"  y="<<y<<endl;
}