Lines and integral

faca87 · March 6, 2023, 4:54pm

Hello, I’ve two problems with this macro

Here the root file

The Vertical and Horizontal lines are not drawn in the plot.

I get different area values (histo integral) depending on the binning.

NBin=1000 results

****************Background Estimation run0610 ****************
Time min= -5000
Time peak min= -2500
Time peak max= 2500
Time max= 5000
Time background external= 5000
Time peak= 5000
Total Area= 449491
Background Area Sx= 2197
Background Area Dx= 2149
Background Area= 4346
Peak Area= 445171
Peak Background Area= 4346
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.976254
True Events area= 440825

Nbins=100 results

****************Background Estimation run0610 ****************
Time min= -5000
Time peak min= -2500
Time peak max= 2500
Time max= 5000
Time background external= 5000
Time peak= 5000
Total Area= 449491
Background Area Sx= 2292
Background Area Dx= 2149
Background Area= 4441
Peak Area= 445258
Peak Background Area= 4441
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.997399
True Events area= 440817

dastudillo · March 7, 2023, 4:52am

For 1, draw the lines after the histogram.

couet · March 7, 2023, 5:58am

For 2
I see that the total is the same, which is good. It should not depend on the binning. Some partial areas are different, but given the bins’ repartition is different that’s not totally surprising.

Wile_E_Coyote · March 7, 2023, 7:54am

For 2

			// ...
			tmin = hgraph->GetXaxis()->GetBinLowEdge(BinMin);
			tmax = hgraph->GetXaxis()->GetBinUpEdge(BinMax);
			tpeakmin = hgraph->GetXaxis()->GetBinLowEdge(BinPeakMin);
			tpeakmax = hgraph->GetXaxis()->GetBinUpEdge(BinPeakMax);
			double TotArea= hgraph->Integral(BinMin,BinMax);
			double BkgAreaSx= hgraph->Integral(BinMin, BinPeakMin - 1);
			double BkgAreaDx= hgraph->Integral(BinPeakMax + 1, BinMax);
			// ...

faca87 · March 7, 2023, 8:48am

Hello, thank you to all

@dastudillo

I moved the draw lines after the histogram, I get the vertical lines, but I still don’t get the horizontal one.

@couet

Thank you, I adde the error hoping at least to trick the problem

Wile_E_Coyote:

For 2

			// ...
			tmin = hgraph->GetXaxis()->GetBinLowEdge(BinMin);
			tmax = hgraph->GetXaxis()->GetBinUpEdge(BinMax);
			tpeakmin = hgraph->GetXaxis()->GetBinLowEdge(BinPeakMin);
			tpeakmax = hgraph->GetXaxis()->GetBinUpEdge(BinPeakMax);
			double TotArea= hgraph->Integral(BinMin,BinMax);
			double BkgAreaSx= hgraph->Integral(BinMin, BinPeakMin - 1);
			double BkgAreaDx= hgraph->Integral(BinPeakMax + 1, BinMax);
			// ...

I followed your suggestion,

double BkgAreaSx= hgraph->IntegralAndError(BinMin,BinPeakMin-1,errorBkgAreaSx,""); 
			double BkgAreaDx= hgraph->IntegralAndError(BinPeakMax+1,BinMax,errorBkgAreaDx,"");

but I still get different partial areas:

Nbin=100

****************Background Estimation run0610 ****************
Time min= -5000 +/- 1
Time peak min= -2500 +/- 1
Time peak max= 2500 +/- 1
Time max= 5000 +/- 1
Time background external Sx= 2500 +/- 1.41421
Time background external Dx= 2500 +/- 1.41421
Time background external = 5000 +/- 2
Time peak= 5000 +/- 1.41421
Total Area= 449491 +/- 670.441
Background Area Sx= 2185 +/- 46.744
Background Area Dx= 2048 +/- 45.2548
Background Area= 4233 +/- 65.0615
Peak Area= 445258 +/- 667.277
Peak Background Area= 4233 +/- 65.0945
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.950685 +/- 1.54507
True Events area= 441025 +/- 670.444

