You are here

How to Read a Tree?

Read a TTree via the Draw method

Use this function if you want to quickly histogram a few variables with a very simple selection mechanism. The loop on the Tree is automatically done by ROOT. The signature of the TTree::Draw method is:

void Draw(const char *varexp, const char *selection, 
          Option_t *option,Int_t nevents, Int_t firstevent);

where the parameters have the following meaning:

  • varexp is an expression with the forms e1, e1:e2 or e1:e2:e3 to generate 1-D, 2-D or 3-D distributions respectively and where e1,e2,e3 are expressions with combinations of the original Tree variables. Examples:
varexp = x           : draw a 1-Dim distribution of column named x
            = sqrt(x)   : draw distribution of sqrt(x)
            = x*y/z
            = y:sqrt(x) : 2-Dim distribution of y versus sqrt(x)
  • selection is an expression with a combination of the Tree variables. In a selection all the C++ operators are authorized. Example:
      selection = "x<y && sqrt(z)>3.2"
  • option is the drawing option. See TH1::Draw for the list of all drawing options. If option contains the string "goff", no graphics is generated.
  • nevents is the number of events to process (default is all).
  • firstevent is the first event to process (default is 0).

By default the temporary histogram created is called htemp. If varexp contains >>hnew, the new histogram created is called hnew and it is kept in the current directory. Example:

    tree.Draw("sqrt(x)>>hsqrt","y>0");

will draw sqrt(x) and save the histogram as "hsqrt" in the current directory. By default, the specified histogram is reset. To continue to append data to an existing histogram, use "+" in front of the histogram name.

     tree.Draw("sqrt(x)>>+hsqrt","y>0");

will not reset hsqrt, but will continue filling.

Look at the tutorial macro ntuple1 for a more detailed example. The kind of expressions you can specify in this function is restricted on purpose (see class TFormula). Note that this method does not require access to the original class code, if you only want to perform the kind of operations described above.

Read a TTree via the TTreeReader

In case the complexity of the analysis which has to be performed on the data contained in the Tree cannot be expressed using the TTree::Draw method, other ways of looping on the Tree can be adopted. For example, this operation can be accomplished with a TTreeReader. In the following example, we'll see how:

  • to loop on a tree.
  • to read a sub-branch only (the branch containing the number of tracks).
  • to read all the branches, when the selection is positive.
  • to histogram the result of a member function of the event object.
TH1F* hnseg(nullptr);
void readEvent(){

   hnseg = new TH1F("hnseg",
                    "Number of segments for selected tracks",
                    4096,0,8192);
   //   Connect file generated in $ROOTSYS/test
   TFile fileIn("Event.root");
   TTree* theTree = nullptr;    
   TTreeReader theReader("T",&fileIn);
   TTreeReaderValue<Event> eventRV(theReader, "event");
   TTreeReaderValue<Int_t> nTracksRV(theReader, "fNtrack");

   while(theReader.Next()){
      if (*nTracksRV < 587) {
         continue; // Check if we don't have too many tracks
      }
      auto event = eventRV.Get();      //Read complete accepted event
                                       //in memory.
      hnseg->Fill(event->GetNseg());   //Fill histogram with number of
                                       //segments.
   }
   hnseg->Draw();
}

It is important to note that not the whole even is read immediately, but the complete deserialization happens only if a certain condition is met.

Read a TTree Manipulating Branches Explicitely

An equivalent result can be achieved manipulating branches explicitly.

TH1F* hnseg(nullptr);
void readEvent2(){
   hnseg = new TH1F("hnseg",
                    "Number of segments for selected tracks",
                    4096,0,8192);

   TFile fileIn("Event.root");

   TTree* theTree = nullptr;
   fileIn.GetObject("T",theTree);

   TClonesArray *tracks = nullptr;
   Event *event = new Event(); //object must be created before
                               //setting the branch address

   auto bntrack = theTree->GetBranch("fNtrack");
   auto branch  = theTree->GetBranch("event");
   branch->SetAddress(&event);
   auto nevent = theTree->GetEntries();
   for (Int_t i=0;i<nevent;i++) {
      bntrack->GetEvent(i);
      if (event->GetNtrack() < 587)continue;
      theTree->GetEvent(i);           //Read complete accepted event
                                      //in memory.
      hnseg->Fill(event->GetNseg());  //Fill histogram with number of
                                      //segments.
      tracks = event->GetTracks();    //Get pointer to the
                                      //TClonesArray object.
      tracks->Clear();                //Clear it.
   }

   hnseg->Draw();
}