I need looping in a root macro

Dear Experts,

I have written a macro in which I am plotting variables from root file with a cut Kpt>0.2. I want to plot the same variable for different Kpt cut values. Kpt values should increase at the interval of 0.05 upto 2.So I have to add a for loop. Could you please tell me how to do it?

Thanks,
Deepak

Hey,

I am usually using a few simple lines of code that are not very efficient, but hey, it works.

TString cutvalue = "";

for (float cut = 0.2; cut<=2.0;cut+=0.05)
{
cutvalue +=cut;
tree->Draw("plot_variable>>HistogramKptcut"+cutvalue,"Kpt>"+cutvalue); 
cutvalue="";
}

I assume you are plotting some variables stored in TTree and that “Kpt>0.2” is the selection string in TTree::Draw. As this selection parameter is a text string you need to do something like:

 for (....) {
    val = .....;
    t->Draw ("...",Form("Kpt > %g", val), .....);
}
1 Like

Hi mwojtas,

Could you please explain “plot_variable>>HistogramKptcut”+cutvalue ?

Thanks,
Deepak

@mwojtas concatenates a variable text string to the selection string to have a variable cut. It works fine but I would suggest you try the method using Form I explained in my previous post. It is more direct as you do not need a temporary text string and more flexible as you can put the varying value anywhere in the cut expression.

Hello,
As @couet said, the form() method is better. I just got used to using the concatenation. The “plot_variable>>HistogramKptcut”+cutvalue
Plots a plot_variable and saves it to a TH1 named “HistogramKptcut”+cutvalue, so for example the last histogram saved would be “HistogramKptcut2.0”. For more details on how to save drawn histograms see
https://root.cern.ch/doc/master/classTTree.html#a73450649dc6e54b5b94516c468523e45

Cheers,
Maks

Hi ,

I am pasting my macro. Please have a look.

{

  TChain *ch = new TChain("tree");
  ch->Add("sel_BsToPhiMuMu_mumuphi_mc_s*.root");
  TTree *tr = ch;
  
  tr->SetBranchStatus("*", kFALSE);
  tr->SetBranchStatus("Bmass", kTRUE);
  tr->SetBranchStatus("Kppt", kTRUE);      
  tr->SetBranchStatus("Kmpt", kTRUE);  
  tr->SetBranchStatus("Triggers", kTRUE);

  TH1D *h_bm   = new TH1D("h_bm"  ,"B_{s} #rightarrow  #mu#mu#phi; M(B_{s}) [GeV/c^{2}] ; Events/0.014GeV", 50, 5.0, 5.7);
  h_bm->GetYaxis()->SetTitleOffset(1.6);

  TCut trig = "Triggers ==1";

  /*for(double i=0.2; i<0.5;i+=0.05)

    {
      TCut kp = Form("Kppt>%d",i);
      //	cout<<kp<<endl;
      TCut km = Form("Kmpt>%d",i);
  */  
      
      TCut selbm = trig && kp && km;
      
      //cout<<i<<endl;

  
      /*
  TString cutvalue = "";

  for (float cut = 0.2; cut<=2.0;cut+=0.05)
    {
      cutvalue +=cut;
     
  */   
      TCanvas *cv1 = new TCanvas("cv1","",600,600);                                                                                                
       tr->Draw("Bmass >> h_bm",selbm);
      //  tr->Draw("Bmass >> h_bm","Kppt>"+cutvalue);
      //cout<<cutvalue<<endl;    
      //  cutvalue="";
   
      double Sraw, S;  
      double scalefactor= 864.467005076;
      
      TAxis *axis = h_bm->GetXaxis();
      int bmin = axis->FindBin(5.287);
      int bmax = axis->FindBin(5.447); 
      double integral = h_bm->Integral(bmin,bmax);
      Sraw = integral;
      S= Sraw/scalefactor;
      cout<<"Sraw="<<Sraw<<endl;
      cout<<"S="<<S<<endl;
      cv1->SaveAs("BsToPhiMuMumc2_5sigma*i.png");
}