Nbin=1000

****************Background Estimation run0610 ****************
Time min= -5000 +/- 1
Time peak min= -2500 +/- 1
Time peak max= 2500 +/- 1
Time max= 5000 +/- 1
Time background external Sx= 2500 +/- 1.41421
Time background external Dx= 2500 +/- 1.41421
Time background external = 5000 +/- 2
Time peak= 5000 +/- 1.41421
Total Area= 449491 +/- 670.441
Background Area Sx= 2185 +/- 46.744
Background Area Dx= 2135 +/- 46.2061
Background Area= 4320 +/- 65.7267
Peak Area= 445171 +/- 667.211
Peak Background Area= 4320 +/- 65.7608
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.970414 +/- 1.5296
True Events area= 440851 +/- 670.444

@all

Moreover,

I added the TPave. Even if I set

	t->SetLineColor(0);
 			t->SetFillStyle(0);

I don’t get transparent shadow, as you can see in the plot.

I would like to highlight the background by drawing some boxes. Then I used

TBox *box1 = new TBox(0.5,0.5,0.8,0.8);
   			box1->SetFillColor(19);
   			box1->Draw();

but I don’t get the box.

To better understand…

Currently I get this plot

but , I should get something like this one

deltatbkgest.cpp (8.9 KB)

Wile_E_Coyote · March 7, 2023, 9:03am

No, you DID NOT. You only took part of it.

You now introduced another problem. You should have:
errort = hgraph->GetXaxis()->GetBinWidth(1); // or divide it by 2.

faca87 · March 7, 2023, 9:13am

Hello,

I did it, but still getting different partial areas

NBin=1000

****************Background Estimation run0610 ****************
Time min= -5000 +/- 10
Time peak min= -2500 +/- 10
Time peak max= 2500 +/- 10
Time max= 5000 +/- 10
Time background external Sx= 2500 +/- 14.1421
Time background external Dx= 2500 +/- 14.1421
Time background external = 5000 +/- 20
Time peak= 5000 +/- 14.1421
Total Area= 449491 +/- 670.441
Background Area Sx= 2185 +/- 46.744
Background Area Dx= 2135 +/- 46.2061
Background Area= 4320 +/- 65.7267
Peak Area= 445171 +/- 667.211
Peak Background Area= 4320 +/- 69.05
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.970414 +/- 1.60539
True Events area= 440851 +/- 670.775

NBin=100

****************Background Estimation run0610 ****************
Time min= -5000 +/- 100
Time peak min= -2500 +/- 100
Time peak max= 2500 +/- 100
Time max= 5000 +/- 100
Time background external Sx= 2500 +/- 141.421
Time background external Dx= 2500 +/- 141.421
Time background external = 5000 +/- 200
Time peak= 5000 +/- 141.421
Total Area= 449491 +/- 670.441
Background Area Sx= 2185 +/- 46.744
Background Area Dx= 2048 +/- 45.2548
Background Area= 4233 +/- 65.0615
Peak Area= 445258 +/- 667.277
Peak Background Area= 4233 +/- 217.341
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.950685 +/- 5.13662
True Events area= 441025 +/- 701.78

deltatbkgest.cpp (9.0 KB)

Wile_E_Coyote · March 7, 2023, 9:15am

Apply the WHOLE set of fixes:

couet · March 7, 2023, 9:22am

Regarding graphics:

double tmin=-5000;
double tpeakmin=-2500;
double tpeakmax=2500;
double tmax=5000;

