ROOT logo
#include "TFile.h"
#include "TTree.h"
#include "TBrowser.h"
#include "TH2.h"
#include "TRandom.h"
#include "TCanvas.h"

// This example illustrates how to make a Tree from variables or arrays
// in a C struct. Use of C structs is strongly discouraged and one should
// use classes instead. However support for C structs is important for 
// legacy applications written in C or Fortran.
//    see tree2a.C for the same example using a class instead of a C-struct.
//
// In this example, we are mapping a C struct to one of the Geant3
// common blocks /gctrak/. In the real life, this common will be filled
// by Geant3 at each step and only the Tree Fill function should be called.
// The example emulates the Geant3 step routines.
//
// to run the example, do:
// .x tree2.C   to execute with the CINT interpreter
// .x tree2.C++ to execute with native compiler
//
//  Author: Rene Brun

const Int_t MAXMEC = 30;
//      PARAMETER (MAXMEC=30) 
//      COMMON/GCTRAK/VECT(7),GETOT,GEKIN,VOUT(7),NMEC,LMEC(MAXMEC) 
//     + ,NAMEC(MAXMEC),NSTEP ,PID,DESTEP,DESTEL,SAFETY,SLENG 
//     + ,STEP  ,SNEXT ,SFIELD,TOFG  ,GEKRAT,UPWGHT
typedef struct { 
  Float_t  vect[7]; 
  Float_t  getot; 
  Float_t  gekin; 
  Float_t  vout[7]; 
  Int_t    nmec; 
  Int_t    lmec[MAXMEC]; 
  Int_t    namec[MAXMEC]; 
  Int_t    nstep; 
  Int_t    pid; 
  Float_t  destep; 
  Float_t  destel; 
  Float_t  safety; 
  Float_t  sleng; 
  Float_t  step; 
  Float_t  snext; 
  Float_t  sfield; 
  Float_t  tofg; 
  Float_t  gekrat; 
  Float_t  upwght; 
} Gctrak_t; 


void helixStep(Float_t step, Float_t *vect, Float_t *vout)
{
  // extrapolate track in constant field
   Float_t field = 20;      //magnetic field in kilogauss
   enum Evect {kX,kY,kZ,kPX,kPY,kPZ,kPP};
   vout[kPP] = vect[kPP];
   Float_t h4    = field*2.99792e-4;
   Float_t rho   = -h4/vect[kPP];
   Float_t tet   = rho*step;
   Float_t tsint = tet*tet/6;
   Float_t sintt = 1 - tsint;
   Float_t sint  = tet*sintt;
   Float_t cos1t = tet/2;
   Float_t f1 = step*sintt;
   Float_t f2 = step*cos1t;
   Float_t f3 = step*tsint*vect[kPZ];
   Float_t f4 = -tet*cos1t;
   Float_t f5 = sint;
   Float_t f6 = tet*cos1t*vect[kPZ];
   vout[kX]   = vect[kX]  + (f1*vect[kPX] - f2*vect[kPY]);
   vout[kY]   = vect[kY]  + (f1*vect[kPY] + f2*vect[kPX]);
   vout[kZ]   = vect[kZ]  + (f1*vect[kPZ] + f3);
   vout[kPX]  = vect[kPX] + (f4*vect[kPX] - f5*vect[kPY]);
   vout[kPY]  = vect[kPY] + (f4*vect[kPY] + f5*vect[kPX]);
   vout[kPZ]  = vect[kPZ] + (f4*vect[kPZ] + f6);   
}

