TGraphErrors; Unsure of x and y errors

ROOT
1

Hello,
Very new beginner here- I’m working on a simple program that does many things, one of which is to plot RMS of Px, Py, Pz versus Pt. There should be 10 points on each graph with error bars, clearly by using TGraphErrors. In the tutorial for TGraphErros, it shows values of error on x and y values as an array containing numbers, but I’m not sure how they got those values, nor how I can get them in my program. Could somebody please help me figure out how to get the errors on x and y values (RMS of Px, Py, Pz and Pt) so that they display as error bars through TGraphErrors in my particular program? I’ve copied my program below, but the lines of interest are at the end. I’ve only included the lines dealing with the graph for RMS of Px vs. Pt, but will add the other two based on what I learn for the Px case. Thank you!


#define merged_cxx
using namespace std;

#include <vector>
#include <iostream>
#include "merged.h"
#include "TH1F.h"
#include "TH2.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TBranch.h"
#include "TF1.h"
#include "THStack.h"
#include "TObjArray.h"
#include "TList.h"
#include "TGraphErrors.h"

void merged::Loop()
{
   if (fChain == 0) return;

   Long64_t nentries = fChain->GetEntries();
   cout << "Number of events = " << nentries << endl;
   Long64_t nbytes = 0, nb = 0;

   //Create a canvas                                                                                                                                                        
   TCanvas*c1 = new TCanvas("c1", "Pions",1000,2000);

   //Declaring histograms for truth table pions                                                                                                                             
   TH1F*hist_px = new TH1F("hist_px", "Px", 50, -0.5, 0.5);
   hist_px->SetMarkerStyle(21);
   hist_px->SetMarkerColor(kBlue);

   TH1F*hist_py = new TH1F("hist_py", "Py", 50, -0.5, 0.5);
   hist_py->SetMarkerStyle(21);
   hist_py->SetMarkerColor(kBlue);

   TH1F*hist_pz = new TH1F("hist_pz", "Pz", 50, -0.5, 0.5);
   hist_pz->SetMarkerStyle(21);
   hist_pz->SetMarkerColor(kBlue);
 TH1F*hist_E = new TH1F("hist_E", "E" , 40, -5., 5.);
   hist_E->SetMarkerStyle(21);
   hist_E->SetMarkerColor(kBlue);


   //Declaring histograms for reco pions                                                                                                                                    
   TH1F*hist_px2 = new TH1F("hist_px2", "Px2", 50, -0.05, 0.05);
   hist_px2->SetMarkerStyle(21);
   hist_px2->SetMarkerColor(kGreen);

   TH1F*hist_py2 = new TH1F("hist_py2", "Py2", 50, -0.05, 0.05);
   hist_py2->SetMarkerStyle(21);
   hist_py2->SetMarkerColor(kGreen);

   TH1F*hist_pz2 = new TH1F("hist_pz2", "Pz2", 50, -0.05, 0.05);
   hist_pz2->SetMarkerStyle(21);
   hist_pz2->SetMarkerColor(kGreen);


   TH1F * hist_log10chisq = new TH1F("hist_log10chisq", "log10chisq", 50, 0., 2.5);
   TH1F * hist_respx = new TH1F("hist_respx","respx", 50, -0.05, 0.05);
   TH1F * hist_respy = new TH1F("hist_respy","respy", 50, -0.05, 0.05);
   TH1F * hist_respz = new TH1F("hist_respz","respz", 50, -0.05, 0.05);
 //histograms for bins of Pt for Px                                                                                                                                       
   TH1F * hrespx[10];
   TH1F * hrespy[10];
   TH1F * hrespz[10];
   for(int i = 0; i < 10; i++){
     string hxname = string("hrespx_") + to_string(i);
     hrespx[i] = new TH1F(hxname.c_str(), "px resolution", 100, -0.05, 0.05);
     string hyname = string("hrespy_") + to_string(i);
     hrespy[i] = new TH1F(hyname.c_str(), "py resolution", 100, -0.05, 0.05);
     string hzname = string("hrespz_") + to_string(i);
     hrespz[i] = new TH1F(hzname.c_str(), "pz resolution", 100, -0.05, 0.05);
   } // End for i = 0 ... 10                                                     
 double errsq = 0.0135 * 0.0135;
   //the following for loop loops over candidates                                                                                                                           
   for (Long64_t jentry=0; jentry<nentries;jentry++) {
      Long64_t ientry = LoadTree(jentry);
      if (ientry < 0) break;
      nb = fChain->GetEntry(jentry);   nbytes += nb;
      //Loops over particles in h16 and fills hist with h16 data                                                                                                            
      for(int i = 0; i < (int)h16_idhep->size(); i++){
        if((*h16_idhep)[i]==211){
          hist_px->Fill((*h16_px)[i]);
          hist_py->Fill((*h16_py)[i]);
          hist_pz->Fill((*h16_pz)[i]);
          hist_E->Fill((*h16_e)[i]);
        }
      }
      //Loops over particles in h19 and fills hists                                                                                                                         
      for(int i19 = 0; i19 < (int)h19_bpcand->size(); i19++){
        hist_px2->Fill((*h19_pipx)[i19]);
        hist_py2->Fill((*h19_pipy)[i19]);
        hist_pz2->Fill((*h19_pipz)[i19]);
        for(int i16 = 0; i16 < (int)h16_idhep->size(); i16++){
          if((*h16_idhep)[i16]==211){
            double chisq =
              ((*h16_px)[i16] - (*h19_pipx)[i19]) * ((*h16_px)[i16] - (*h19_pipx)[i19]) / errsq +
              ((*h16_py)[i16] - (*h19_pipy)[i19]) * ((*h16_py)[i16] - (*h19_pipy)[i19]) / errsq +
              ((*h16_pz)[i16] - (*h19_pipz)[i19]) * ((*h16_pz)[i16] - (*h19_pipz)[i19]) / errsq;
            double log10chisq = log10(chisq);
            double Pt = sqrt(((*h16_px)[i16]) * ((*h16_px)[i16]) + ((*h16_py)[\
                                                                               i16]) * ((*h16_py)[i16]));
                                                      
hist_log10chisq->Fill(log10chisq);
            if(chisq < 10.){
              double res_px = ((*h16_px)[i16] - (*h19_pipx)[i19]);
              double res_py = ((*h16_py)[i16] - (*h19_pipy)[i19]);
              double res_pz = ((*h16_pz)[i16] - (*h19_pipz)[i19]);

              hist_respx->Fill(res_px);
              hist_respy->Fill(res_py);
              hist_respz->Fill(res_pz);

              if(0 < Pt && Pt < 2.5){
                int ipt = Pt/0.25;
                hrespx[ipt]->Fill(res_px);
                hrespy[ipt]->Fill(res_py);
                hrespz[ipt]->Fill(res_pz);
              }
            }
          } // End if((*h16_idhep)[i]==211)                                                                                                                                 
        } // End for i16 = 0 ... (int)h16_idhep->size()                                                                                                                     
      }
   }
    hist_px->Draw("P PLC PMC");
   hist_px2->Draw("SAME PLC PMC");
   c1->SaveAs("Px_of_pion.png");
   c1->Clear();

   hist_py->Draw("P PLC PMC");
   hist_py2->Draw("SAME P PLC PMC");
   c1->SaveAs("Py_of_pion.png");
   c1->Clear();

   hist_pz->Draw("P PLC PMC");
   hist_pz2->Draw("SAME PLC PMC");
   c1->SaveAs("Pz_of_pion.png");
   c1->Clear();

   hist_E->Draw("P PLC PMC");
   // hist_E2->Draw("SAME PLC PMC");                                                                                                                                        
   c1->SaveAs("E_of_pion.png");
   c1->Clear();

   hist_log10chisq->Draw("log10chisq");
   c1->SaveAs("log10chisq.png");
   c1->Clear();

   hist_respx->Draw();
   c1->SaveAs("res_px.png");
   c1->Clear();

   hist_respy->Draw();
   c1->SaveAs("res_py.png");
   c1->Clear();

   hist_respz->Draw();
   c1->SaveAs("res_pz.png");
   c1->Clear();

   double ptval[10], rmspx[10], rmspy[10], rmspz[10];
   for(int i = 0; i < 10; i++){
     ptval[i] = 0.125 + (double)i * 0.25;

     hrespx[i]->Draw();
     rmspx[i] = hrespx[i]->GetRMS();
     string hxname = string("hrespx_") + to_string(i) + string(".png");
     c1->SaveAs(hxname.c_str());
     c1->Clear();

     hrespy[i]->Draw();
     rmspy[i] = hrespy[i]->GetRMS();
     string hyname = string("hrespy_") + to_string(i) + string(".png");
     c1->SaveAs(hyname.c_str());
     c1->Clear();

     hrespz[i]->Draw();
     rmspz[i] = hrespz[i]->GetRMS();
     string hzname = string("hrespz_") + to_string(i) + string(".png");
     c1->SaveAs(hzname.c_str());
     c1->Clear();
 } // End for i = 0 ... 10                                                                                                                                                
   const Int_t i = 10;
   double ex1[i] = {};//??                                                                                                                                                  
   double ey1[i] = {};//??                                                                                                                                                  

   TGraphErrors *gr1 = new TGraphErrors(i,ptval,rmspx,ex1,ey1);
   gr1->SetTitle("Px RMS vs. Pt");
   gr1->SetMarkerColor(4);
   gr1->SetMarkerStyle(21);
   gr1->Draw();
   c1->SaveAs("tgraph.png");
   c1->Update;
   // Plot RMS values vs pt using TGraphErrors                                                                                                                              
}
2