void deltatbkgest()
{

   float fillcolor=0.35;
   gStyle->SetFitFormat("8.6g");
   gStyle->SetOptStat(111110);
   ofstream results;
   results.open("ResultsBackgroundRun0610.txt",  ios::out);
   TFile *fin = TFile::Open("run0610.root");
   fin->ls ();
   TGaxis::SetMaxDigits(3);
   if (fin == 0) {
      printf("Error: cannot open the file!\n");
   } else {
      TTree *te=0;
      fin->GetObject("tree",te);
      gStyle->SetOptFit();
      gStyle->SetOptStat(111110);
      te->SetLineColor(kBlue);
      TH1::StatOverflows(true);
      TCanvas *c01 = new TCanvas("c01","c01",1280,1024);
      c01->SetLogy();
      TString string = TString::Format("t1-t0 >> htemp(1000, -5000, 5000)");
      TCut coincidence = TString::Format("abs(t1-t0)<5000").Data();
      TCut pulser0 = TString::Format("e2==0").Data();
      TCut pulser1 = TString::Format("e3==0").Data();
      TCut GePD = TString::Format("e0>0").Data();
      TCut GeDD = TString::Format("e1>0").Data();
      te->Draw(string,coincidence && GePD && GeDD && pulser0 && pulser1);
      te->GetHistogram()->SetTitle("");
      te->SetScanField(0);
      TH1F *hgraph= (TH1F*)gPad->GetPrimitive("htemp");
      hgraph->GetXaxis()->SetTitle("#DeltaT (ns)");
      hgraph->GetYaxis()->SetTitle("Counts");
      hgraph->GetYaxis()->SetTitleSize(40);
      hgraph->GetYaxis()->SetTitleFont(43);
      hgraph->GetYaxis()->SetTitleOffset(1.1);
      hgraph->GetYaxis()->SetLabelFont(43);
      hgraph->GetYaxis()->SetLabelSize(40);
      hgraph->GetXaxis()->SetTitleSize(40);
      hgraph->GetXaxis()->SetTitleFont(43);
      hgraph->GetXaxis()->SetTitleOffset(1.1);
      hgraph->GetXaxis()->SetLabelFont(43);
      hgraph->GetXaxis()->SetLabelSize(40);
      hgraph->SetFillStyle(3005);
      hgraph->SetFillColorAlpha(kBlue, fillcolor);


      hgraph->SetStats(0);
      hgraph->Draw();
      hgraph->SetName("");

      auto lh = new TLine(tpeakmin,10,tpeakmax,10);
      lh->SetLineColor(kRed);
      lh->SetLineWidth(3.);
      lh->Draw();

      auto lv1 = new TLine(tpeakmin,1,tpeakmin,10000);
      lv1->SetLineColor(kRed);
      lv1->SetLineWidth(3.);
      lv1->Draw();

      auto lv2 = new TLine(tpeakmax,1,tpeakmax,10000);
      lv2->SetLineColor(kRed);
      lv2->SetLineWidth(3.);
      lv2->Draw();

      auto t1 = new TText(0,5,"Random coinc");
      t1->SetTextAlign(23);
      t1->SetTextColor(kRed);
      t1->SetTextFont(132);
      t1->Draw();

      auto t2 = new TText(0,100,"True coinc");
      t2->SetTextAlign(23);
      t2->SetTextColor(kCyan);
      t2->SetTextFont(132);
      t2->SetTextAngle(45);
      t2->Draw();

      int BinMin = hgraph->GetXaxis()-> FindBin(tmin);
      int BinMax = hgraph->GetXaxis()-> FindBin(tmax);
      int BinPeakMin = hgraph->GetXaxis()-> FindBin(tpeakmin);
      int BinPeakMax = hgraph->GetXaxis()-> FindBin(tpeakmax);
      double TotArea= hgraph->Integral(BinMin,BinMax);
      double BkgAreaSx= hgraph->Integral(BinMin,BinPeakMin);
      double BkgAreaDx= hgraph->Integral(BinPeakMax,BinMax);
      double BkgArea= BkgAreaSx+BkgAreaDx;
      double PeakArea= hgraph->Integral(BinPeakMin,BinPeakMax);
      double texternal=abs(tpeakmin-tmin)+abs(tmax-tpeakmax);
      double tpeak=abs(tpeakmax-tpeakmin);
      double PeakBkgArea= (BkgArea/texternal)*tpeak;
      double PeakBkgPercentage= 100*PeakBkgArea/PeakArea;
      double TrueEventsArea= PeakArea-PeakBkgArea;
      results << "****************Background Estimation run0610 ****************"<<  endl;
      results << "Time min= " << tmin <<  endl;
      results << "Time peak min= " << tpeakmin <<  endl;
      results << "Time peak max= " << tpeakmax <<  endl;
      results << "Time max= " << tmax <<  endl;
      results << "Time background external= " << texternal <<  endl;
      results << "Time peak= " << tpeak <<  endl;
      results << "Total Area= " << TotArea <<  endl;
      results << "Background Area Sx= " << BkgAreaSx <<  endl;
      results << "Background Area Dx= " << BkgAreaDx <<  endl;
      results << "Background Area= " << BkgArea <<  endl;
      results << "Peak Area= " << PeakArea <<  endl;
      results << "Peak Background Area= " << PeakBkgArea<<  endl;
      results << "Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  " << PeakBkgPercentage <<  endl;
      results << "True Events area= " << TrueEventsArea <<  endl;
   }
}

