Filling a TH1F Conditionally in a "For" Loop

Newbie
1

Hello ROOTers,
NB: I’ve looked extensively through the forums for this topic and while there are many answers that hint at what I need to accomplish, there are still gaps. Please be gentle.

I have a TTree with branches “tof” (time of flight) and “energy”. I’d like to set up a nested for loop or some other script that iterates through varying energy ranges, makes a histogram of the tof data, does some analysis on the tof TH1F, and goes to the next energy range iteration.

void test_macro(string filename) 
{
    //open files to compare
	TFile *_file0 = TFile::Open(filename.c_str());


    // use nested for loop to iterate over energy ranges where i 
    // is the lower bound and j the upper
    
    for (double i = 0.; i < 1000.; i = i + 200.) {
      
      for (double j = 0.; j < 1000.; j = j + 200.) {
        
        if (j > i) {

            TTree *tree=(TTree*)_file0->GetObjectChecked("tree","TTree"); 
	        TH1D *NoTgt_tof=new TH1D("NoTgt_tof","EJ-315 ToF Gammas",200,8,200);

            TCanvas *c=new TCanvas();
	        c->cd(1);

            tree->Draw("tof>>NoTgt_tof", "energy > i && energy < j ","");


            // do analysis on the tof histogram             

        }

      }
      
    }

}

Some things that I’ve tried without success:

 tree->Draw("tof>>NoTgt_tof", Form("Energy > %d",i),"");

