Logo ROOT   6.10/09
Reference Guide
tree2.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_tree
3 /// \notebook
4 /// This example illustrates how to make a Tree from variables or arrays
5 /// in a C struct. **Use of C structs is strongly discouraged and one should
6 /// use classes instead**. However support for C structs is important for
7 /// legacy applications written in C or Fortran.
8 /// See tree2a.C for the same example using a class instead of a C-struct.
9 ///
10 /// In this example, we are mapping a C struct to one of the Geant3
11 /// common blocks /gctrak/. In the real life, this common will be filled
12 /// by Geant3 at each step and only the Tree Fill function should be called.
13 /// The example emulates the Geant3 step routines.
14 ///
15 /// to run the example, do:
16 /// ~~~
17 /// .x tree2.C to execute with the Cling interpreter
18 /// .x tree2.C++ to execute with native compiler
19 /// ~~~
20 /// \macro_code
21 ///
22 /// \author Rene Brun
23 
24 #include "TFile.h"
25 #include "TTree.h"
26 #include "TBrowser.h"
27 #include "TH2.h"
28 #include "TRandom.h"
29 #include "TCanvas.h"
30 #include "TMath.h"
31 #include "TROOT.h"
32 
33 const Int_t MAXMEC = 30;
34 
35 typedef struct {
36  Float_t vect[7];
37  Float_t getot;
38  Float_t gekin;
39  Float_t vout[7];
40  Int_t nmec;
41  Int_t lmec[MAXMEC];
42  Int_t namec[MAXMEC];
43  Int_t nstep;
44  Int_t pid;
45  Float_t destep;
46  Float_t destel;
47  Float_t safety;
48  Float_t sleng;
49  Float_t step;
50  Float_t snext;
51  Float_t sfield;
52  Float_t tofg;
53  Float_t gekrat;
54  Float_t upwght;
55 } Gctrak_t;
56 
57 
58 void helixStep(Float_t step, Float_t *vect, Float_t *vout)
59 {
60  // extrapolate track in constant field
61  Float_t field = 20; //magnetic field in kilogauss
62  enum Evect {kX,kY,kZ,kPX,kPY,kPZ,kPP};
63  vout[kPP] = vect[kPP];
64  Float_t h4 = field*2.99792e-4;
65  Float_t rho = -h4/vect[kPP];
66  Float_t tet = rho*step;
67  Float_t tsint = tet*tet/6;
68  Float_t sintt = 1 - tsint;
69  Float_t sint = tet*sintt;
70  Float_t cos1t = tet/2;
71  Float_t f1 = step*sintt;
72  Float_t f2 = step*cos1t;
73  Float_t f3 = step*tsint*vect[kPZ];
74  Float_t f4 = -tet*cos1t;
75  Float_t f5 = sint;
76  Float_t f6 = tet*cos1t*vect[kPZ];
77  vout[kX] = vect[kX] + (f1*vect[kPX] - f2*vect[kPY]);
78  vout[kY] = vect[kY] + (f1*vect[kPY] + f2*vect[kPX]);
79  vout[kZ] = vect[kZ] + (f1*vect[kPZ] + f3);
80  vout[kPX] = vect[kPX] + (f4*vect[kPX] - f5*vect[kPY]);
81  vout[kPY] = vect[kPY] + (f4*vect[kPY] + f5*vect[kPX]);
82  vout[kPZ] = vect[kPZ] + (f4*vect[kPZ] + f6);
83 }
84 
85 void tree2w()
86 {
87  //create a Tree file tree2.root
88 
89  //create the file, the Tree and a few branches with
90  //a subset of gctrak
91  TFile f("tree2.root","recreate");
92  TTree t2("t2","a Tree with data from a fake Geant3");
93  Gctrak_t gstep;
94  t2.Branch("vect",gstep.vect,"vect[7]/F");
95  t2.Branch("getot",&gstep.getot,"getot/F");
96  t2.Branch("gekin",&gstep.gekin,"gekin/F");
97  t2.Branch("nmec",&gstep.nmec,"nmec/I");
98  t2.Branch("lmec",gstep.lmec,"lmec[nmec]/I");
99  t2.Branch("destep",&gstep.destep,"destep/F");
100  t2.Branch("pid",&gstep.pid,"pid/I");
101 
102  //Initialize particle parameters at first point
103  Float_t px,py,pz,p,charge=0;
104  Float_t vout[7];
105  Float_t mass = 0.137;
106  Bool_t newParticle = kTRUE;
107  gstep.step = 0.1;
108  gstep.destep = 0;
109  gstep.nmec = 0;
110  gstep.pid = 0;
111 
112  //transport particles
113  for (Int_t i=0;i<10000;i++) {
114  //generate a new particle if necessary
115  if (newParticle) {
116  px = gRandom->Gaus(0,.02);
117  py = gRandom->Gaus(0,.02);
118  pz = gRandom->Gaus(0,.02);
119  p = TMath::Sqrt(px*px+py*py+pz*pz);
120  charge = 1; if (gRandom->Rndm() < 0.5) charge = -1;
121  gstep.pid += 1;
122  gstep.vect[0] = 0;
123  gstep.vect[1] = 0;
124  gstep.vect[2] = 0;
125  gstep.vect[3] = px/p;
126  gstep.vect[4] = py/p;
127  gstep.vect[5] = pz/p;
128  gstep.vect[6] = p*charge;
129  gstep.getot = TMath::Sqrt(p*p + mass*mass);
130  gstep.gekin = gstep.getot - mass;
131  newParticle = kFALSE;
132  }
133 
134  // fill the Tree with current step parameters
135  t2.Fill();
136 
137  //transport particle in magnetic field
138  helixStep(gstep.step, gstep.vect, vout); //make one step
139 
140  //apply energy loss
141  gstep.destep = gstep.step*gRandom->Gaus(0.0002,0.00001);
142  gstep.gekin -= gstep.destep;
143  gstep.getot = gstep.gekin + mass;
144  gstep.vect[6] = charge*TMath::Sqrt(gstep.getot*gstep.getot - mass*mass);
145  gstep.vect[0] = vout[0];
146  gstep.vect[1] = vout[1];
147  gstep.vect[2] = vout[2];
148  gstep.vect[3] = vout[3];
149  gstep.vect[4] = vout[4];
150  gstep.vect[5] = vout[5];
151  gstep.nmec = (Int_t)(5*gRandom->Rndm());
152  for (Int_t l=0;l<gstep.nmec;l++) gstep.lmec[l] = l;
153  if (gstep.gekin < 0.001) newParticle = kTRUE;
154  if (TMath::Abs(gstep.vect[2]) > 30) newParticle = kTRUE;
155  }
156 
157  //save the Tree header. The file will be automatically closed
158  //when going out of the function scope
159  t2.Write();
160 }
161 
162 void tree2r()
163 {
164  //read the Tree generated by tree2w and fill one histogram
165  //we are only interested by the destep branch.
166 
167  //note that we use "new" to create the TFile and TTree objects !
168  //because we want to keep these objects alive when we leave
169  //this function.
170  TFile *f = new TFile("tree2.root");
171  TTree *t2 = (TTree*)f->Get("t2");
172  static Float_t destep;
173  TBranch *b_destep = t2->GetBranch("destep");
174  b_destep->SetAddress(&destep);
175 
176  //create one histogram
177  TH1F *hdestep = new TH1F("hdestep","destep in Mev",100,1e-5,3e-5);
178 
179  //read only the destep branch for all entries
180  Long64_t nentries = t2->GetEntries();
181  for (Long64_t i=0;i<nentries;i++) {
182  b_destep->GetEntry(i);
183  hdestep->Fill(destep);
184  }
185 
186  //we do not close the file.
187  //We want to keep the generated histograms
188  //We fill a 3-d scatter plot with the particle step coordinates
189  TCanvas *c1 = new TCanvas("c1","c1",600,800);
190  c1->SetFillColor(42);
191  c1->Divide(1,2);
192  c1->cd(1);
193  hdestep->SetFillColor(45);
194  hdestep->Fit("gaus");
195  c1->cd(2);
196  gPad->SetFillColor(37);
197  t2->SetMarkerColor(kRed);
198  t2->Draw("vect[0]:vect[1]:vect[2]");
199  if (gROOT->IsBatch()) return;
200 
201  // invoke the x3d viewer
202  gPad->GetViewer3D("x3d");
203 }
204 
205 void tree2() {
206  tree2w();
207  tree2r();
208 }
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3126
virtual void SetAddress(void *add)
Set address of this branch.
Definition: TBranch.cxx:2148
#define snext(osub1, osub2)
Definition: triangle.c:1167
long long Long64_t
Definition: RtypesCore.h:69
float Float_t
Definition: RtypesCore.h:53
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:235
return c1
Definition: legend1.C:41
Definition: Rtypes.h:56
THist< 1, float, THistStatContent, THistStatUncertainty > TH1F
Definition: THist.hxx:311
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format...
Definition: TFile.h:46
TVirtualPad * cd(Int_t subpadnumber=0)
Set current canvas & pad.
Definition: TCanvas.cxx:679
virtual TObject * Get(const char *namecycle)
Return pointer to object identified by namecycle.
#define gROOT
Definition: TROOT.h:375
tomato 1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:551
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
Short_t Abs(Short_t d)
Definition: TMathBase.h:108
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
virtual Double_t Rndm()
Machine independent random number generator.
Definition: TRandom.cxx:512
virtual TBranch * GetBranch(const char *name)
Return pointer to the branch with the given name in this tree or its friends.
Definition: TTree.cxx:4979
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
TLine * l
Definition: textangle.C:4
virtual Int_t GetEntry(Long64_t entry=0, Int_t getall=0)
Read all leaves of entry and return total number of bytes read.
Definition: TBranch.cxx:1291
R__EXTERN TRandom * gRandom
Definition: TRandom.h:62
const Bool_t kFALSE
Definition: RtypesCore.h:92
The Canvas class.
Definition: TCanvas.h:31
double f(double x)
virtual void Draw(Option_t *opt)
Default Draw method for all objects.
Definition: TTree.h:355
int nentries
Definition: THbookFile.cxx:89
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
virtual Long64_t GetEntries() const
Definition: TTree.h:381
virtual void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0)
Automatic pad generation by division.
Definition: TPad.cxx:1135
double f2(const double *x)
TF1 * f1
Definition: legend1.C:11
#define gPad
Definition: TVirtualPad.h:284
A TTree object has a header with a name and a title.
Definition: TTree.h:78
Double_t Sqrt(Double_t x)
Definition: TMath.h:591
A TTree is a list of TBranches.
Definition: TBranch.h:57
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
Definition: TH1.cxx:3564
const Bool_t kTRUE
Definition: RtypesCore.h:91