faca87 · March 7, 2023, 9:50am

Hello

@wile_E_Coyote

Now I think I I apply the whole set of fixes

int BinMin = hgraph->GetXaxis()-> FindBin(tmin);
   			int BinMax = hgraph->GetXaxis()-> FindBin(tmax);
			int BinPeakMin = hgraph->GetXaxis()-> FindBin(tpeakmin);
   			int BinPeakMax = hgraph->GetXaxis()-> FindBin(tpeakmax);
   			errort=hgraph->GetXaxis()->GetBinWidth(1);
   			tmin = hgraph->GetXaxis()->GetBinLowEdge(BinMin);
			tmax = hgraph->GetXaxis()->GetBinUpEdge(BinMax);
			tpeakmin = hgraph->GetXaxis()->GetBinLowEdge(BinPeakMin);
			tpeakmax = hgraph->GetXaxis()->GetBinUpEdge(BinPeakMax);
			double TotArea= hgraph->IntegralAndError(BinMin,BinMax,errorTotArea,""); 
			double BkgAreaSx= hgraph->IntegralAndError(BinMin,BinPeakMin-1,errorBkgAreaSx,""); 
			double BkgAreaDx= hgraph->IntegralAndError(BinPeakMax+1,BinMax,errorBkgAreaDx,"");

but I still get different partial areas

NBin=1000

****************Background Estimation run0610 ****************
Time min= -5000 +/- 10
Time peak min= -2500 +/- 10
Time peak max= 2510 +/- 10
Time max= 5010 +/- 10
Time background external Sx= 2500 +/- 14.1421
Time background external Dx= 2500 +/- 14.1421
Time background external = 5000 +/- 20
Time peak= 5010 +/- 14.1421
Total Area= 449491 +/- 670.441
Background Area Sx= 2185 +/- 46.744
Background Area Dx= 2135 +/- 46.2061
Background Area= 4320 +/- 65.7267
Peak Area= 445171 +/- 667.211
Peak Background Area= 4328.64 +/- 69.1837
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.972354 +/- 1.60529
True Events area= 440842 +/- 670.789

NBin=100

****************Background Estimation run0610 ****************
Time min= -5000 +/- 100
Time peak min= -2500 +/- 100
Time peak max= 2600 +/- 100
Time max= 5100 +/- 100
Time background external Sx= 2500 +/- 141.421
Time background external Dx= 2500 +/- 141.421
Time background external = 5000 +/- 200
Time peak= 5100 +/- 141.421
Total Area= 449491 +/- 670.441
Background Area Sx= 2185 +/- 46.744
Background Area Dx= 2048 +/- 45.2548
Background Area= 4233 +/- 65.0615
Peak Area= 445258 +/- 667.277
Peak Background Area= 4317.66 +/- 220.377
Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  0.969698 +/- 5.10629
True Events area= 440940 +/- 702.726

