Re: [ROOT] ROOT: Speed problem

From: Rene Brun (Rene.Brun@cern.ch)
Date: Tue Jul 02 2002 - 14:33:48 MEST


Hi jan,

You should use twice the same file as shown in a modified
version of your example below, otherwise you destroy
completly the Tree cache.

Rene Brun

  UShort_t nchar=0;
  Float_t px[MAX_TRACK];

  TFile *fin = new TFile("oxygen.root","READ");
  if (!(fin.IsOpen())) break;
  TTree *oxy = (TTree*)fin->Get("oxy");
  TFile *fin2 = new TFile("oxygen.root","READ");
  TTree *oxy2 = (TTree*)fin2->Get("oxy");

  oxy->SetBranchAddress("nchar",&nchar);
  oxy->SetBranchAddress("px",px);
  oxy2->SetBranchAddress("nchar",&nchar);
  oxy2->SetBranchAddress("px",px);

  UInt_t nevents=oxy->GetEntries();
  for (UInt_t i=0; i<nevents; i++) {
      oxy->GetEntry(i);
      ....
      using nchar
      using px[...]
      ....
   
      for (UInt_t i2=1; i2<500; i2++) {
           UShort_t event2 = (gRandom->Rndm())*nevents;
           if (event2 == i) event2 = (gRandom->Rndm())*nevents;
           oxy2->GetEntry(event2);
           ...
           using nchar
           using px[...]
           ....
       } // end i2
   } // end i


On Tue, 2 Jul 2002, Jan Musinsky wrote:

>  Hi,
> I have met a stupid (maybe) problem.
> I use tree to analyze ...:
> 
> 
>   UShort_t nchar=0;
>   Float_t px[MAX_TRACK];
> 
>   TFile *fin = new TFile("oxygen.root","READ");
>   if (!(fin.IsOpen())) break;
>   TTree *oxy = (TTree*)fin->Get("oxy");
> 
>   oxy->SetBranchAddress("nchar",&nchar);
>   oxy->SetBranchAddress("px",px);
> 
>   UInt_t nevents=oxy->GetEntries();
>   for (UInt_t i=0; i<nevents; i++) {
>       oxy->GetEntry(i);
>       ....
>       using nchar
>       using px[...]
>       ....
>    
>       for (UInt_t i2=1; i2<500; i2++) {
>            UShort_t event2 = (gRandom->Rndm())*nevents;
>            if (event2 == i) event2 = (gRandom->Rndm())*nevents;
>            oxy->GetEntry(event2);
>            ...
>            using nchar
>            using px[...]
>            ....
>        } // end i2
>    } // end i
> 
> 
> 
> The result is that the analysis of the root takes days, days, ..
> Can I somehow speed it up ?
> Thanks in advance,
> Jan
> 
> 



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