[ROOT] Problem to sum large TTree

From: Andrei Daniel (daniel@suntimpx.jinr.ru)
Date: Tue Jun 13 2000 - 07:41:45 MEST


Hi Rene:

You was right. Compression  is huge. But it was one of
the main reason to use TFile. It allows to work easy
with these type of data. An Original uncompressed file
occupied ~16 Gbytes and we have 3 similar files.
Results of Cube.Print() are here.

Let me continue the discussion. On the local PC
working under W2000 I have only one TFile and
cannot make full check. TFile was copied two times
with different names and code was applied to these
three files. On PC it worked without errors.
The same code applied to the 3 real TFiles give
error on Alpha (Digital UNIX), which was asked in the
previous mail.

       Fatal in <CustomReAlloc2>: storage exhausted
aborting

What may be the difference. I see only that compression
of last TFile is much more higher then compression of two
first TFiles. Can it be a reason of error?

Andrei

FILE1:
****************************************************************************
**
*Tree    :Cube      : Cf252 data of 1995, 1st week
*
*Entries :   993345 : Total  Size = -831707011 bytes  File  Size =
1161634262 *
*        :          : Tree compression factor =  14.15
*
****************************************************************************
**
*Branch  :Spectrum  : sbuf[8192]/s
*
*Entries :   993345 : Total  Size = -838350589 bytes  File Size  =
1154990684 *
*Baskets :   993345 : Basket Size =     16384 bytes  Compression=  14.15
*
*...........................................................................
.*

FILE2:
****************************************************************************
**
*Tree    :Cube      : Cf252 data of 1995
*
*Entries :   993345 : Total  Size = -832558858 bytes  File  Size =
933632868 *
*        :          : Tree compression factor =  17.61
*
****************************************************************************
**
*Branch  :Spectrum  : sbuf[8192]/s
*
*Entries :   993345 : Total  Size = -838350589 bytes  File Size  =
927841137 *
*Baskets :   993345 : Basket Size =     32000 bytes  Compression=  17.61
*
*...........................................................................
.*

FILE3:
****************************************************************************
**
*Tree    :Cube      : Cf252 data of 1995
*
*Entries :   993345 : Total  Size = -833245747 bytes  File  Size =
140854651 *
*        :          : Tree compression factor = 120.38
*
****************************************************************************
**
*Branch  :Spectrum  : sbuf[8192]/s
*
*Entries :   993345 : Total  Size = -838350589 bytes  File Size  =
135749809 *
*Baskets :   993345 : Basket Size =     32000 bytes  Compression= 120.38
*
*...........................................................................
.*

> -----Original Message-----
> From: brun@pcbrun.cern.ch [mailto:brun@pcbrun.cern.ch]On Behalf
> Of Rene Brun
> Sent: Thursday, June 08, 2000 12:00 PM
> To: Andrei Daniel
> Cc: Roottalk
> Subject: Re: [ROOT] Problem to sum 3 large TTree
>
>
> Hi Andrei,
> I do not understand how you have been able to create the original 3 files.
> With 991936 entries, each file should be around 16.2 GBytes
> (uncompressed).
> Probably, you had a huge compression factor such that each file
> was just below
> the 2 GByte size.
> Could you send me the result of FCube.Print() for each original file ?
>
> Rene Brun
>
> Andrei Daniel wrote:
> >
> > Hi Rooters!
> >
> > After 3 ours of work program gave me the next message.
> > I used Alpha working under Digital Unix and
> > ROOT version 2.23/12.
> >
> >       Fatal in <CustomReAlloc2>: storage exhausted
> > aborting
> >
> > I have 3 TFiles. Each file includes the same TTree created
> > from the different data. Now I want to sum these data and
> > write result in the new TFile. What may be wrong or can be
> > changed in the next part of code? Let me note that each TFile
> > can be read separately without errors.
> >
> > Thanks in advance.
> >
> > Andrei
> >
> > const int nx = 1408;
> > const int n  = nx * (nx + 1) / 2;
> > const int nz = 8192;
> >
> >    TROOT sumofcubes ("sumofcubes","Sum Data of 3 Cubes");
> >
> >    struct
> >    {
> >       unsigned short sbuf[nz];
> >    } buf1,buf2,buf3;
> >    struct
> >    {
> >       float fbuf[nz];
> >    } buf4;
> >
> >    int i, n1, n2, ix, iy, iz;
> >
> >    TFile *F1 = new TFile("a1.root", "read");
> >    TTree *Cube1 = (TTree*) F1->Get("Cube");
> >    Cube1->SetBranchAddress("Spectrum", buf1.sbuf);
> >
> >    TFile *F2 = new TFile("a2.root", "read");
> >    TTree *Cube2 = (TTree*) F2->Get("Cube");
> >    Cube2->SetBranchAddress("Spectrum", buf2.sbuf);
> >
> >    TFile *F3 = new TFile("a3.root", "read");
> >    TTree *Cube3 = (TTree*) F3->Get("Cube");
> >    Cube3->SetBranchAddress("Spectrum", buf3.sbuf);
> >
> >
> >    TFile *FOUT = new TFile("as.root","recreate");
> >    TTree *FCube = new TTree("FCube","SUM OF 3 Weeks");
> >    TBranch *Spectrum =
> FCube->Branch("Spectrum",buf4.fbuf,"fbuf[8192]/F");
> >
> >    for(i = 0; i < n; ++i)
> >    {
> >       Cube1->GetEntry(i);
> >       Cube2->GetEntry(i);
> >       Cube3->GetEntry(i);
> >
> >       for(iz = 0; iz < nz; ++iz)
> >       {
> >          buf4.fbuf[iz] = (float) buf1.sbuf[iz]
> >                        + (float) buf2.sbuf[iz]
> >                        + (float) buf3.sbuf[iz];
> >       }
> >
> >       FCube->Fill();
> >    }
> >
> >    FCube->Write();
> >
> >    F1->Close();
> >    F2->Close();
> >    F3->Close();
> >
> >    FOUT->Close();
> > }
> >

----------------------------------------
Andrei Daniel
FLNR, JINR, Dubna 141980, Russia
Tel: (7-09621)64568  Fax: (7-09621)65083
daniel@cv.jinr.ru



This archive was generated by hypermail 2b29 : Tue Jan 02 2001 - 11:50:27 MET