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"
40 ClassImp(TTreeDataElement)
42 ClassImp(TTreeSRedirector)
46 void TTreeStream::Test()
50 TFile *ftest = new TFile("teststreamer.root","recreate");
51 if (!ftest) ftest = new TFile("teststreamer.root","new");
53 //create to streems Tree1 and Tree2
54 TTreeStream stream1("Tree1");
55 TTreeStream stream2("Tree2");
60 TObject *po = new TObject;
61 TObject *po2 = new TObject;
62 for (Int_t i=0;i<100000;i++) {
65 po2->SetUniqueID(i*100);
69 // The data layout of stream is defined during first invocation of streamer.
70 // Endl is the trigger which define the end of structure.
72 // The name of branch can be specified using strings with = at the the end
73 // if string is not specified automatic convention is u (sed B0, B1, ...Bn)
74 stream1<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
77 //3.) just another example - we can fill the same tree with different objects
80 stream2<<f2<<po2<<"\n";
83 //4.) Close the streeamers (Write the streamed tree's to the file) and close the corresponding file.
90 //5.) and now see results in file tteststreamer.root
93 void TTreeSRedirector::Test2()
96 //Example test function to show functionality of TTreeSRedirector
99 //1.)create the redirector associated with file (testredirector.root)
102 TFile* file = new TFile("test.root","recreate");
103 TTreeSRedirector *pmistream= new TTreeSRedirector();
104 TTreeSRedirector &mistream = *pmistream;
108 TObject *po = new TObject;
109 TObject *po2 = new TObject;
110 for (Int_t i=0;i<100000;i++) {
113 po2->SetUniqueID(i*100);
116 //2.) create the tree with identifier specified by first argument
117 // layout specified by sequence of arguments
118 // Tree identifier has to be specified as first argument !!!
119 // if the tree and layout was already defined the consistency if layout is checked
120 // if the data are consisten fill given tree
121 // the name of branch can be specified using strings with = at the the end
122 // if string is not specified use automatic convention B0, B1, ...Bn
123 mistream<<"TreeIdentifier"<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
127 //3.) just another example - we can fill the same tree with different objects
129 mistream<<"TreeK"<<f<<po<<"\n";
130 mistream<<"TreeK"<<f2<<po2<<"\n";
133 //4.) write the streamed tree's to the file and close the corresponding file in destructor
138 //5.) and now see results in file testredirector.root
141 void TTreeSRedirector::Test()
144 //Example test function to show functionality of TTreeSRedirector
147 //1.)create the redirector associated with file (testredirector.root)
150 TTreeSRedirector *pmistream= new TTreeSRedirector("testredirector.root");
151 TTreeSRedirector &mistream = *pmistream;
155 TObject *po = new TObject;
156 TObject *po2 = new TObject;
157 for (Int_t i=0;i<100000;i++) {
160 po2->SetUniqueID(i*100);
163 //2.) create the tree with identifier specified by first argument
164 // layout specified by sequence of arguments
165 // Tree identifier has to be specified as first argument !!!
166 // if the tree and layout was already defined the consistency if layout is checked
167 // if the data are consisten fill given tree
168 // the name of branch can be specified using strings with = at the the end
169 // if string is not specified use automatic convention B0, B1, ...Bn
170 mistream<<"TreeIdentifier"<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
174 //3.) just another example - we can fill the same tree with different objects
176 mistream<<"TreeK"<<f<<po<<"\n";
177 mistream<<"TreeK"<<f2<<po2<<"\n";
180 //4.) write the streamed tree's to the file and close the corresponding file in destructor
184 //5.) and now see results in file testredirector.root
188 TTreeSRedirector::TTreeSRedirector(const char *fname,const char * option) :
190 fDirectoryOwner(kTRUE),
198 fDirectory = new TFile(fname,option);
202 fDirectory = gDirectory;
203 fDirectoryOwner = kFALSE;
207 TTreeSRedirector::~TTreeSRedirector()
212 Close(); //write the tree to the selected file
219 void TTreeSRedirector::StoreObject(TObject* object){
223 TDirectory * backup = gDirectory;
226 if (backup) backup->cd();
229 void TTreeSRedirector::SetDirectory(TDirectory *sfile){
231 // Set the external file
232 // In case other file already attached old file is closed before
233 // Redirector will be the owner of file ?
234 if (fDirectory && fDirectoryOwner) {
241 TTreeStream & TTreeSRedirector::operator<<(Int_t id)
244 // return reference to the data layout with given identifier
245 // if not existing - creates new
246 if (!fDataLayouts) fDataLayouts = new TObjArray(10000);
247 TTreeStream *clayout=0;
248 Int_t entries = fDataLayouts->GetEntriesFast();
249 for (Int_t i=0;i<entries;i++){
250 TTreeStream * layout = (TTreeStream*)fDataLayouts->At(i);
251 if (!layout) continue;
252 if (layout->fId==id) {
258 TDirectory * backup = gDirectory;
261 snprintf(chname,100,"Tree%d",id);
262 clayout = new TTreeStream(chname);
264 fDataLayouts->AddAt(clayout,entries);
265 if (backup) backup->cd();
270 void TTreeSRedirector::SetExternalTree(const char* name, TTree* externalTree)
272 TTreeStream *clayout=(TTreeStream*)fDataLayouts->FindObject(name);
275 TDirectory * backup = gDirectory;
277 clayout = new TTreeStream(name,externalTree);
279 clayout->SetName(name);
280 Int_t entries = fDataLayouts->GetEntriesFast();
281 fDataLayouts->AddAt(clayout,entries);
282 if (backup) backup->cd();
285 // AliError(Form("identifier %s already associated",name));
289 TTreeStream & TTreeSRedirector::operator<<(const char* name)
292 // return reference to the data layout with given identifier
293 // if not existing - creates new
294 if (!fDataLayouts) fDataLayouts = new TObjArray(10000);
295 TTreeStream *clayout=(TTreeStream*)fDataLayouts->FindObject(name);
296 Int_t entries = fDataLayouts->GetEntriesFast();
299 TDirectory * backup = gDirectory;
301 clayout = new TTreeStream(name);
303 clayout->SetName(name);
304 fDataLayouts->AddAt(clayout,entries);
305 if (backup) backup->cd();
313 void TTreeSRedirector::Close(){
316 TDirectory * backup = gDirectory;
319 Int_t entries = fDataLayouts->GetEntriesFast();
320 for (Int_t i=0;i<entries;i++){
321 TTreeStream * layout = (TTreeStream*)fDataLayouts->At(i);
323 if (layout->fTree) layout->fTree->Write(layout->GetName());
329 if (backup) backup->cd();
332 //-------------------------------------------------------------
333 TTreeDataElement:: TTreeDataElement(Char_t type) :
345 TTreeDataElement:: TTreeDataElement(TDataType* type) :
357 TTreeDataElement:: TTreeDataElement(TClass* cl) :
369 //-------------------------------------------------------------------
370 TTreeStream::TTreeStream(const char *treename, TTree* externalTree):
371 TNamed(treename,treename),
384 if (!fTree) fTree = new TTree(treename, treename);
387 TTreeStream::~TTreeStream()
398 void TTreeStream::Close()
401 // Flush data to disk and close
406 Int_t TTreeStream::CheckIn(Char_t type, void *pointer)
409 // Insert object of given type
411 if (!fElements) fElements = new TObjArray(10000);
412 if (fElements->GetSize()<=fCurrentIndex) fElements->Expand(fCurrentIndex*2);
413 TTreeDataElement* element = (TTreeDataElement*)fElements->At(fCurrentIndex);
415 element = new TTreeDataElement(type);
418 if (fNextName.Length()>0){
419 if (fNextNameCounter==0){
420 snprintf(name,1000,"%s",(const char*)fNextName);
422 if (fNextNameCounter>0){
423 snprintf(name,1000,"%s%d",(const char*)fNextName,fNextNameCounter);
427 snprintf(name,1000,"B%d.",fCurrentIndex);
429 element->SetName(name);
431 element->SetPointer(pointer);
432 fElements->AddAt(element,fCurrentIndex);
434 return 0; //new element added
436 if (element->GetType()!=type){
438 return 1; //mismatched data element
440 element->SetPointer(pointer);
445 Int_t TTreeStream::CheckIn(TObject *o){
450 if (!fElements) fElements = new TObjArray(1000);
451 TTreeDataElement* element = (TTreeDataElement*)fElements->At(fCurrentIndex);
453 element = new TTreeDataElement(o->IsA());
456 if (fNextName.Length()>0){
457 if (fNextNameCounter==0){
458 snprintf(name,1000,"%s",(const char*)fNextName);
460 if (fNextNameCounter>0){
461 snprintf(name,1000,"%s%d",(const char*)fNextName,fNextNameCounter);
465 snprintf(name,1000,"B%d",fCurrentIndex);
467 element->SetName(name);
469 element->SetPointer(o);
470 fElements->AddAt(element,fCurrentIndex);
472 return 0; //new element added
474 if (element->fClass!=o->IsA()){
476 return 1; //mismatched data element
478 element->SetPointer(o);
483 void TTreeStream::BuildTree(){
487 if (fTree && fTree->GetEntries()>0) return;
488 if (!fTree) fTree = new TTree(GetName(),GetName());
489 Int_t entries = fElements->GetEntriesFast();
490 fBranches = new TObjArray(entries);
492 for (Int_t i=0;i<entries;i++){
494 TTreeDataElement* element = (TTreeDataElement*)fElements->At(i);
496 if (element->GetName()[0]==0){
497 snprintf(bname1,1000,"B%d",i);
500 snprintf(bname1,1000,"%s",element->GetName());
502 if (element->fClass){
503 if (element->fClass->GetBaseClass("TClonesArray")){
504 TBranch * br = fTree->Branch(bname1,element->fClass->GetName(),&(element->fPointer));
505 fBranches->AddAt(br,i);
508 TBranch * br = fTree->Branch(bname1,element->fClass->GetName(),&(element->fPointer));
509 fBranches->AddAt(br,i);
512 if (element->GetType()>0){
514 snprintf(bname2,1000,"B%d/%c",i,element->GetType());
515 TBranch * br = fTree->Branch(bname1,element->fPointer,bname2);
516 fBranches->AddAt(br,i);
521 void TTreeStream::Fill(){
526 Int_t entries=fElements->GetEntriesFast();
527 if (entries>fTree->GetNbranches()) BuildTree();
528 for (Int_t i=0;i<entries;i++){
529 TTreeDataElement* el = (TTreeDataElement*)fElements->At(i);
531 if (!el->GetType()) continue;
532 TBranch * br = (TBranch*)fBranches->At(i);
534 if (el->GetType()) br->SetAddress(el->fPointer);
537 if (fStatus==0) fTree->Fill(); //fill only in case of non conflicts
542 TTreeStream & TTreeStream::Endl()
545 // Perform pseudo endl operation
547 if (fTree->GetNbranches()==0) BuildTree();
555 TTreeStream &TTreeStream::operator<<(const Char_t *name)
564 //if tree was already defined ignore
565 if (fTree->GetEntries()>0) return *this;
566 //check branch name if tree was not
569 for (last=0;;last++){
570 if (name[last]==0) break;
573 if (last>0&&name[last-1]=='='){