It’s not working .

Hi,

After a brief look I can tell one main mistake, the h_bm histogram is not retrieved in a way you described. The correct way to retrieve it is:

tree.Draw("sqrt(x)>>hsqrt","y>0")
TH1F *hsqrt = (TH1F*)gDirectory->Get("hsqrt");

All of that can be found in the link I send you before. What I am not sure about is whether or not you need to have an TFile open do write to store the histograms, so you could try adding

 TFile * MockFile = new TFile("mock.root ", "recreate");

in the very beginning.

Cheers,
Maks

Your macro is almost un-readable … The loops are commented… what should we keep ? the syntax is incorrect etc …

can you provide "the* real “non working macro” with a few explanations saying why it does not work.

The best would be that you provide something we could run i.e: macro + data file

Hi Couet,

I am pasting the non working macro. The loop is running but it is not changing the Kppt and Kmpt values.

{

  TChain *ch = new TChain("tree");
  ch->Add("sel_BsToPhiMuMu_mumuphi_mc_s0.root");
  TTree *tr = ch;
  
  tr->SetBranchStatus("*", kFALSE);
  tr->SetBranchStatus("Bmass", kTRUE);
  tr->SetBranchStatus("Kppt", kTRUE);      
  tr->SetBranchStatus("Kmpt", kTRUE);  
  tr->SetBranchStatus("Triggers", kTRUE);

  TH1D *h_bm   = new TH1D("h_bm"  ,"B_{s} #rightarrow  #mu#mu#phi; M(B_{s}) [GeV/c^{2}] ; Events/0.014GeV", 50, 5.0, 5.7);
  h_bm->GetYaxis()->SetTitleOffset(1.6);

  //TCut trig = "Triggers ==1";

  for(double i=0.2; i<0.5;i+=0.05)

    {
      TCut trig = "Triggers ==1";
      
      TCut kp = Form("Kppt >%d",i);
      TCut km = Form("Kmpt > %d",i);
      
      TCut selbm = trig && kp && km;
      
      //cout<<i<<endl;

  
      /*
	TString cutvalue = "";

	for (float cut = 0.2; cut<=2.0;cut+=0.05)
	{
	cutvalue +=cut;
	
      */
   
      TCanvas *cv1 = new TCanvas "cv1","",600,600);                                                                                                
      tr->Draw("Bmass >> h_bm",selbm);
      //  tr->Draw("Bmass >> h_bm","Kppt>"+cutvalue);
      //cout<<cutvalue<<endl;    
      //  cutvalue="";
   
      double Sraw, S;  
      double scalefactor= 864.467005076;
      
      TAxis *axis = h_bm->GetXaxis();
      int bmin = axis->FindBin(5.287);
      int bmax = axis->FindBin(5.447); 
      double integral = h_bm->Integral(bmin,bmax);
      Sraw = integral;
      S= Sraw/scalefactor;
      cout<<"Sraw="<<Sraw<<endl;
      cout<<"S="<<S<<endl;
      cv1->SaveAs("BsToPhiMuMumc2_5sigma*i.png");
}
   
}

Thanks,
Deepak

i is a double not an int. do:

{
   for (double i=0.2; i<0.5;i+=0.05) {
      TCut trig = "Triggers ==1";
      TCut kp = Form("Kppt >%g",i); // %g !... not %d !...
      TCut km = Form("Kmpt > %g",i);
      TCut selbm = trig && kp && km;
      printf("%s\n",selbm.GetTitle());
   }
}

Hi,

Here’s what I’d do: create a 2 dimensional histogram, with Kpt as second dimension! You can always create slices using ProjectX(). This is much more efficient than re-running over all tree entries for every cut value!

Axel.

1 Like

Hi Couet,

Thanks,It’s working. But I have added another root file and similarly I want to calculate. But it is not working properly. I am attaching the macro .Please have a look on it. It is giving the same “S” value for all cuts where as different values for B .