deltatbkgest.cpp (9.2 KB)

@couet

I see in your plot (as well in the example that I uploaded) that the vertical lines stop before than the end of the canvas. I know it is a well known problem wehn using log scale…that’s because I tried to use the VerticalLine and HorizontalLine functions instead of the TLine (I found the code for VerticalLine in your old post and I tried to adapt it also ffor the Horizontal one). But it isnt’ a big problem…I also can use the TLine…but I see it gives problems depending on the binning. Indeed, if i set nbin=100, I get

Moreover, do you know why I couldn’t get the boxes? (i.e. drawing the green and orange box as in the example)?

Thanks
deltatbkgest.cpp (9.4 KB)

Wile_E_Coyote · March 7, 2023, 9:54am

BTW. Try with:

double tpeakmax = 2499.999;
double tmax = 4999.999;

couet · March 7, 2023, 10:21am

#include "TFile.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TH1.h"
#include "TH2.h"
#include "TH3.h"
#include "TGaxis.h"
#include "TRandom.h"
#include "TLegend.h"
#include "TPaveStats.h"
#include "TGraph.h"
#include "TSystem.h"
#include "TTree.h"
#include "TTreePlayer.h"
#include "TF1.h"
#include "TCut.h"
#include "TPaletteAxis.h"
#include "TFitResult.h"
#include "TMath.h"

double tmin=-5000;
double tpeakmin=-2500;
double tpeakmax=2500;
double tmax=5000;