This would work except that I need two iterators (eventually I’ll need four!) and I can’t get the Form( workaround to work with both i and j.

I’ve also tried:

            TCut cut1 = Form("Energy > %d",i);
            TCut cut2 = Form("Energy < %d",j);

            tree->Draw("tof>>NoTgt_tof", cut1 && cut2,"");

This also doesn’t work because the i and j in the TCut never iterate.
I understand I can’t pass a local variable to the Draw() method, but I truly don’t need to draw the histograms to do the required analysis. So if I could either append the iterators to the TTree every loop, or if there were some analysis method that allowed me to conditionally fill a TH1F where I pass it local variables – that’d be swell.
Thanks for your help everyone.
_
Please read tips for efficient and successful posting and posting code

ROOT Version: 6.26/02
Platform: Windows
Compiler: VS

2

Hi Rahon,
if this line is changed

with the one below

tree->Draw("tof>>NoTgt_tof", Form("Energy > %d &&  Energy < %d",i,j));

You should be able to use both i and j iterator.

I add two suggestions

  • you could put the TTree declaration outside the for loop.
  • if you do not need to draw the histogram you could use the option "goff" after the cuts selection.

Best,
Stefano

3

Hello Stefano,
Thank you for the response!

I tried that line but it doesn’t appear that the values from the Form() statement are passing to the Draw() method. When I run the code the only output I see is that of the full-energy tof histogram (ie: without any cuts applied). My print statements tell me that i and j are iterating just fine but the ‘significance’ output tells me the Form() statement isn’t operating.

The full(er) code is below, I actually apply the cuts to two separate TH1Fs from two separate input files and then compare them.

As for the goff suggestion, is the below correct?

tree->Draw("tof>>NoTgt_tof", Form("Energy > %d &&  Energy < %d",i,j),"goff");

void peak_comparison(string filename) 
{
    //open files to compare
	TFile *_file0 = TFile::Open(filename.c_str());
    TFile *_file1 = TFile::Open("EJ315_081722_Ag_full_shielding.root");

    //use nested for loop to iterate over energy bounds  

    for (double i = 1.; i < 602.; i = i + 200.) {
      
      for (double j = 1.; j < 602.; j = j + 200.) {
        
        if (j > i) {

          TTree *tree=(TTree*)_file0->GetObjectChecked("tree","TTree"); 
	        TH1D *NoTgt_tof=new TH1D("NoTgt_tof","EJ-315 ToF Gammas",200,8,200);

          TCanvas *c=new TCanvas();
	        c->cd(1);

          tree->Draw("tof>>NoTgt_tof", Form("calibEnergy > %d && calibEnergy < %d",i,j),"");
          std::cout << "NoTgt applied from " << i << " to " << j << " keV:" << "\n";


          TTree *tree1=(TTree*)_file1->GetObjectChecked("tree","TTree");	
	        TH1D *Ag_tof=new TH1D("Ag_tof","ToF Ag",200,8,200);

	        //Draw onto same canvas
	        tree1->Draw("tof>>Ag_tof",Form("calibEnergy > %d && calibEnergy < %d",i,j),"same");
	

          //Scale histograms
          Ag_tof->Scale((NoTgt_tof->Integral(20,200))/(Ag_tof->Integral(20,200)));

          // determine t-test
          std::cout << "For the energy group of " << i << " to " << j << " keV:" << "\n";
	    
                int counts_NoTgt = NoTgt_tof->Integral(60,72);
	        int counts_Ag = Ag_tof->Integral(60,72);


          std::cout << "Significance: " << (counts_Ag - counts_NoTgt)/sqrt(counts_Ag+counts_NoTgt) << "\n";

            

        }

      }
      
    }

Thank you again for the help!

4

I just noticed that i and j are double. So in the Form you need %lf instead of %d.

Pay attention to this line

If the counts of the histograms are integers the division inside scale is rounded. maybe recast one of the integral as a float

 //Scale histograms
 Ag_tof->Scale(float(NoTgt_tof->Integral(20,200))/(Ag_tof->Integral(20,200)));

The goff is used correctly

5

Thank you! Yes, that fixed it. Good catch on the float!

7

Of course Integral returns a Double_t, but if the TH1 fill is done without any weight and no scale has been applied the result is an integer number.

8

@Dilicus could you please provide an example? I don’t see it (and would be extremely surprised if I did) with something super simple like

root [0] TH1F h("h", "my histogram", 5, 0., 5.);
root [1] h.Fill(2);
root [2] h.Integral()
(double) 1.0000000
10

Appologies for my quick answer above - what I changed was not the variable in the ‘Scale()’ line but rather the variable type of the iterators:

    for (float i = 50.; i < 600.; i = i + 200.) {
      
      for (float j = 50.; j < 600.; j = j + 200.) {
        
        if (j > i) {
	        
          _file0->cd();
          tree->Draw("tof>>NoTgt_tof", Form("calibEnergy > %lf && calibEnergy < %lf",i,j),"");
          std::cout << "NoTgt applied from " << i << " to " << j << " keV:" << "\n";

	        _file1->cd();
	        tree1->Draw("tof>>Ag_tof",Form("calibEnergy > %lf && calibEnergy < %lf",i,j),"same");
	        
          //Scale histograms
          Ag_tof->Scale((NoTgt_tof->Integral(20,200.))/(Ag_tof->Integral(20,200.)));

          // determine t-test
          std::cout << "For the energy group of " << i << " to " << j << " keV:" << "\n";
	        int counts_NoTgt = NoTgt_tof->Integral(60,72);
	        std::cout << "Counts in No Tgt " << counts_NoTgt << "\n";

	        int counts_Ag = Ag_tof->Integral(60,72);
	        std::cout << "Counts in Ag Tgt " << counts_Ag << "\n";

          std::cout << "Significance: " << (counts_Ag - counts_NoTgt)/sqrt(counts_Ag+counts_NoTgt) << "\n";
            
        }

      }
      
    }

And this gives the same results that I get when manually calculating energy bounds

11

Double_t counts_NoTgt
Double_t counts_Ag

12

@yus When I said integer number I was referring to the set of Natural number and not to the computer representation of an integer number (e.g int or Int_t )

@Wile_E_Coyote obviously I never seen a division between two double leading to an “integer rounding”
When I gave the advice regarding the cast inside the Scale, I gave it by heart and I was clearly wrong.
In the past, I used Scale(1/N) where N was an int and I had the wrong result because I should have used Scale(1./N), and for this reason when I see a Scale() I always thinking of a rounding error.
I should have double checked and have a look to the documentation of TH1::Integral().

I am sorry, I generated a lot of confusion.