Set TRatioPlot ratio error

bethlong06 · May 27, 2023, 4:50pm

I want to set the errors on the ratio of my TRatioPlot by hand. I’ve calculated the errors for each point and they’re saved in hRatioErrorPerPoint. With the code below, I think the errors I create are drawn in a separate graph on top of the automatic graph that’s generated by the TRatioPlot, but since the automatic ones are up to 13 orders of magnitude greater than the ones I calculate, I can only see the automatic ones.

How can I force the errors shown to be those from hRatioErrorPerPoint?

void TRatioPlotMaker(TH1* hD, TH1* hM, TCanvas* cRatio, TH1F* hRatioErrorPerPoint)
{
  //Plot ratio Data/MC with TRatioPlot
  TRatioPlot *hRatioPlot = new TRatioPlot(hD,hM);
  hRatioPlot->SetH1DrawOpt("E");
  hRatioPlot->SetH2DrawOpt("HIST");

  cRatio->cd();
  hRatioPlot->Draw();
  
  TGraphErrors* gRatioWithErrors = (TGraphErrors*)hRatioPlot->GetLowerRefGraph();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    gRatioWithErrors->SetPointError(ii,1,hRatioErrorPerPoint->GetBinContent(ii+1));

  hRatioPlot->GetLowerPad()->cd();
  gRatioWithErrors->Draw("P");

Wile_E_Coyote · May 27, 2023, 5:02pm

Try:

hRatioPlot->GetLowerPad()->cd(); gPad->Clear();
delete gRatioWithErrors->GetHistogram(); gRatioWithErrors->SetHistogram(0);
gRatioWithErrors->Draw("AP");
gPad->Modified(); gPad->Update();

bethlong06 · May 27, 2023, 5:49pm

This still doesn’t work, it draws the axis twice, like I asked about here.

On further investigation, it seems that the error isn’t being set properly:

  TGraphErrors* gRatioWithErrors = (TGraphErrors*)hRatioPlot->GetLowerRefGraph();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      gRatioWithErrors->SetPointError(ii,1,hRatioErrorPerPoint->GetBinContent(ii+1));
      std::cout<<ii<<" "<<hRatioErrorPerPoint->GetBinContent(ii+1)<<" "<<gRatioWithErrors->GetErrorY(ii)<<std::endl;
    }

0 2.24397 68.6574
1 0.962936 84.5099
2 2.11517 13009.7
3 26.4041 577946
4 12.1803 5266.68
5 7.8174 3485.87
6 1.92447 3681.33
7 1.86803 139.855
8 7.75246 7.48251e+08
9 2.12146 18.201
10 115.378 2.35133e+17
...

Wile_E_Coyote · May 27, 2023, 6:35pm

bethlong06 · May 27, 2023, 9:13pm

I’ll try

bethlong06 · May 28, 2023, 8:58am

Ok, for example this shows the behaviour:

void MWE(){

  TH1F* h1 = new TH1F("h1","h1",100,-3,3);
  h1->FillRandom("gaus",10000000);

  TH1F* h2 = new TH1F("h2","h2",100,-3,3);
  h2->FillRandom("gaus",100);

  TCanvas* cRatio = new TCanvas("cRatio","cRatio",900,700);
  TRatioPlot* rRatio = new TRatioPlot(h2,h1);
  rRatio->SetH1DrawOpt("E");
  rRatio->SetH2DrawOpt("HIST");
  rRatio->Draw();

  TGraphErrors* gRatioWithErrors = (TGraphErrors*)rRatio->GetLowerRefGraph();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      gRatioWithErrors->SetPointError(ii,1,5);
      //      std::cout<<ii<<" "<<hRatioErrorPerPoint->GetBinContent(ii+1)<<" "<<gRatioWithErrors->GetErrorY(ii)<<std::endl;
    }
  rRatio->GetLowerPad()->cd();
  gPad->Clear();
  gRatioWithErrors->Draw("AP");
  gPad->Modified(); gPad->Update();
}

bethlong06 · May 28, 2023, 9:07am

Here, nothing is drawn on the bottom pad

