Re: [ROOT] sub-branch SetAddress/GetEvent behavior question

From: Rene Brun (Rene.Brun@cern.ch)
Date: Thu Aug 23 2001 - 11:10:53 MEST


Hi Kate,

Looking at your file, I see that the branch contents are wrong.
I do not understand how your Tree had been created. I do not see any problems
with branches having the same attribute names in different classes.
Could you run the example below on your system and let me know if you get
correct results? if results are correct (this should work with versions
2.25,3.00 or 3.01), we have to understand what is different in your
configuration. Could you reproduce the program writing the file
with a small test program ?

void dots() {
   // test of branches with the same class
   
   TFile f("dots.root","recreate");
   TTree *T = new TTree("T","dots tree");
   TLine *l1 = new TLine();
   TLine *l2 = new TLine();
   TBox *l3  = new TBox();
   T->Branch("line1.","TLine",&l1,32000,99);
   T->Branch("line2.","TLine",&l2,32000,99);
   T->Branch("box3.","TBox", &l3,32000,99);
   for (Int_t i=0;i<1000;i++) {
      l1->SetX1(100000+i);
      l1->SetY1(110000+i);
      l1->SetX2(120000+i);
      l1->SetY2(130000+i);
      l2->SetX1(200000+i);
      l2->SetY1(210000+i);
      l2->SetX2(220000+i);
      l2->SetY2(230000+i);
      l3->SetX1(300000+i);
      l3->SetY1(310000+i);
      l3->SetX2(320000+i);
      l3->SetY2(330000+i);
      l1->SetLineColor(1);
      l2->SetLineColor(2);
      l3->SetLineColor(3);
      T->Fill();
   }
   T->Print();
   T->Write();
}
void dotstest() {
   TFile f("dots.root");
   TTree *T = (TTree*)f.Get("T");
   T->Show(123);
}   

Rene Brun

Kate Scholberg wrote:
> 
> ear rooters,
> 
> I'm still having problems reading Root tree contents by directly setting
> sub-branch addresses.  It behaves inconsistently-- sometimes seems to
> work, but doesn't work under other conditions.
> 
> The Event class example I sent on July 30 (macros testread2 and
> testread3) did NOT work.  Philippe Canal suggested a workaround, which
> works for Root 3.01 with Event classes.  But, it doesn't work for Root 3.00
> (because GetValue is apparently not yet implemented for
> TBranchElement).  It also segfaults when I try to implement the
> workaround using my own classes (see macro and link below).  The
> workaround is also rather awkward to implement for my application.
> 
> However, the more obvious approach, e.g. my original code and Rene's
> macro of Jul 6, actually *does* work for Root 3.00 (and 3.01), but
> only under some conditions.  One necessary (but not sufficient)
> condition for working is to *not* load my shared lib.
> It also seems to work fine if the sub-branch names are "plain",
> i.e.  created with a split branch function like
> 
> TBranch *tw2=_tree->Branch("tweedle", "Tweedle",  &twd, 64000,1);
> 
> However I'm dealing with several split branches
> made from classes which have degenerate data member names,
> 
> e.g.
> 
> class Tweedle
> {
>   int dee;
>   int dum;
> }
> 
> class Slithery
> {
>   int dee;
>   int urgghl;
> }
> 
> leading to ambiguities in finding the branch "dee" from its name...
> The solution to this (I think), which I found on Thu Apr 22 1999
> roottalk, is to add "." to the branch name, i.e. create with TBranch
> *tw2=_tree->Branch("tweedle.", "Tweedle", &twd, 64000,1);
> 
> Then I try to access the sub-branch directly using e.g.
> 
> int blah=0;
> amstree->SetBranchAddress("tweedle.dee", &blah);
> 
> or
> 
> amstree->SetBranchAddress("slithery.dee", &blah);
> 
> Well, this works, but only for the *first* super-branch... GetEvent seems
> to read the contents correctly for the first super-branch, but nothing
> is read for subsequent branches (although GetBranch seems to
> find the branches successfully).
> 
> I guess I could rewrite the classes to have non-degenerate member names,
> but I'm not the only user of the classes and I'd rather not.
> 
> Any ideas?  Is such behavior a bug or a feature?  If it's a bug, will
> it be fixed?  Anything I should check about my classes? (but I
> strongly suspect that it's not only a problem with my classes, since
> the "obvious approach" doesn't work at all with the Event test class--
> see earlier in this thread).  Any other workaround ideas?  (It would
> be nice if I could get this to work with Root 3.00).
> 
> Thanks,
> 
> Kate Scholberg
> schol@mit.edu
> 
> ============================================================
> 
> Following are some macros.   Macros and files can be found at
> http://www-lns.mit.edu/~schol/roottest/
> 
> Same behavior for Root 3.00, 3.01, Linux, OSF1.
> 
> testread7 ("obvious approach") reads the first super-branch (event02)
> correctly, but get nothing for the next ones (e.g. beta02)
> 
> testread8 is Philippe's workaround (for 3.01).  Segfaults.
> 
> Please let me know if more info needed (e.g. class definitions).
> 
> ------------------------------------------------------------
> 
> void testread7()
> {
>   TFile file("prv3.root");
>   TTree *amstree = (TTree*)file->Get("AMSRoot");
>   int nevent = amstree->GetEntries();
>   printf("nevent %d\n",nevent);
> 
>   int evtno = 0; int run=0; int time[2];
>   int nbeta = 0;
>   int npart = 0;
>   amstree->SetBranchAddress("event02.Eventno", &evtno);
>   amstree->SetBranchAddress("event02.Run", &run);
>   amstree->SetBranchAddress("event02.Time[2]", time);
> 
>   amstree->SetBranchAddress("beta02.Nbeta", &nbeta);
> 
>   amstree->SetBranchAddress("part02.Npart", &npart);
> 
>   nevent = 10;
>   for(int i=0;i<nevent;i++) {
>     amstree->GetEvent(i);
>     printf("i %d, evtno=%d, run=%d,time %d %d, nbeta=%d npart=%d \n",i,evtno,run
> ,time[0],time[1],nbeta,npart);
> 
>  }
> }
> 
> ------------------------------------------------------------
> 
> void testread8()
> {
>   TFile file("prv3.root");
>   TTree *amstree = (TTree*)file->Get("AMSRoot");
>   int nevent = amstree->GetEntries();
>   printf("nevent %d\n",nevent);
> 
>   Int_t evtno = 0; Int_t run=0; Int_t time[2];
>   Int_t nbeta = 0;
>   Int_t npart = 0;
> 
>   TBranchElement *evbranch = (TBranchElement*)amstree->GetBranch("event02.Event\
> no");
>   TBranchElement *betabranch = (TBranchElement*)amstree->GetBranch("beta02.Nbet\
> a");
> 
>   nevent = 10;
>   for(int i=0;i<nevent;i++) {
>     evbranch->GetEntry(i);
>     betabranch->GetEntry(i);
>     evtno = evbranch->GetValue(0,0);
> //    evtno = Int_t(evbranch->GetValue(0,0,0));
>     nbeta = betabranch->GetValue(0,0);
>     printf("i %d, evtno=%d, run=%d,time %d %d, nbeta=%d npart=%d \n",i,evtno,ru\
> n,time[0],time[1],nbeta,npart);
>   }
> }



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