1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
19 // marian.ivanov@cern.ch
21 // ------------------------------------------------------------------------------------------------
23 // Standard stream (cout) like input for the tree
24 // Run and see TTreeStreamer::Test() - to see TTreeStreamer functionality
25 // ------------------------------------------------------------------------------------------------
27 // -------------------------------------------------------------------------------------------------
29 // Redirect file to different TTreeStreams
30 // Run and see TTreeSRedirector::Test() as an example of TTreeSRedirectorer functionality
35 #include <TDirectory.h>
36 #include <TObjArray.h>
38 #include "TTreeStream.h"
39 // includes for test procedures
43 ClassImp(TTreeDataElement)
45 ClassImp(TTreeSRedirector)
49 void TTreeStream::Test()
53 TFile *ftest = new TFile("teststreamer.root","recreate");
54 if (!ftest) ftest = new TFile("teststreamer.root","new");
56 //create to streems Tree1 and Tree2
57 TTreeStream stream1("Tree1");
58 TTreeStream stream2("Tree2");
63 TObject *po = new TObject;
64 TObject *po2 = new TObject;
65 for (Int_t i=0;i<100000;i++) {
68 po2->SetUniqueID(i*100);
72 // The data layout of stream is defined during first invocation of streamer.
73 // Endl is the trigger which define the end of structure.
75 // The name of branch can be specified using strings with = at the the end
76 // if string is not specified automatic convention is u (sed B0, B1, ...Bn)
77 stream1<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
80 //3.) just another example - we can fill the same tree with different objects
83 stream2<<f2<<po2<<"\n";
86 //4.) Close the streeamers (Write the streamed tree's to the file) and close the corresponding file.
93 //5.) and now see results in file tteststreamer.root
96 void TTreeSRedirector::Test2()
99 //Example test function to show functionality of TTreeSRedirector
102 //1.)create the redirector associated with file (testredirector.root)
105 TFile* file = new TFile("test.root","recreate");
106 TTreeSRedirector *pmistream= new TTreeSRedirector();
107 TTreeSRedirector &mistream = *pmistream;
111 TObject *po = new TObject;
112 TObject *po2 = new TObject;
113 for (Int_t i=0;i<100000;i++) {
116 po2->SetUniqueID(i*100);
119 //2.) create the tree with identifier specified by first argument
120 // layout specified by sequence of arguments
121 // Tree identifier has to be specified as first argument !!!
122 // if the tree and layout was already defined the consistency if layout is checked
123 // if the data are consisten fill given tree
124 // the name of branch can be specified using strings with = at the the end
125 // if string is not specified use automatic convention B0, B1, ...Bn
126 mistream<<"TreeIdentifier"<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
130 //3.) just another example - we can fill the same tree with different objects
132 mistream<<"TreeK"<<f<<po<<"\n";
133 mistream<<"TreeK"<<f2<<po2<<"\n";
136 //4.) write the streamed tree's to the file and close the corresponding file in destructor
141 //5.) and now see results in file testredirector.root
144 void TTreeSRedirector::Test()
147 //Example test function to show functionality of TTreeSRedirector
150 //1.)create the redirector associated with file (testredirector.root)
153 TTreeSRedirector *pmistream= new TTreeSRedirector("testredirector.root");
154 TTreeSRedirector &mistream = *pmistream;
158 TObject *po = new TObject;
159 TObject *po2 = new TObject;
160 for (Int_t i=0;i<100000;i++) {
163 po2->SetUniqueID(i*100);
166 //2.) create the tree with identifier specified by first argument
167 // layout specified by sequence of arguments
168 // Tree identifier has to be specified as first argument !!!
169 // if the tree and layout was already defined the consistency if layout is checked
170 // if the data are consisten fill given tree
171 // the name of branch can be specified using strings with = at the the end
172 // if string is not specified use automatic convention B0, B1, ...Bn
173 mistream<<"TreeIdentifier"<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
177 //3.) just another example - we can fill the same tree with different objects
179 mistream<<"TreeK"<<f<<po<<"\n";
180 mistream<<"TreeK"<<f2<<po2<<"\n";
183 //4.) write the streamed tree's to the file and close the corresponding file in destructor
187 //5.) and now see results in file testredirector.root
190 void TTreeSRedirector::UnitTest(Int_t testEntries){
194 UnitTestSparse(0.5,testEntries);
195 UnitTestSparse(0.1,testEntries);
196 UnitTestSparse(0.01,testEntries);
199 void TTreeSRedirector::UnitTestSparse(Double_t scale, Int_t testEntries){
201 // Unit test for the TTreeSRedirector
202 // 1.) Test TTreeRedirector
203 // a.) Fill tree with random vectors
204 // b.) Fill downscaled version of vectors
205 // c.) The same skipping first entry
206 // 2.) Check results wtitten to terminale
207 // a.) Disk consumption
208 // skip data should be scale time smaller than full
209 // zerro replaced ata should be compresed time smaller than full
210 // b.) Test invariants
211 // Input parameter scale => downscaling of sprse element
213 if (scale<=0) scale=1;
214 if (scale>1) scale=1;
215 TTreeSRedirector *pcstream = new TTreeSRedirector("testpcstreamSparse.root","recreate");
216 for (Int_t ientry=0; ientry<testEntries; ientry++){
217 TVectorD vecRandom(200);
218 TVectorD vecZerro(200); // zerro vector
219 for (Int_t j=0; j<200; j++) vecRandom[j]=j+ientry+0.1*gRandom->Rndm();
220 Bool_t isSelected= (gRandom->Rndm()<scale);
221 TVectorD *pvecFull = &vecRandom;
222 TVectorD *pvecSparse = isSelected ? &vecRandom:0;
223 TVectorD *pvecSparse0 = isSelected ? &vecRandom:0;
224 TVectorD *pvecSparse1 = isSelected ? &vecRandom:&vecZerro;
228 pvecSparse=&vecRandom;
230 (*pcstream)<<"Full"<< // stored all vectors
234 (*pcstream)<<"SparseSkip"<< // fraction of vectors stored
236 "vec.="<<pvecSparse<<
238 (*pcstream)<<"SparseSkip0"<< // fraction with -pointer
240 "vec.="<<pvecSparse0<<
242 (*pcstream)<<"SparseZerro"<< // all vectors filled, franction filled with 0
244 "vec.="<<pvecSparse1<<
251 TFile* f = TFile::Open("testpcstreamSparse.root");
252 TTree * treeFull = (TTree*)f->Get("Full");
253 TTree * treeSparseSkip = (TTree*)f->Get("SparseSkip");
254 TTree * treeSparseSkip0 = (TTree*)f->Get("SparseSkip0");
255 TTree * treeSparseZerro = (TTree*)f->Get("SparseZerro");
258 Double_t ratio=(1./scale)*treeSparseSkip->GetZipBytes()/Double_t(treeFull->GetZipBytes());
259 Double_t ratio0=(1./scale)*treeSparseSkip0->GetZipBytes()/Double_t(treeFull->GetZipBytes());
260 Double_t ratio1=(1./scale)*treeSparseZerro->GetZipBytes()/Double_t(treeFull->GetZipBytes());
261 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tRatioSkip\t%f\n",scale,ratio);
262 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tRatioSkip0\t%f\n",scale,ratio0);
263 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tRatioZerro\t%f\n",scale,ratio1);
265 Int_t outlyersSparseSkip=treeSparseSkip->Draw("1","(vec.fElements-ientry-Iteration$-0.5)>0.5","goff");
266 Int_t outlyersSparseSkip0=treeSparseSkip0->Draw("1","(vec.fElements-ientry-Iteration$-0.5)>0.5","goff");
267 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tOutlyersSkip\t%d\n",scale,outlyersSparseSkip!=0);
268 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tOutlyersSkip0\t%d\n",scale,outlyersSparseSkip0!=0);
269 // c.) Number of entries
271 Int_t entries=treeFull->GetEntries();
272 Int_t entries0=treeSparseSkip0->GetEntries();
273 Bool_t isOKStat =(entries==entries0);
274 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tEntries\t%d\n",scale,isOKStat);
277 TVectorD *pvecRead = 0;
278 treeSparseSkip0->SetBranchAddress("vec.",&pvecRead);
280 for (Int_t ientry=0; ientry<testEntries; ientry++){
281 if (!pvecRead) continue;
282 if (pvecRead->GetNrows()==0) continue;
283 if (TMath::Abs((*pvecRead)[0]-ientry)>0.5) readOK=kFALSE;
285 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tReadOK\t%d\n",scale,readOK);
288 Bool_t isOK=(outlyersSparseSkip0==0)&&isOKStat&&readOK;
289 printf("#UnitTest:\tTTreeSRedirector::TestSparse(%f)\tisOk\t%d\n",scale,isOK);
293 TTreeSRedirector::TTreeSRedirector(const char *fname,const char * option) :
295 fDirectoryOwner(kTRUE),
303 fDirectory = new TFile(fname,option);
307 fDirectory = gDirectory;
308 fDirectoryOwner = kFALSE;
312 TTreeSRedirector::~TTreeSRedirector()
317 Close(); //write the tree to the selected file
324 void TTreeSRedirector::StoreObject(TObject* object){
328 TDirectory * backup = gDirectory;
331 if (backup) backup->cd();
334 void TTreeSRedirector::SetDirectory(TDirectory *sfile){
336 // Set the external file
337 // In case other file already attached old file is closed before
338 // Redirector will be the owner of file ?
339 if (fDirectory && fDirectoryOwner) {
346 TTreeStream & TTreeSRedirector::operator<<(Int_t id)
349 // return reference to the data layout with given identifier
350 // if not existing - creates new
351 if (!fDataLayouts) fDataLayouts = new TObjArray(10000);
352 TTreeStream *clayout=0;
353 Int_t entries = fDataLayouts->GetEntriesFast();
354 for (Int_t i=0;i<entries;i++){
355 TTreeStream * layout = (TTreeStream*)fDataLayouts->At(i);
356 if (!layout) continue;
357 if (layout->fId==id) {
363 TDirectory * backup = gDirectory;
366 snprintf(chname,100,"Tree%d",id);
367 clayout = new TTreeStream(chname);
369 fDataLayouts->AddAt(clayout,entries);
370 if (backup) backup->cd();
375 void TTreeSRedirector::SetExternalTree(const char* name, TTree* externalTree)
377 TTreeStream *clayout=(TTreeStream*)fDataLayouts->FindObject(name);
380 TDirectory * backup = gDirectory;
382 clayout = new TTreeStream(name,externalTree);
384 clayout->SetName(name);
385 Int_t entries = fDataLayouts->GetEntriesFast();
386 fDataLayouts->AddAt(clayout,entries);
387 if (backup) backup->cd();
390 // AliError(Form("identifier %s already associated",name));
394 TTreeStream & TTreeSRedirector::operator<<(const char* name)
397 // return reference to the data layout with given identifier
398 // if not existing - creates new
399 if (!fDataLayouts) fDataLayouts = new TObjArray(10000);
400 TTreeStream *clayout=(TTreeStream*)fDataLayouts->FindObject(name);
401 Int_t entries = fDataLayouts->GetEntriesFast();
404 TDirectory * backup = gDirectory;
406 clayout = new TTreeStream(name);
408 clayout->SetName(name);
409 fDataLayouts->AddAt(clayout,entries);
410 if (backup) backup->cd();
418 void TTreeSRedirector::Close(){
421 TDirectory * backup = gDirectory;
424 Int_t entries = fDataLayouts->GetEntriesFast();
425 for (Int_t i=0;i<entries;i++){
426 TTreeStream * layout = (TTreeStream*)fDataLayouts->At(i);
428 if (layout->fTree) layout->fTree->Write(layout->GetName());
434 if (backup) backup->cd();
437 //-------------------------------------------------------------
438 TTreeDataElement:: TTreeDataElement(Char_t type) :
450 TTreeDataElement:: TTreeDataElement(TDataType* type) :
462 TTreeDataElement:: TTreeDataElement(TClass* cl) :
474 //-------------------------------------------------------------------
475 TTreeStream::TTreeStream(const char *treename, TTree* externalTree):
476 TNamed(treename,treename),
489 if (!fTree) fTree = new TTree(treename, treename);
492 TTreeStream::~TTreeStream()
503 void TTreeStream::Close()
506 // Flush data to disk and close
511 Int_t TTreeStream::CheckIn(Char_t type, void *pointer)
514 // Insert object of given type
516 if (!fElements) fElements = new TObjArray(10000);
517 if (fElements->GetSize()<=fCurrentIndex) fElements->Expand(fCurrentIndex*2);
518 TTreeDataElement* element = (TTreeDataElement*)fElements->At(fCurrentIndex);
520 element = new TTreeDataElement(type);
523 if (fNextName.Length()>0){
524 if (fNextNameCounter==0){
525 snprintf(name,1000,"%s",(const char*)fNextName);
527 if (fNextNameCounter>0){
528 snprintf(name,1000,"%s%d",(const char*)fNextName,fNextNameCounter);
532 snprintf(name,1000,"B%d.",fCurrentIndex);
534 element->SetName(name);
536 element->SetPointer(pointer);
537 fElements->AddAt(element,fCurrentIndex);
539 return 0; //new element added
541 if (element->GetType()!=type){
543 return 1; //mismatched data element
545 element->SetPointer(pointer);
550 Int_t TTreeStream::CheckIn(TObject *pObject){
555 if (pObject) pClass=pObject->IsA();
556 if (!fElements) fElements = new TObjArray(1000);
557 TTreeDataElement* element = (TTreeDataElement*)fElements->At(fCurrentIndex);
559 element = new TTreeDataElement(pClass);
562 if (fNextName.Length()>0){
563 if (fNextNameCounter==0){
564 snprintf(name,1000,"%s",(const char*)fNextName);
566 if (fNextNameCounter>0){
567 snprintf(name,1000,"%s%d",(const char*)fNextName,fNextNameCounter);
571 snprintf(name,1000,"B%d",fCurrentIndex);
573 element->SetName(name);
575 element->SetPointer(pObject);
576 fElements->AddAt(element,fCurrentIndex);
578 return 0; //new element added
580 if (element->fClass==0) {
581 element->fClass=pClass;
583 if (element->fClass!=pClass && pClass!=0){
585 return 1; //mismatched data element
588 element->SetPointer(pObject);
593 void TTreeStream::BuildTree(){
597 //if (fTree && fTree->GetEntries()>0) return;
598 Int_t entriesFilled=0;
600 fTree = new TTree(GetName(),GetName());
602 entriesFilled=fTree->GetEntries();
604 Int_t entries = fElements->GetEntriesFast();
605 if (!fBranches) fBranches = new TObjArray(entries);
607 for (Int_t i=0;i<entries;i++){
609 TTreeDataElement* element = (TTreeDataElement*)fElements->At(i);
610 if (fBranches->At(i)) continue;
612 if (element->GetName()[0]==0){
613 snprintf(bname1,1000,"B%d",i);
616 snprintf(bname1,1000,"%s",element->GetName());
618 if (element->fClass){
619 if (element->fClass->GetBaseClass("TClonesArray")){
620 TBranch * br = fTree->Branch(bname1,element->fClass->GetName(),&(element->fPointer));
621 if (entriesFilled!=0) {
623 for (Int_t ientry=0; ientry<entriesFilled;ientry++) br->Fill();
624 br->SetAddress(&(element->fPointer));
626 fBranches->AddAt(br,i);
629 TBranch * br = fTree->Branch(bname1,element->fClass->GetName(),&(element->fPointer));
630 if (entriesFilled!=0) {
632 for (Int_t ientry=0; ientry<entriesFilled;ientry++) br->Fill();
633 br->SetAddress(&(element->fPointer));
635 fBranches->AddAt(br,i);
638 if (element->GetType()>0){
640 snprintf(bname2,1000,"B%d/%c",i,element->GetType());
641 TBranch * br = fTree->Branch(bname1,element->fPointer,bname2);
642 if (entriesFilled!=0) {
644 for (Int_t ientry=0; ientry<entriesFilled;ientry++) br->Fill();
645 br->SetAddress(element->fPointer);
648 fBranches->AddAt(br,i);
653 void TTreeStream::Fill(){
658 Int_t entries=fElements->GetEntriesFast();
659 if (entries>fTree->GetNbranches()) BuildTree();
660 for (Int_t i=0;i<entries;i++){
661 TTreeDataElement* el = (TTreeDataElement*)fElements->At(i);
663 if (!el->GetType()) continue;
664 TBranch * br = (TBranch*)fBranches->At(i);
666 if (el->GetType()) br->SetAddress(el->fPointer);
669 if (fStatus==0) fTree->Fill(); //fill only in case of non conflicts
674 TTreeStream & TTreeStream::Endl()
677 // Perform pseudo endl operation
679 if (fTree->GetNbranches()==0) BuildTree();
687 TTreeStream &TTreeStream::operator<<(const Char_t *name)
696 //if tree was already defined ignore
697 if (fTree->GetEntries()>0) return *this;
698 //check branch name if tree was not
701 for (last=0;;last++){
702 if (name[last]==0) break;
705 if (last>0&&name[last-1]=='='){