void tree2w()
{
   //create a Tree file tree2.root
   
   //create the file, the Tree and a few branches with 
   //a subset of gctrak
   TFile f("tree2.root","recreate");
   TTree t2("t2","a Tree with data from a fake Geant3");
   Gctrak_t gstep;
   t2.Branch("vect",gstep.vect,"vect[7]/F");
   t2.Branch("getot",&gstep.getot,"getot/F");
   t2.Branch("gekin",&gstep.gekin,"gekin/F");
   t2.Branch("nmec",&gstep.nmec,"nmec/I");
   t2.Branch("lmec",gstep.lmec,"lmec[nmec]/I");
   t2.Branch("destep",&gstep.destep,"destep/F");
   t2.Branch("pid",&gstep.pid,"pid/I");
   
   //Initialize particle parameters at first point
   Float_t px,py,pz,p,charge=0;
   Float_t vout[7];
   Float_t mass  = 0.137;
   Bool_t newParticle = kTRUE;
   gstep.step    = 0.1;
   gstep.destep  = 0;
   gstep.nmec    = 0;
   gstep.pid     = 0;
      
   //transport particles 
   for (Int_t i=0;i<10000;i++) {
      //generate a new particle if necessary
      if (newParticle) {
         px = gRandom->Gaus(0,.02);
         py = gRandom->Gaus(0,.02);
         pz = gRandom->Gaus(0,.02);
         p  = TMath::Sqrt(px*px+py*py+pz*pz);
         charge = 1; if (gRandom->Rndm() < 0.5) charge = -1;
         gstep.pid    += 1;
         gstep.vect[0] = 0;
         gstep.vect[1] = 0;
         gstep.vect[2] = 0;
         gstep.vect[3] = px/p;
         gstep.vect[4] = py/p;
         gstep.vect[5] = pz/p;
         gstep.vect[6] = p*charge;
         gstep.getot   = TMath::Sqrt(p*p + mass*mass);
         gstep.gekin   = gstep.getot - mass;
         newParticle = kFALSE;
      }
      
      // fill the Tree with current step parameters
      t2.Fill();
      
      //transport particle in magnetic field
      helixStep(gstep.step, gstep.vect, vout); //make one step
      
      //apply energy loss
      gstep.destep = gstep.step*gRandom->Gaus(0.0002,0.00001);
      gstep.gekin -= gstep.destep;
      gstep.getot   = gstep.gekin + mass;
      gstep.vect[6] = charge*TMath::Sqrt(gstep.getot*gstep.getot - mass*mass);
      gstep.vect[0] = vout[0];
      gstep.vect[1] = vout[1];
      gstep.vect[2] = vout[2];
      gstep.vect[3] = vout[3];
      gstep.vect[4] = vout[4];
      gstep.vect[5] = vout[5];
      gstep.nmec    = (Int_t)(5*gRandom->Rndm());
      for (Int_t l=0;l<gstep.nmec;l++) gstep.lmec[l] = l;
      if (gstep.gekin < 0.001)            newParticle = kTRUE;
      if (TMath::Abs(gstep.vect[2]) > 30) newParticle = kTRUE;
   }
  
   //save the Tree header. The file will be automatically closed
   //when going out of the function scope
   t2.Write();
}

void tree2r()
{
   //read the Tree generated by tree2w and fill one histogram
   //we are only interested by the destep branch.
     
   //note that we use "new" to create the TFile and TTree objects !
   //because we want to keep these objects alive when we leave 
   //this function.
   TFile *f = new TFile("tree2.root");
   TTree *t2 = (TTree*)f->Get("t2");
   static Float_t destep;
   TBranch *b_destep = t2->GetBranch("destep");
   b_destep->SetAddress(&destep);
   
   //create one histogram
   TH1F *hdestep   = new TH1F("hdestep","destep in Mev",100,1e-5,3e-5);
   
   //read only the destep branch for all entries
   Int_t nentries = (Int_t)t2->GetEntries();
   for (Int_t i=0;i<nentries;i++) {
      b_destep->GetEntry(i); 
      hdestep->Fill(destep);
   }
  
   //we do not close the file. 
   //We want to keep the generated histograms
   //We fill a 3-d scatter plot with the particle step coordinates
   TCanvas *c1 = new TCanvas("c1","c1",600,800);
   c1->SetFillColor(42);
   c1->Divide(1,2);
   c1->cd(1);
   hdestep->SetFillColor(45);
   hdestep->Fit("gaus");
   c1->cd(2);
   gPad->SetFillColor(37);
   t2->SetMarkerColor(kRed);
   t2->Draw("vect[0]:vect[1]:vect[2]");
   if (gROOT->IsBatch()) return;
   
   // invoke the x3d viewer
   gPad->GetViewer3D("x3d");
}   