{
  TChain *ch = new TChain("tree");
  ch->Add("sel_BsToPhiMuMu_mumuphi_mc_s0.root");
  TTree *tr = ch;
  
  tr->SetBranchStatus("*", kFALSE);
  tr->SetBranchStatus("Bmass", kTRUE);
  tr->SetBranchStatus("Kppt", kTRUE);      
  tr->SetBranchStatus("Kmpt", kTRUE);  
  tr->SetBranchStatus("Triggers", kTRUE);
  
    TH1D *h_bms   = new TH1D("h_bms"  ,"B_{s} #rightarrow  #mu#mu#phi; M(B_{s}) [GeV/c^{2}] ; Events/0.014GeV", 50, 5.0, 5.7);
  h_bms->GetYaxis()->SetTitleOffset(1.6);

  for (double i=0.2; i<2.0;i+=0.05) {
    TCut trig = "Triggers ==1";
    TCut kp = Form("Kppt >%g",i); // %g !... not %d !...
    TCut km = Form("Kmpt > %g",i);
    TCut selbm = trig && kp && km;
    printf("%s\n",selbm.GetTitle());
    
     
    TCanvas *cv1 = new TCanvas("cv1","",600,600);
    tr->Draw("Bmass >> h_bms",selbm);

    double Sraw, S;  
    double scalefactor= 864.467005076;
      
    TAxis *axis = h_bms->GetXaxis();
    int bmin = axis->FindBin(5.287);
    int bmax = axis->FindBin(5.447); 
    double integral = h_bms->Integral(bmin,bmax);
    Sraw = integral;
    S= Sraw/scalefactor;
    cout<<"Sraw="<<Sraw<<endl;
    cout<<"S="<<S<<endl;
    cv1->SaveAs("BsToPhiMuMumc2_5sigma%g,i.png");
     
   
    TFile *f = new TFile("sel_BsToPhiMuMu_2016BCDEF_data_s0.root");
    TTree *t1  =(TTree*)f->Get("tree");

    TH1D *h_bmb   = new TH1D("h_bmb"  ,"B_{s} #rightarrow  #mu#mu#phi Data; M(B_{s}) [GeV/c^{2}] ; Events/0.001", 700, 5.0, 5.7);

    TCut bm = "Bmass < 5.223 || Bmass>5.511";
    TCut selb = trig && kp && km && bm ;
    printf("%s\n",selb.GetTitle());
    
    TCanvas *cv2 = new TCanvas("cv2","",600,600);
    t1->Draw("Bmass >> h_bmb",selb);
    h_bmb->GetYaxis()->SetTitleOffset(1.6);
    
    double N1,N2,B;

    TAxis *axis = h_bmb->GetXaxis();
    int bmin1 = axis->FindBin(5.143); //in your case xmin=-1.5
    int bmax1 = axis->FindBin(5.223); 
    double integral1 = h_bmb->Integral(bmin1,bmax1);

    N1=integral1;
    cout<<"lower sideband entries"<<N1<<endl;
    
    int bmin2 = axis->FindBin(5.511); //in your case xmin=-1.5                   
    int bmax2 = axis->FindBin(5.591);
    double integral2 = h_bmb->Integral(bmin2,bmax2);
      
    N2=integral2;
    cout<<"upper sideband entries"<<N2<<endl;
    
    B= N1+N2;
    cout<<"B="<<B<<endl;
    
     
    }

}

To run your macro we are missing the file sel_BsToPhiMuMu_mumuphi_mc_s0.root

Hi Couet,

Here is the file. Please replace the second file as same as first file.

Sorry, I tried to upload the file but it is showing error in uploading for several times.

Thanks,
Deepak

If you have a CERN account, you can put it into your public folder there and post the path here.

Dear Couet,

The loop is working fine. Now I want to save only S and B values in a text file. Could you please tell me how to do it?

Thanks,
Deepak

That’s standard C++. You can easily file tutorials on the web explaining this. For instance this one.

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