How to get errors and make a plot with it

Roofit and RooStats
1

hDiff is the subtraction of Electron from Pion. That’s working fine now, thanks to you all. Right now, I’m trying to get the actual values of hDiff and its errors, then loop over it and make a plot of momentum on the x-axis versus value/error of each point on the y-axis. This is the code I’m trying to make work:

 histo2 ->hDiff->clone();
  for (i=0;GetNBins;i++){
    value = hDiff->GetBin(i);
    error->GetBinError(i);
    histo2->SetBin(value/error,i);
  }

This is the larger code.
Thanks in anticipation for a favorable response

{

  gStyle->SetTitleSize(0.05,"XY");
  gStyle->SetTitleOffset(.9, "XY");
  gStyle->SetStatY(0.9);
  gStyle->SetStatX(0.9);
  
  hDummy = new TH2D("hDummy","hDummy",1000,0,80,1000,0,20);
  hDummy->Draw("axis");
  hDummy->GetXaxis()->SetTitle("TOF [ns]");
  hDummy->GetYaxis()->SetTitle("Momentum [MeV/c]");
  hDummy->GetXaxis()->CenterTitle();
  hDummy->GetYaxis()->CenterTitle();
  TCanvas *c1 = new TCanvas("c1", "c1",18,41,700,500);
  c1->Range(24,50,40,70);
  TH2D *Electron = new TH2D("Electron","",100,0.,2.0,250,30,1000);
  for(int i=0; i<10000; i++){
    double d =10;
    double Me=0.000511;
    double p =gRandom->Uniform(0., 2.0);
    double x = TMath::Power((p/Me),2);
    double v1 = sqrt(x/(1+x));
    double v = sqrt(x/(1+x))*(3*1e8);
    double t = (d/v)/ 0.000000001;
    double Nrandomt= gRandom->Gaus(0,0.2);
    double tof= t+Nrandomt;
    
    
    Electron->Fill(p,tof);
    
  }

  TH2D *Pion = new TH2D("Pion","",100,0.,2.0,250,30,1000);
  for(int i=0; i<10000; i++){
    double d =10;
    double Mp=0.136;
    double p =gRandom->Uniform(0., 2.0);
    double x = TMath::Power((p/Mp),2);
    double v1 = sqrt(x/(1+x));
    double v = sqrt(x/(1+x))*(3*1e8);
    double t = (d/v)/ 0.000000001;
    double Nrandomt= gRandom->Gaus(0,0.2);
    double tof= t+Nrandomt;
    
    
    Pion->Fill(p,tof);
    
  }
 
  TH1D *h1Electron = Electron->ProfileX();
  TH1D *h1Pion = Pion->ProfileX();
  h1Pion->Print("all");
  h1Electron->Print("all");
  TH1D* hDiff = (TH1D*)h1Pion->Clone("hDiff");
  hDiff->Add(h1Electron, -1);
  hDiff->Print("all");
  hDiff->Draw();




  histo2 ->hDiff->clone();
  for (i=0;GetNBins;i++){
    value = hDiff->GetBin(i);
    error->GetBinError(i);
    histo2->SetBin(value/error,i);
  }

  histo2->Draw("SAME");
  // hDiff->Print("all");
  
}
2

Try:

  TH1D *histo2 = ((TH1D*)(hDiff->Clone("histo2")));
  histo2->SetLineColor(kRed);
  for (int i = 0; i <= histo2->GetNbinsX() + 1; i++) {
    Double_t content = histo2->GetBinContent(i);
    Double_t error = histo2->GetBinError(i);
    if (error != 0.) histo2->SetBinContent(i, content / error);
    else if (content == 0.) histo2->SetBinContent(i, 1.); // assume "0/0 = 1"
    else histo2->SetBinContent(i, 0.); // should be set to "infinity"?
  }
  histo2->Draw("SAME"); // "SAME" or "SAME HIST"
  // histo2->Print("all");

BTW. When you post “source code” or “output” here, do remember to enclose them into two lines which contain just three characters ``` (see how your post has been edited above).

3

Thanks,
The GetBin seems to be workingThe SetBinContent doesn’t seems to be working.

For instance the bin content when i = 3 is 33.7865, and the bin error at i = 3 is 4.57392. histo2->SetBinContent(i, content / error); suppose to give a result of i = 3, and content/error = 8.3, but I was getting 0.015.

4

I guess it’s because you again mistake TProfile histograms for TH1D histograms (my example source code expects that the “hDiff” is a TH1D histogram).

5

How do you think I can fix this error?

6

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