void tree2() {
   tree2w();
   tree2r();
}
 tree2.C:1
 tree2.C:2
 tree2.C:3
 tree2.C:4
 tree2.C:5
 tree2.C:6
 tree2.C:7
 tree2.C:8
 tree2.C:9
 tree2.C:10
 tree2.C:11
 tree2.C:12
 tree2.C:13
 tree2.C:14
 tree2.C:15
 tree2.C:16
 tree2.C:17
 tree2.C:18
 tree2.C:19
 tree2.C:20
 tree2.C:21
 tree2.C:22
 tree2.C:23
 tree2.C:24
 tree2.C:25
 tree2.C:26
 tree2.C:27
 tree2.C:28
 tree2.C:29
 tree2.C:30
 tree2.C:31
 tree2.C:32
 tree2.C:33
 tree2.C:34
 tree2.C:35
 tree2.C:36
 tree2.C:37
 tree2.C:38
 tree2.C:39
 tree2.C:40
 tree2.C:41
 tree2.C:42
 tree2.C:43
 tree2.C:44
 tree2.C:45
 tree2.C:46
 tree2.C:47
 tree2.C:48
 tree2.C:49
 tree2.C:50
 tree2.C:51
 tree2.C:52
 tree2.C:53
 tree2.C:54
 tree2.C:55
 tree2.C:56
 tree2.C:57
 tree2.C:58
 tree2.C:59
 tree2.C:60
 tree2.C:61
 tree2.C:62
 tree2.C:63
 tree2.C:64
 tree2.C:65
 tree2.C:66
 tree2.C:67
 tree2.C:68
 tree2.C:69
 tree2.C:70
 tree2.C:71
 tree2.C:72
 tree2.C:73
 tree2.C:74
 tree2.C:75
 tree2.C:76
 tree2.C:77
 tree2.C:78
 tree2.C:79
 tree2.C:80
 tree2.C:81
 tree2.C:82
 tree2.C:83
 tree2.C:84
 tree2.C:85
 tree2.C:86
 tree2.C:87
 tree2.C:88
 tree2.C:89
 tree2.C:90
 tree2.C:91
 tree2.C:92
 tree2.C:93
 tree2.C:94
 tree2.C:95
 tree2.C:96
 tree2.C:97
 tree2.C:98
 tree2.C:99
 tree2.C:100
 tree2.C:101
 tree2.C:102
 tree2.C:103
 tree2.C:104
 tree2.C:105
 tree2.C:106
 tree2.C:107
 tree2.C:108
 tree2.C:109
 tree2.C:110
 tree2.C:111
 tree2.C:112
 tree2.C:113
 tree2.C:114
 tree2.C:115
 tree2.C:116
 tree2.C:117
 tree2.C:118
 tree2.C:119
 tree2.C:120
 tree2.C:121
 tree2.C:122
 tree2.C:123
 tree2.C:124
 tree2.C:125
 tree2.C:126
 tree2.C:127
 tree2.C:128
 tree2.C:129
 tree2.C:130
 tree2.C:131
 tree2.C:132
 tree2.C:133
 tree2.C:134
 tree2.C:135
 tree2.C:136
 tree2.C:137
 tree2.C:138
 tree2.C:139
 tree2.C:140
 tree2.C:141
 tree2.C:142
 tree2.C:143
 tree2.C:144
 tree2.C:145
 tree2.C:146
 tree2.C:147
 tree2.C:148
 tree2.C:149
 tree2.C:150
 tree2.C:151
 tree2.C:152
 tree2.C:153
 tree2.C:154
 tree2.C:155
 tree2.C:156
 tree2.C:157
 tree2.C:158
 tree2.C:159
 tree2.C:160
 tree2.C:161
 tree2.C:162
 tree2.C:163
 tree2.C:164
 tree2.C:165
 tree2.C:166
 tree2.C:167
 tree2.C:168
 tree2.C:169
 tree2.C:170
 tree2.C:171
 tree2.C:172
 tree2.C:173
 tree2.C:174
 tree2.C:175
 tree2.C:176
 tree2.C:177
 tree2.C:178
 tree2.C:179
 tree2.C:180
 tree2.C:181
 tree2.C:182
 tree2.C:183
 tree2.C:184
 tree2.C:185
 tree2.C:186
 tree2.C:187
 tree2.C:188
 tree2.C:189
 tree2.C:190
 tree2.C:191
 tree2.C:192
 tree2.C:193
 tree2.C:194
 tree2.C:195
 tree2.C:196
 tree2.C:197
 tree2.C:198
 tree2.C:199
 tree2.C:200
 tree2.C:201
 tree2.C:202
 tree2.C:203
 tree2.C:204