Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
tcl.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_tree
3/// \notebook
4/// How to write a TClonesArray to a TTree
5///
6/// The following tests can be run
7/// Interactive tests
8/// ~~~
9/// Root > .x tcl.C //no-split interpreted
10/// Root > .x tcl.C(1) //split interpreted
11/// Root > .x tcl.C++ //no-split compiled
12/// Root > .x tcl.C++(1) //split compiled
13/// ~~~
14/// Batch tests: same as above but with no graphics
15/// ~~~{.bash}
16/// root -b -q tcl.C
17/// root -b -q tcl.C++
18/// root -b -q "tcl.C(1)"
19/// root -b -q "tcl.C++(1)"
20/// ~~~
21/// \macro_code
22///
23/// \author Rene Brun
24
25#include "TFile.h"
26#include "TClonesArray.h"
27#include "TH2.h"
28#include "TLine.h"
29#include "TTree.h"
30#include "TBenchmark.h"
31#include "TRandom.h"
32
33void tclwrite(Int_t split)
34{
35// Generate a Tree with a TClonesArray
36// The array can be split or not
37 TFile f("tcl.root","recreate");
38 f.SetCompressionLevel(1); //try level 2 also
39 TTree T("T","test tcl");
40 TClonesArray *arr = new TClonesArray("TLine");
41 TClonesArray &ar = *arr;
42 T.Branch("tcl",&arr,256000,split);
43 //By default a TClonesArray is created with its BypassStreamer bit set.
44 //However, because TLine has a custom Streamer, this bit was reset
45 //by TTree::Branch above. We set again this bit because the current
46 //version of TLine uses the automatic Streamer.
47 //BypassingStreamer saves space and time.
48 arr->BypassStreamer();
49 for (Int_t ev=0;ev<10000;ev++) {
50 ar.Clear();
51 Int_t nlines = Int_t(gRandom->Gaus(50,10));
52 if(nlines < 0) nlines = 1;
53 for (Int_t i=0;i<nlines;i++) {
55 Float_t y1 = gRandom->Rndm();
57 Float_t y2 = gRandom->Rndm();
58 new(ar[i]) TLine(x1,y1,x2,y2);
59 }
60 T.Fill();
61 }
62 T.Print();
63 T.Write();
64}
65
66void tclread()
67{
68// read file generated by tclwrite
69// loop on all entries.
70// histogram center of lines
71 TFile *f = new TFile("tcl.root");
72 TTree *T = (TTree*)f->Get("T");
73 TH2F *h2 = new TH2F("h2","center of lines",40,0,1,40,0,1);
74
75 TClonesArray *arr = new TClonesArray("TLine");
76 T->GetBranch("tcl")->SetAutoDelete(kFALSE);
77 T->SetBranchAddress("tcl",&arr);
78 Long64_t nentries = T->GetEntries();
79 for (Long64_t ev=0;ev<nentries;ev++) {
80 arr->Clear();
81 T->GetEntry(ev);
82 Int_t nlines = arr->GetEntriesFast();
83 for (Int_t i=0;i<nlines;i++) {
84 TLine *line = (TLine*)arr->At(i);
85 h2->Fill(0.5*(line->GetX1()+line->GetX2()), 0.5*(line->GetY1()+line->GetY2()));
86 }
87 }
88 h2->Draw("lego");
89}
90
91void tcl(Int_t split=0)
92{
93 gBenchmark->Start("tcl");
94 tclwrite(split);
95 tclread();
96 gBenchmark->Show("tcl");
97}
#define f(i)
Definition RSha256.hxx:104
static const double x2[5]
static const double x1[5]
int Int_t
Definition RtypesCore.h:45
const Bool_t kFALSE
Definition RtypesCore.h:101
long long Long64_t
Definition RtypesCore.h:80
float Float_t
Definition RtypesCore.h:57
R__EXTERN TBenchmark * gBenchmark
Definition TBenchmark.h:59
int nentries
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
virtual void Start(const char *name)
Starts Benchmark with the specified name.
virtual void Show(const char *name)
Stops Benchmark name and Prints results.
An array of clone (identical) objects.
void BypassStreamer(Bool_t bypass=kTRUE)
When the kBypassStreamer bit is set, the automatically generated Streamer can call directly TClass::W...
virtual void Clear(Option_t *option="")
Clear the clones array.
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition TFile.h:54
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition TH1.cxx:3074
2-D histogram with a float per channel (see TH1 documentation)}
Definition TH2.h:251
Int_t Fill(Double_t)
Invalid Fill method.
Definition TH2.cxx:358
Use the TLine constructor to create a simple line.
Definition TLine.h:22
Double_t GetY1() const
Definition TLine.h:52
Double_t GetX2() const
Definition TLine.h:51
Double_t GetX1() const
Definition TLine.h:50
Double_t GetY2() const
Definition TLine.h:53
Int_t GetEntriesFast() const
Definition TObjArray.h:58
TObject * At(Int_t idx) const
Definition TObjArray.h:164
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition TRandom.cxx:274
virtual Double_t Rndm()
Machine independent random number generator.
Definition TRandom.cxx:552
A TTree represents a columnar dataset.
Definition TTree.h:79
TLine * line
double T(double x)