Try:

   double ptval[10],   rmspx[10],   rmspy[10],   rmspz[10];
   double ptval_e[10], rmspx_e[10], rmspy_e[10], rmspz_e[10];
   for(int i = 0; i < 10; i++) {
     ptval[i] = 0.125 + (double)i * 0.25;
     ptval_e[i] = 0.25 / 2.0; // half of the "step"
     rmspx[i] = hrespx[i]->GetRMS();
     rmspx_e[i] = hrespx[i]->GetRMSError();
     rmspy[i] = hrespy[i]->GetRMS();
     rmspy_e[i] = hrespy[i]->GetRMSError();
     rmspz[i] = hrespz[i]->GetRMS();
     rmspz_e[i] = hrespz[i]->GetRMSError();
   }
   TGraphErrors *gr1 = new TGraphErrors(10, ptval, rmspx, ptval_e, rmspx_e);
   gr1->SetTitle("Px RMS vs. Pt");
   TGraphErrors *gr2 = new TGraphErrors(10, ptval, rmspy, ptval_e, rmspy_e);
   gr2->SetTitle("Py RMS vs. Pt");
   TGraphErrors *gr3 = new TGraphErrors(10, ptval, rmspz, ptval_e, rmspz_e);
   gr3->SetTitle("Pz RMS vs. Pt");

or maybe even better:

   TGraphErrors *gr1 = new TGraphErrors(10); gr1->SetTitle("Px RMS vs. Pt");
   TGraphErrors *gr2 = new TGraphErrors(10); gr2->SetTitle("Py RMS vs. Pt");
   TGraphErrors *gr3 = new TGraphErrors(10); gr3->SetTitle("Pz RMS vs. Pt");
   for(int i = 0; i < 10; i++) {
     Double_t pt = 0.125 + (double)i * 0.25;
     Double_t pt_e = 0.25 / 2.0; // half of the "step"
     gr1->SetPoint(i, pt, hrespx[i]->GetRMS());
     gr1->SetPointError(i, pt_e, hrespx[i]->GetRMSError());
     gr2->SetPoint(i, pt, hrespy[i]->GetRMS());
     gr2->SetPointError(i, pt_e, hrespy[i]->GetRMSError());
     gr3->SetPoint(i, pt, hrespz[i]->GetRMS());
     gr3->SetPointError(i, pt_e, hrespz[i]->GetRMSError());
   }
1 Like
3

That works wonderfully, thank you so much!!

4

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.