Re: [ROOT] How to write lists of objects?

From: Rene Brun (Rene.Brun@cern.ch)
Date: Mon Nov 12 2001 - 16:52:56 MET


Hi Daniel,

On request from H1, I implemented about one year ago a special branch
constructor
  TTree::Branch(TList*list, Int_t bufsize, Int_t splitlevel)

I could extend this function in the way shown below. Note that this will also
require a special new function to set the corresponding branch address
when reading a Tree created this way.
Let me know what you think of this proposal

Rene Brun


Int_t TTree::Branch(TCollection *list, Int_t bufsize, Int_t splitlevel)
{
//   This function creates one branch for each element in the collection.
//   Each entry in the collection becomes a top level branch if the 
//   corresponding class is not a collection. If it is a collection, the entry
//   in the collection becomes in turn top level branches, etc.
//   The splitlevel is decreased by 1 everytime a new collection is found.
//   For example if list is a TObjArray* 
//     - if splitlevel = 0 only one top level branch is created
//     - if splitlevel = 1, one top level branch is created for each element
//        of the TObjArray.
//     - if splitlevel = 2, one top level branch is created for each array
element.
//       if, in turn, one of the array elements is a TCollection, one top level
//       branch will be created for each element of this collection.
//
//   In case a collection element is a TClonesArray, the special Tree
constructor
//   for TClonesArray is called.
//   The collection itself cannot be a TClonesArray.
//         
//   The function returns the total number of branches created.

   if (list == 0) return 0;
   TObject *obj;
   Int_t nbranches = GetListOfBranches()->GetEntries();
   if (list->InheritsFrom(TClonesArray::Class())) {
         Error("Branch", "Cannot call this constructor for a TClonesArray");
         return 0;
   }
   TIter next(list);

   while ((obj = next())) {
      if (obj->InheritsFrom(TClonesArray::Class())) {
         TClonesArray *clones = (TClonesArray*)obj;
         Branch(clones->GetName(),clones->ClassName(),
                clones->GetObjectRef(obj),bufsize,splitlevel-1);
      } else if (obj->InheritsFrom(TCollection::Class())) {
         TCollection *col = (TCollection*)obj;
         Branch(col,bufsize,splitlevel-1);
      } else {
         Branch(obj->GetName(),obj->ClassName(),
                list->GetObjectRef(obj),bufsize,splitlevel-1);
      }
   }
   return GetListOfBranches()->GetEntries() - nbranches;
}


Magestro Daniel wrote:
> 
> > Hi Thomas,
> >
> > Except for a TClonesArray, it does not make sense to use the
> > split mode for a TCollection. Replace the line :
> >
> >     t->Branch("TList", "TList", &list, 32000, 9);
> > by
> >     t->Branch("TList", "TList", &list, 32000, 0);
> >
> > Rene Brun
> 
> Hi - a long response to Rene's answer... there may be situations where a
> non-zero split mode would be useful for TList's or TObjArray's.
> 
> In our analysis (HADES) the event structure is defined by a single TObjArray
> which contains several other TObjArray's (one for each detector component),
> each of which contains objects with TClonesArray's of the different data
> containers for that detector.  Schematically:
> 
>    TObjArray             TObjArray's
>    ---------             -----------
>    event data  CONTAINS  objArray1 (drift chambers (DC))
>    structure             objArray2 (ring cherenkov)
>                          objArray3 (tof wall)
>                          .....
> 
>                                    Categories
>                                    --------------
>    Then, e.g., objArray1 CONTAINS  category1 (TClonesArray of DC raw data)
>    'category' objects which each   category2 (TClonesArray of DC cal data)
>    contain a TClonesArray          category3 (TClonesArray of DC hit data)
>                                    .....
> 
> The TObjArray's are necessary because there are different types of
> 'category' objects for optimization reasons.
> 
> Currently, we split the two TObjArray 'levels' by hand, creating new
> branches for each object.  We implemented this 3 years ago.  However, our
> manual splitting essentially overrides much of Root's powerful automatic
> splitting implemented in Root v3.xx.
> 
> If TObjArray's were automatically split by Root, then a new branch
> containing its objects as sub-branches (which may be split further) would be
> created.  This would be very useful for us, and I can think of other uses,
> too.  Is this technically difficult?
> 
> Otherwise, to fully implement Root's automatic splitting scheme in an
> elegant way, we may have to eliminate our different types of data
> 'categories' and replace all of our TObjArray's by TClonesArrays.  This is
> not a big deal, but it involves changing our concept somewhat...
> 
> Best regards,
> Dan
>                  --------------------------------------------
>                 | Dr. Daniel Magestro   +49-6159-71-2147     |
>                 | magestro@gsi.de       GSI/Kernphysik I     |
>                 | www.gsi.de/~magestro  Planckstr. 1         |
>                 | (last updated Nov.8)  64291 Darmstadt (DE) |
>                  --------------------------------------------



This archive was generated by hypermail 2b29 : Tue Jan 01 2002 - 17:51:08 MET