void MWE(){

  TH1F* h1 = new TH1F("h1","h1",100,-3,3);
  h1->FillRandom("gaus",10000000);

  TH1F* h2 = new TH1F("h2","h2",100,-3,3);
  h2->FillRandom("gaus",100);

  TCanvas* cRatio = new TCanvas("cRatio","cRatio",900,700);
  TRatioPlot* rRatio = new TRatioPlot(h2,h1);
  rRatio->SetH1DrawOpt("E");
  rRatio->SetH2DrawOpt("HIST");
  rRatio->Draw();

  TGraphErrors* gRatioWithErrors = (TGraphErrors*)rRatio->GetLowerRefGraph();
  TGraphErrors* gNewRatio = new TGraphErrors();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      gNewRatio->SetPointError(ii,1,5);
      gNewRatio->SetPoint(ii,ii+60,gRatioWithErrors->GetPointY(ii));
      //      std::cout<<ii<<" "<<hRatioErrorPerPoint->GetBinContent(ii+1)<<" "<<gRatioWithErrors->GetErrorY(ii)<<std::endl;
    }
  rRatio->GetLowerPad()->cd();
  rRatio->GetLowerPad()->Clear();
  gNewRatio->Draw("P");
  gPad->Modified(); gPad->Update();

}

Wile_E_Coyote · May 28, 2023, 9:14am

I think the problem is that you actually have a TGraphAsymmErrors:
std::cout << rRatio->GetLowerRefGraph()->ClassName() << std::endl;

bethlong06 · May 28, 2023, 10:02am

Making gRatioWithErrors and gNewRatio into TGraphAsymmErrors still gives the empty bottom plot

  TGraphAsymmErrors* gRatioWithErrors = (TGraphAsymmErrors*)rRatio->GetLowerRefGraph();
  TGraphAsymmErrors* gNewRatio = new TGraphAsymmErrors();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      gNewRatio->SetPointError(ii,1,1,5,5);
      gNewRatio->SetPoint(ii,ii,gRatioWithErrors->GetPointY(ii));
      //      std::cout<<ii<<" "<<hRatioErrorPerPoint->GetBinContent(ii+1)<<" "<<gRatioWithErrors->GetErrorY(ii)<<std::endl;
    }

bethlong06 · May 28, 2023, 10:07am

Ok, there were a few bugs. This finally works:

void MWE(){

  TH1F* h1 = new TH1F("h1","h1",100,-3,3);
  h1->FillRandom("gaus",10000000);

  TH1F* h2 = new TH1F("h2","h2",100,-3,3);
  h2->FillRandom("gaus",100);

  TCanvas* cRatio = new TCanvas("cRatio","cRatio",900,700);
  TRatioPlot* rRatio = new TRatioPlot(h2,h1);
  rRatio->SetH1DrawOpt("E");
  rRatio->SetH2DrawOpt("HIST");
  rRatio->Draw();

  TGraphErrors* gRatioWithErrors = (TGraphErrors*)rRatio->GetLowerRefGraph();
  TGraphErrors* gNewRatio = new TGraphErrors();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      //   gNewRatio->SetPointError(ii,1,1,5,5);
      gNewRatio->SetPointError(ii,1,5);
      gNewRatio->SetPoint(ii,gRatioWithErrors->GetPointX(ii),gRatioWithErrors->GetPointY(ii));
      //      std::cout<<ii<<" "<<hRatioErrorPerPoint->GetBinContent(ii+1)<<" "<<gRatioWithErrors->GetErrorY(ii)<<std::endl;
    }
  rRatio->GetLowerPad()->cd();
  rRatio->GetLowerPad()->Clear();
  gNewRatio->Draw("AP");
  gPad->Modified(); gPad->Update();

}

bethlong06 · May 28, 2023, 10:13am

