2 // marian.ivanov@cern.ch
4 // ------------------------------------------------------------------------------------------------
6 // Standard stream (cout) like input for the tree
7 // Run and see TTreeStreamer::Test() - to see TTreeStreamer functionality
8 // ------------------------------------------------------------------------------------------------
10 // -------------------------------------------------------------------------------------------------
12 // Redirect file to different TTreeStreams
13 // Run and see TTreeSRedirector::Test() as an example of TTreeSRedirectorer functionality
19 #include "TObjArray.h"
21 #include "TTreeStream.h"
24 ClassImp(TTreeDataElement)
26 ClassImp(TTreeSRedirector)
30 void TTreeStream::Test()
34 TFile *ftest = new TFile("teststreamer.root","recreate");
35 if (!ftest) ftest = new TFile("teststreamer.root","new");
37 //create to streems Tree1 and Tree2
38 TTreeStream stream1("Tree1");
39 TTreeStream stream2("Tree2");
44 TObject *po = new TObject;
45 TObject *po2 = new TObject;
46 for (Int_t i=0;i<100000;i++) {
49 po2->SetUniqueID(i*100);
53 // The data layout of stream is defined during first invocation of streamer.
54 // Endl is the trigger which define the end of structure.
56 // The name of branch can be specified using strings with = at the the end
57 // if string is not specified automatic convention is u (sed B0, B1, ...Bn)
58 stream1<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
61 //3.) just another example - we can fill the same tree with different objects
64 stream2<<f2<<po2<<"\n";
67 //4.) Close the streeamers (Write the streamed tree's to the file) and close the corresponding file.
74 //5.) and now see results in file tteststreamer.root
79 void TTreeSRedirector::Test()
82 //Example test function to show functionality of TTreeSRedirector
85 //1.)create the redirector associated with file (testredirector.root)
88 TTreeSRedirector *pmistream= new TTreeSRedirector("testredirector.root");
89 TTreeSRedirector &mistream = *pmistream;
93 TObject *po = new TObject;
94 TObject *po2 = new TObject;
95 for (Int_t i=0;i<100000;i++) {
98 po2->SetUniqueID(i*100);
101 //2.) create the tree with identifier specified by first argument
102 // layout specified by sequence of arguments
103 // Tree identifier has to be specified as first argument !!!
104 // if the tree and layout was already defined the consistency if layout is checked
105 // if the data are consisten fill given tree
106 // the name of branch can be specified using strings with = at the the end
107 // if string is not specified use automatic convention B0, B1, ...Bn
108 mistream<<"TreeIdentifier"<<"i="<<i<<"ch="<<ch<<"f="<<f<<"po="<<po<<"\n";
112 //3.) just another example - we can fill the same tree with different objects
114 mistream<<"TreeK"<<f<<po<<"\n";
115 mistream<<"TreeK"<<f2<<po2<<"\n";
118 //4.) write the streamed tree's to the file and close the corresponding file in destructor
122 //5.) and now see results in file testredirector.root
126 TTreeSRedirector::TTreeSRedirector(const char *fname) :
127 fFile(new TFile(fname,"recreate")),
134 fFile = new TFile(fname,"new");
138 TTreeSRedirector::~TTreeSRedirector()
143 Close(); //write the tree to the selected file
147 void TTreeSRedirector::StoreObject(TObject* object){
151 TFile * backup = gFile;
154 if (backup) backup->cd();
159 TTreeStream & TTreeSRedirector::operator<<(Int_t id)
162 // return reference to the data layout with given identifier
163 // if not existing - creates new
164 if (!fDataLayouts) fDataLayouts = new TObjArray(10000);
165 TTreeStream *clayout=0;
166 Int_t entries = fDataLayouts->GetEntriesFast();
167 for (Int_t i=0;i<entries;i++){
168 TTreeStream * layout = (TTreeStream*)fDataLayouts->At(i);
169 if (!layout) continue;
170 if (layout->fId==id) {
178 sprintf(chname,"Tree%d",id);
179 clayout = new TTreeStream(chname);
181 fDataLayouts->AddAt(clayout,entries);
187 TTreeStream & TTreeSRedirector::operator<<(const char* name)
190 // return reference to the data layout with given identifier
191 // if not existing - creates new
192 if (!fDataLayouts) fDataLayouts = new TObjArray(10000);
193 TTreeStream *clayout=(TTreeStream*)fDataLayouts->FindObject(name);
194 Int_t entries = fDataLayouts->GetEntriesFast();
198 clayout = new TTreeStream(name);
200 clayout->SetName(name);
201 fDataLayouts->AddAt(clayout,entries);
209 void TTreeSRedirector::Close(){
212 TFile * backup = gFile;
215 Int_t entries = fDataLayouts->GetEntriesFast();
216 for (Int_t i=0;i<entries;i++){
217 TTreeStream * layout = (TTreeStream*)fDataLayouts->At(i);
219 if (layout->fTree) layout->fTree->Write(layout->GetName());
225 if (backup) backup->cd();
230 //-------------------------------------------------------------
231 TTreeDataElement:: TTreeDataElement(Char_t type) :
243 TTreeDataElement:: TTreeDataElement(TDataType* type) :
255 TTreeDataElement:: TTreeDataElement(TClass* cl) :
267 //-------------------------------------------------------------------
268 TTreeStream::TTreeStream(const char *treename):
269 TNamed(treename,treename),
272 fTree(new TTree(treename, treename)),
284 TTreeStream::~TTreeStream()
295 void TTreeStream::Close()
298 // Flush data to disk and close
303 Int_t TTreeStream::CheckIn(Char_t type, void *pointer)
306 // Insert object of given type
308 if (!fElements) fElements = new TObjArray(1000);
309 TTreeDataElement* element = (TTreeDataElement*)fElements->At(fCurrentIndex);
311 element = new TTreeDataElement(type);
314 if (fNextName.Length()>0){
315 if (fNextNameCounter==0){
316 sprintf(name,"%s",(const char*)fNextName);
318 if (fNextNameCounter>0){
319 sprintf(name,"%s%d",(const char*)fNextName,fNextNameCounter);
323 sprintf(name,"B%d.",fCurrentIndex);
325 element->SetName(name);
327 element->SetPointer(pointer);
328 fElements->AddAt(element,fCurrentIndex);
330 return 0; //new element added
332 if (element->GetType()!=type){
334 return 1; //mismatched data element
336 element->SetPointer(pointer);
341 Int_t TTreeStream::CheckIn(TObject *o){
346 if (!fElements) fElements = new TObjArray(1000);
347 TTreeDataElement* element = (TTreeDataElement*)fElements->At(fCurrentIndex);
349 element = new TTreeDataElement(o->IsA());
352 if (fNextName.Length()>0){
353 if (fNextNameCounter==0){
354 sprintf(name,"%s",(const char*)fNextName);
356 if (fNextNameCounter>0){
357 sprintf(name,"%s%d",(const char*)fNextName,fNextNameCounter);
361 sprintf(name,"B%d",fCurrentIndex);
363 element->SetName(name);
365 element->SetPointer(o);
366 fElements->AddAt(element,fCurrentIndex);
368 return 0; //new element added
370 if (element->fClass!=o->IsA()){
372 return 1; //mismatched data element
374 element->SetPointer(o);
379 void TTreeStream::BuildTree(){
383 if (fTree->GetEntries()>0) return;
384 fTree = new TTree(GetName(),GetName());
385 Int_t entries = fElements->GetEntriesFast();
386 fBranches = new TObjArray(entries);
388 for (Int_t i=0;i<entries;i++){
390 TTreeDataElement* element = (TTreeDataElement*)fElements->At(i);
392 if (element->GetName()[0]==0){
393 sprintf(bname1,"B%d",i);
396 sprintf(bname1,element->GetName());
398 if (element->fClass){
399 if (element->fClass->GetBaseClass("TClonesArray")){
400 TBranch * br = fTree->Branch(bname1,element->fClass->GetName(),&(element->fPointer));
401 fBranches->AddAt(br,i);
404 TBranch * br = fTree->Branch(bname1,element->fClass->GetName(),&(element->fPointer));
405 fBranches->AddAt(br,i);
408 if (element->GetType()>0){
410 sprintf(bname2,"B%d/%c",i,element->GetType());
411 TBranch * br = fTree->Branch(bname1,element->fPointer,bname2);
412 fBranches->AddAt(br,i);
417 void TTreeStream::Fill(){
422 Int_t entries=fElements->GetEntriesFast();
423 if (entries>fTree->GetNbranches()) BuildTree();
424 for (Int_t i=0;i<entries;i++){
425 TTreeDataElement* el = (TTreeDataElement*)fElements->At(i);
427 if (!el->GetType()) continue;
428 TBranch * br = (TBranch*)fBranches->At(i);
430 if (el->GetType()) br->SetAddress(el->fPointer);
433 if (fStatus==0) fTree->Fill(); //fill only in case of non conflicts
438 TTreeStream & TTreeStream::Endl()
441 // Perform pseudo endl operation
443 if (fTree->GetNbranches()==0) BuildTree();
451 TTreeStream &TTreeStream::operator<<(Char_t *name)
460 //if tree was already defined ignore
461 if (fTree->GetEntries()>0) return *this;
462 //check branch name if tree was not
465 for (last=0;;last++){
466 if (name[last]==0) break;
469 if (last>0&&name[last-1]=='='){