void deltatbkgest()
{

   float fillcolor=0.35;
   gStyle->SetFitFormat("8.6g");
   gStyle->SetOptStat(111110);
   ofstream results;
   results.open("ResultsBackgroundRun0610.txt",  ios::out);
   TFile *fin = TFile::Open("run0610.root");
   fin->ls ();
   TGaxis::SetMaxDigits(3);
   if (fin == 0) {
      printf("Error: cannot open the file!\n");
   } else {
      TTree *te=0;
      fin->GetObject("tree",te);
      gStyle->SetOptFit();
      gStyle->SetOptStat(111110);
      te->SetLineColor(kBlue);
      TH1::StatOverflows(true);
      TCanvas *c01 = new TCanvas("c01","c01",1280,1024);
      c01->SetLogy();
      TString string = TString::Format("t1-t0 >> htemp(1000, -5000, 5000)");
      TCut coincidence = TString::Format("abs(t1-t0)<5000").Data();
      TCut pulser0 = TString::Format("e2==0").Data();
      TCut pulser1 = TString::Format("e3==0").Data();
      TCut GePD = TString::Format("e0>0").Data();
      TCut GeDD = TString::Format("e1>0").Data();
      te->Draw(string,coincidence && GePD && GeDD && pulser0 && pulser1);
      te->GetHistogram()->SetTitle("");
      te->SetScanField(0);
      TH1F *hgraph= (TH1F*)gPad->GetPrimitive("htemp");
      hgraph->GetXaxis()->SetTitle("#DeltaT (ns)");
      hgraph->GetYaxis()->SetTitle("Counts");
      hgraph->GetYaxis()->SetTitleSize(40);
      hgraph->GetYaxis()->SetTitleFont(43);
      hgraph->GetYaxis()->SetTitleOffset(1.1);
      hgraph->GetYaxis()->SetLabelFont(43);
      hgraph->GetYaxis()->SetLabelSize(40);
      hgraph->GetXaxis()->SetTitleSize(40);
      hgraph->GetXaxis()->SetTitleFont(43);
      hgraph->GetXaxis()->SetTitleOffset(1.1);
      hgraph->GetXaxis()->SetLabelFont(43);
      hgraph->GetXaxis()->SetLabelSize(40);
      hgraph->SetFillStyle(3005);
      hgraph->SetFillColor(kBlue);


      hgraph->SetStats(0);
      hgraph->Draw();
      hgraph->SetName("");

      auto b1 = new TBox(-5000,0,5000,10);
      b1->SetFillColorAlpha(kGreen-4,0.5);
      b1->Draw();

      auto b2 = new TBox(tpeakmin,0,tpeakmax,10);
      b2->SetFillColorAlpha(kOrange,0.5);
      b2->Draw();

      auto lh = new TLine(tpeakmin,10,tpeakmax,10);
      lh->SetLineColor(kRed);
      lh->SetLineWidth(3.);
      lh->Draw();

      auto lv1 = new TLine(tpeakmin,1,tpeakmin,13765);
      lv1->SetLineColor(kRed);
      lv1->SetLineWidth(3.);
      lv1->Draw();

      auto lv2 = new TLine(tpeakmax,1,tpeakmax,13765);
      lv2->SetLineColor(kRed);
      lv2->SetLineWidth(3.);
      lv2->Draw();

      auto t1 = new TText(0,5,"Random coinc");
      t1->SetTextAlign(23);
      t1->SetTextColor(kRed);
      t1->SetTextFont(132);
      t1->Draw();

      auto t2 = new TText(0,100,"True coinc");
      t2->SetTextAlign(23);
      t2->SetTextColor(kCyan);
      t2->SetTextFont(132);
      t2->SetTextAngle(45);
      t2->Draw();

      int BinMin = hgraph->GetXaxis()-> FindBin(tmin);
      int BinMax = hgraph->GetXaxis()-> FindBin(tmax);
      int BinPeakMin = hgraph->GetXaxis()-> FindBin(tpeakmin);
      int BinPeakMax = hgraph->GetXaxis()-> FindBin(tpeakmax);
      double TotArea= hgraph->Integral(BinMin,BinMax);
      double BkgAreaSx= hgraph->Integral(BinMin,BinPeakMin);
      double BkgAreaDx= hgraph->Integral(BinPeakMax,BinMax);
      double BkgArea= BkgAreaSx+BkgAreaDx;
      double PeakArea= hgraph->Integral(BinPeakMin,BinPeakMax);
      double texternal=abs(tpeakmin-tmin)+abs(tmax-tpeakmax);
      double tpeak=abs(tpeakmax-tpeakmin);
      double PeakBkgArea= (BkgArea/texternal)*tpeak;
      double PeakBkgPercentage= 100*PeakBkgArea/PeakArea;
      double TrueEventsArea= PeakArea-PeakBkgArea;
      results << "****************Background Estimation run0610 ****************"<<  endl;
      results << "Time min= " << tmin <<  endl;
      results << "Time peak min= " << tpeakmin <<  endl;
      results << "Time peak max= " << tpeakmax <<  endl;
      results << "Time max= " << tmax <<  endl;
      results << "Time background external= " << texternal <<  endl;
      results << "Time peak= " << tpeak <<  endl;
      results << "Total Area= " << TotArea <<  endl;
      results << "Background Area Sx= " << BkgAreaSx <<  endl;
      results << "Background Area Dx= " << BkgAreaDx <<  endl;
      results << "Background Area= " << BkgArea <<  endl;
      results << "Peak Area= " << PeakArea <<  endl;
      results << "Peak Background Area= " << PeakBkgArea<<  endl;
      results << "Peak Background Area (%)= 100*PeakBkgArea/PeakArea=  " << PeakBkgPercentage <<  endl;
      results << "True Events area= " << TrueEventsArea <<  endl;
   }
}

faca87 · March 7, 2023, 1:09pm

Hello

Yes I read @couet’s message…by you reply it looke like that maybe there was a solution…but don’t worry, it isn’t so important (I hope).

@couet

Thank you, I got it

system · March 21, 2023, 1:10pm

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