This still redraws the axis labels:

  TH1F* h1 = new TH1F("h1","h1",100,-3,3);
  h1->FillRandom("gaus",10000000);

  TH1F* h2 = new TH1F("h2","h2",100,-3,3);
  h2->FillRandom("gaus",100);

  TCanvas* cRatio = new TCanvas("cRatio","cRatio",900,700);
  TRatioPlot* rRatio = new TRatioPlot(h1,h2);
  rRatio->SetH1DrawOpt("E");
  rRatio->SetH2DrawOpt("HIST");
  rRatio->Draw();

  TGraphErrors* gRatioWithErrors = (TGraphErrors*)rRatio->GetLowerRefGraph();
  TGraphErrors* gNewRatio = new TGraphErrors();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      //   gNewRatio->SetPointError(ii,1,1,5,5);
      gNewRatio->SetPointError(ii,1,5);
      gNewRatio->SetPoint(ii,gRatioWithErrors->GetPointX(ii),gRatioWithErrors->GetPointY(ii));
      //      std::cout<<ii<<" "<<hRatioErrorPerPoint->GetBinContent(ii+1)<<" "<<gRatioWithErrors->GetErrorY(ii)<<std::endl;
    }
  rRatio->GetLowerPad()->cd();
  rRatio->GetLowerPad()->Clear();
  gNewRatio->Draw("AP");
  gPad->Modified(); gPad->Update();
}

Wile_E_Coyote · May 28, 2023, 10:15am

TGraph *gRatioWithErrors = rRatio->GetLowerRefGraph();
or:
TGraphAsymmErrors *gRatioWithErrors = (TGraphAsymmErrors*)rRatio->GetLowerRefGraph();

bethlong06 · May 28, 2023, 10:29am

error: cannot initialize a variable of type 'TGraphErrors *' with an rvalue of type 'TGraph *'
  TGraphErrors* gRatioWithErrors = rRatio->GetLowerRefGraph();

I’ve tried with TGraphAsymmErrors but the axis labels still get redrawn

Wile_E_Coyote · May 28, 2023, 10:48am

I’m afraid the TRatioPlot always tries to draw its own TGraphAsymmErrors in the bottom pad (and that’s why you get two y-axes).

If you want to “completely replace” the bottom pad, try to play with: ${ROOTSYS}/tutorials/hist/ratioplotOld.C

bethlong06 · May 28, 2023, 10:53am

I wouldn’t mind it drawing its own TGraphAsymmErrors as long as I could replace the errors it was trying to draw. But I’ll have a look at that macro, thank you for your help

Wile_E_Coyote · May 28, 2023, 10:54am

You should be able to modify the “original” TGraphAsymmErrors errors.

bethlong06 · May 28, 2023, 11:01am

Ok yes, I’ve finally managed. For some reason the Y axis range doesn’t get updated but I can set it by hand. This finally works. Thank you for the help!

void MWE(){

  TH1F* h1 = new TH1F("h1","h1",100,-3,3);
  h1->FillRandom("gaus",10000000);

  TH1F* h2 = new TH1F("h2","h2",100,-3,3);
  h2->FillRandom("gaus",100);

  TCanvas* cRatio = new TCanvas("cRatio","cRatio",900,700);
  TRatioPlot* rRatio = new TRatioPlot(h1,h2);
  rRatio->SetH1DrawOpt("E");
  rRatio->SetH2DrawOpt("HIST");
  rRatio->Draw();

  TGraphAsymmErrors* gRatioWithErrors = (TGraphAsymmErrors*)rRatio->GetLowerRefGraph();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      gRatioWithErrors->SetPointError(ii,1,1,5,5);
    }
}

Wile_E_Coyote · May 28, 2023, 11:02am

You need:
delete gRatioWithErrors->GetHistogram(); gRatioWithErrors->SetHistogram(0);

bethlong06 · May 28, 2023, 11:05am

works a treat

bethlong06 · May 28, 2023, 11:06am

Final version:

void MWE(){

  TH1F* h1 = new TH1F("h1","h1",100,-3,3);
  h1->FillRandom("gaus",10000000);

  TH1F* h2 = new TH1F("h2","h2",100,-3,3);
  h2->FillRandom("gaus",100);

  TCanvas* cRatio = new TCanvas("cRatio","cRatio",900,700);
  TRatioPlot* rRatio = new TRatioPlot(h1,h2);
  rRatio->SetH1DrawOpt("E");
  rRatio->SetH2DrawOpt("HIST");
  rRatio->Draw();

  TGraphAsymmErrors* gRatioWithErrors = (TGraphAsymmErrors*)rRatio->GetLowerRefGraph();
  for(int ii=0;ii<gRatioWithErrors->GetN();ii++)
    {
      gRatioWithErrors->SetPointError(ii,1,1,5,5);
    }
}
delete gRatioWithErrors->GetHistogram(); gRatioWithErrors->SetHistogram(0);