Re: [ROOT] setting the split mode of a TParticle in a tree

From: Rene Brun (Rene.Brun@cern.ch)
Date: Fri Nov 22 2002 - 18:33:52 MET


Hello Anne-Sylvie,

Yes, you will get a warning with TParticle because this class has a custom
Streamer, but it should work. For example, the small root session:
root > gSystem->Load("libEG")
root > TClonesArray* genphotons = new TClonesArray("TParticle");
root > TTree *tree=new TTree("tree","test");
root > tree->Branch("genphotons",&genphotons,16000,10);
root > tree->Print();

will produce the output below. As you can see the tree is split.

Rene brun

******************************************************************************
*Tree    :tree      : test
*
*Entries :        0 : Total =          275907 bytes  File  Size =
0 *
*        :          : Tree compression factor =   1.00
*
******************************************************************************
*Br    0 :genphotons : genphotons_
*
*Entries :        0 : Total  Size=     287893 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    1 :genphotons.fUniqueID : fUniqueID[genphotons_]
*
*Entries :        0 : Total  Size=      12576 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    2 :genphotons.fBits : fBits[genphotons_]
*
*Entries :        0 : Total  Size=      12552 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    3 :genphotons.fLineColor : fLineColor[genphotons_]
*
*Entries :        0 : Total  Size=      12582 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    4 :genphotons.fLineStyle : fLineStyle[genphotons_]
*
*Entries :        0 : Total  Size=      12582 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    5 :genphotons.fLineWidth : fLineWidth[genphotons_]
*
*Entries :        0 : Total  Size=      12582 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    6 :genphotons.fPdgCode : fPdgCode[genphotons_]
*
*Entries :        0 : Total  Size=      12570 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    7 :genphotons.fStatusCode : fStatusCode[genphotons_]
*
*Entries :        0 : Total  Size=      12588 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    8 :genphotons.fMother[2] : fMother[genphotons_]
*
*Entries :        0 : Total  Size=      12573 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br    9 :genphotons.fDaughter[2] : fDaughter[genphotons_]
*
*Entries :        0 : Total  Size=      12585 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   10 :genphotons.fWeight : fWeight[genphotons_]
*
*Entries :        0 : Total  Size=      12564 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   11 :genphotons.fCalcMass : fCalcMass[genphotons_]
*
*Entries :        0 : Total  Size=      12576 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   12 :genphotons.fPx : fPx[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   13 :genphotons.fPy : fPy[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   14 :genphotons.fPz : fPz[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   15 :genphotons.fE : fE[genphotons_]
*
*Entries :        0 : Total  Size=      12534 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   16 :genphotons.fVx : fVx[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   17 :genphotons.fVy : fVy[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   18 :genphotons.fVz : fVz[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   19 :genphotons.fVt : fVt[genphotons_]
*
*Entries :        0 : Total  Size=      12540 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   20 :genphotons.fPolarTheta : fPolarTheta[genphotons_]
*
*Entries :        0 : Total  Size=      12588 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*
*Br   21 :genphotons.fPolarPhi : fPolarPhi[genphotons_]
*
*Entries :        0 : Total  Size=      12576 bytes  File Size  =
0 *
*Baskets :        0 : Basket Size=      16000 bytes  Compression=   1.00
*
*............................................................................*

On Fri, 22 Nov 2002, Anne-Sylvie Nicollerat wrote:

> Hello,
> 
> I am trying to write a root tree containing TClonesArrays of TParticles.
> 
> It seems like I cannot set the split mode for the branch correctly.
> When I try to plot for instance the energy distribution of the particles
> in the array, I am getting an histogram with always the same value
> (62.3). However when I read the tree with a macro and fill a histogram,
> I get the correct distribution.
> 
> When I run the executable, I get the following error message:
> Warning in <TTree::Bronch>: Using split mode on a class: TParticle with
> a custom Streamer
> 
> Any idea ?
> Thanks a lot
> 
> Anne-Sylvie
> 
> 
> Here is the way I am defining the tree:
> 
>  TClonesArray* genphotons = new TClonesArray("TParticle");
> 
>  tree->Branch("genphotons",&genphotons,16000,10);
> 
> ...................
> 
>   for(ParticleTree::const_iterator iloop=Genphotons.begin();
>       iloop!=Genphotons.end();iloop++)
>     { (....)
>           (*genphotons)[genphotons->GetEntries()] = new
>           TParticle((*iloop)->HepId(),(*iloop)->Status(),
> 
> first_mother,last_mother,first_daughter,last_daughter,
> 
> (*iloop)->px(),(*iloop)->py(),(*iloop)->pz(),(*iloop)->Energy(),
> 
> (*iloop)->x(),(*iloop)->y(),(*iloop)->z(),(*iloop)->t()
>                    );
>     }
> 
>     tree->Fill();
> 
> 
> --
> ---------------------------------------------------------
> Anne-Sylvie Nicollerat
> Swiss Federal Insitute of technology (ETHZ)
> Building 32, Office 3-C21, Phone +41 22 767 3585
> CERN, CH-1211 Geneva 23
> Switzerland
> 
> "Pour le croyant Dieu se trouve au début, pour le physicien au terme de
> toute pensée." (Planck)
> 
> 



This archive was generated by hypermail 2b29 : Sat Jan 04 2003 - 23:51:20 MET