Re: selection of ntuple and plot in a Graph

From: Rene Brun (Rene.Brun@cern.ch)
Date: Tue Apr 20 1999 - 16:41:31 MEST


Hi Ulrich.
The number of selectedRows is only computed by TTree::Draw, not by
GetEvent.
If you implement the selection logic yourself, there is no need to call
TTree::Draw. You should modify your program below to set the addresses
for the branches "sy" and "rep" and fill directly the TGraph object.
Below is a sketch of a possible macro.

Rene Brun

{
gROOT->Reset();

TFile* file = new TFile("data/test_repVI_17.4.1999_15:36:54.root",
"read");
TNtuple* ntuple = (TNtuple*) file->Get("ntuple");
TCanvas *canvas = new TCanvas("canvas", "ROOT Canvas", 500, 500);

static Float_t rep;
static Float_t sy;

ntuple->SetBranchAddress("sy",&sy);
ntuple->SetBranchAddress("rep",&rep);
TBranch *brep = ntuple->GetBranch("rep");

Int_t nevent = ntuple->GetEntries();
Int_t nselected = 0;
Float_t repold = -1;
Float_t *gsy  = new Float_t[nevent];
Float_t *grep = new Float_t[nevent];

for (Int_t i=0; i<nevent; i++) {
  brep->GetEvent(i);    // read branch 'rep' only
  if (r == repold)    // make some cuts, whatever
    continue;
  else
    repold = r;   //here put the brackets at teh right place!!!
  ntuple->GetEvent(i);  //read complete accepted event in memory. This
                        // seems not to mark the event as 'selected'.
  {
     gsy[nselected]  = sy;
     grep[nselected] = rep;
     nselected++;
  }
  
}
TGraph *graph = new TGraph(nselected, gsy, grep);
graph->Draw("ac*");

}

Ulrich Jost wrote:
> 
> Hi,
> 
> I try to make a selection in an ntuple and plot the surviving points in a
> graph, following some examples on the web. In the macro below the entire
> ntuple is plotted, ntuple->GetSelectedRows() gives 0. So how do I make the
> events 'selected'?
> We still have version 2.00/13.
> 
> Thanks, Ulrich
> 
> {
> gROOT->Reset();
> 
> TFile* file = new TFile("data/test_repVI_17.4.1999_15:36:54.root",
>                         "read");
> TNtuple* ntuple = (TNtuple*) file->Get("ntuple");
> TCanvas *canvas = new TCanvas("canvas", "ROOT Canvas", 500, 500);
> 
> static Float_t r;
> 
> TBranch *brep = ntuple->GetBranch("rep");
> brep->SetAddress(&r);
> 
> Int_t nevent = ntuple->GetEntries();
> Int_t nselected = 0;
> Float_t repold = -1;
> 
> for (Int_t i=0; i<nevent; i++) {
>   brep->GetEvent(i);    // read branch 'rep' only
>   if (r == repold)    // make some cuts, whatever
>     continue;
>   else
>     repold = r;
>   ntuple->GetEvent(i);  //read complete accepted event in memory. This
>                         // seems not to mark the event as 'selected'.
>   nselected++;
> }
> 
> cout<<"GetSelectedRows: "<<ntuple->GetSelectedRows()<<endl;
> ntuple->SetEstimate(ntuple->GetSelectedRows());
> ntuple->Draw("sy:rep", "", "goff");
> TGraph *graph = new TGraph(ntuple->GetSelectedRows(),
>                            ntuple->GetV2(), ntuple->GetV1());
> graph->Draw("ac*");
> 
> }



This archive was generated by hypermail 2b29 : Tue Jan 04 2000 - 00:43:32 MET