staff.C: to create cernstaff.root, execute tutorial $ROOTSYS/tree/cernbuild.C | Trees I/O, Queries, Graphics | tree.C: This macro displays the Tree data structures |
// Example of macro illustrating how to write a TClonesArray to a TTree // the following tests can be run // Interactive tests // Root > .x tcl.C //no-split interpreted // Root > .x tcl.C(1) //split interpreted // Root > .x tcl.C++ //no-split compiled // Root > .x tcl.C++(1) //split compiled // batch tests: same as above but with no graphics // root -b -q tcl.C // root -b -q tcl.C++ // root -b -q "tcl.C(1)" // root -b -q "tcl.C++(1)" // //Author: Rene Brun #include "TFile.h" #include "TClonesArray.h" #include "TH2.h" #include "TLine.h" #include "TTree.h" #include "TBenchmark.h" #include "TRandom.h" void tclwrite(Int_t split) { // Generate a Tree with a TClonesArray // The array can be split or not TFile f("tcl.root","recreate"); f.SetCompressionLevel(1); //try level 2 also TTree T("T","test tcl"); TClonesArray *arr = new TClonesArray("TLine"); TClonesArray &ar = *arr; T.Branch("tcl",&arr,256000,split); //By default a TClonesArray is created with its BypassStreamer bit set. //However, because TLine has a custom Streamer, this bit was reset //by TTree::Branch above. We set again this bit because the current //version of TLine uses the automatic Streamer. //BypassingStreamer saves space and time. arr->BypassStreamer(); for (Int_t ev=0;ev<10000;ev++) { ar.Clear(); Int_t nlines = Int_t(gRandom->Gaus(50,10)); if(nlines < 0) nlines = 1; for (Int_t i=0;i<nlines;i++) { Float_t x1 = gRandom->Rndm(); Float_t y1 = gRandom->Rndm(); Float_t x2 = gRandom->Rndm(); Float_t y2 = gRandom->Rndm(); new(ar[i]) TLine(x1,y1,x2,y2); } T.Fill(); } T.Print(); T.Write(); } void tclread() { // read file generated by tclwrite // loop on all entries. // histogram center of lines TFile *f = new TFile("tcl.root"); TTree *T = (TTree*)f->Get("T"); TH2F *h2 = new TH2F("h2","center of lines",40,0,1,40,0,1); TClonesArray *arr = new TClonesArray("TLine"); T->GetBranch("tcl")->SetAutoDelete(kFALSE); T->SetBranchAddress("tcl",&arr); Long64_t nentries = T->GetEntries(); for (Long64_t ev=0;ev<nentries;ev++) { arr->Clear(); T->GetEntry(ev); Int_t nlines = arr->GetEntriesFast(); for (Int_t i=0;i<nlines;i++) { TLine *line = (TLine*)arr->At(i); h2->Fill(0.5*(line->GetX1()+line->GetX2()), 0.5*(line->GetY1()+line->GetY2())); } } h2->Draw("lego"); } void tcl(Int_t split=0) { gBenchmark->Start("tcl"); tclwrite(split); tclread(); gBenchmark->Show("tcl"); } tcl.C:1 tcl.C:2 tcl.C:3 tcl.C:4 tcl.C:5 tcl.C:6 tcl.C:7 tcl.C:8 tcl.C:9 tcl.C:10 tcl.C:11 tcl.C:12 tcl.C:13 tcl.C:14 tcl.C:15 tcl.C:16 tcl.C:17 tcl.C:18 tcl.C:19 tcl.C:20 tcl.C:21 tcl.C:22 tcl.C:23 tcl.C:24 tcl.C:25 tcl.C:26 tcl.C:27 tcl.C:28 tcl.C:29 tcl.C:30 tcl.C:31 tcl.C:32 tcl.C:33 tcl.C:34 tcl.C:35 tcl.C:36 tcl.C:37 tcl.C:38 tcl.C:39 tcl.C:40 tcl.C:41 tcl.C:42 tcl.C:43 tcl.C:44 tcl.C:45 tcl.C:46 tcl.C:47 tcl.C:48 tcl.C:49 tcl.C:50 tcl.C:51 tcl.C:52 tcl.C:53 tcl.C:54 tcl.C:55 tcl.C:56 tcl.C:57 tcl.C:58 tcl.C:59 tcl.C:60 tcl.C:61 tcl.C:62 tcl.C:63 tcl.C:64 tcl.C:65 tcl.C:66 tcl.C:67 tcl.C:68 tcl.C:69 tcl.C:70 tcl.C:71 tcl.C:72 tcl.C:73 tcl.C:74 tcl.C:75 tcl.C:76 tcl.C:77 tcl.C:78 tcl.C:79 tcl.C:80 tcl.C:81 tcl.C:82 tcl.C:83 tcl.C:84 tcl.C:85 tcl.C:86 tcl.C:87 tcl.C:88 tcl.C:89 |
|