ROOT  6.06/09
Reference Guide
TGenerator.cxx
Go to the documentation of this file.
1 // @(#)root/eg:$Id$
2 // Author: Ola Nordmann 21/09/95
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 //////////////////////////////////////////////////////////////////////////
13 // //
14 // TGenerator //
15 // //
16 // Is an base class, that defines the interface of ROOT to various //
17 // event generators. Every event generator should inherit from //
18 // TGenerator or its subclasses. //
19 // //
20 // Derived class can overload the member function GenerateEvent //
21 // to do the actual event generation (e.g., call PYEVNT or similar). //
22 // //
23 // The derived class should overload the member function //
24 // ImportParticles (both types) to read the internal storage of the //
25 // generated event into either the internal TObjArray or the passed //
26 // TClonesArray of TParticles. //
27 // //
28 // If the generator code stores event data in the /HEPEVT/ common block //
29 // Then the default implementation of ImportParticles should suffice. //
30 // The common block /HEPEVT/ is structed like //
31 // //
32 // /* C */ //
33 // typedef struct { //
34 // Int_t nevhep; // Event number //
35 // Int_t nhep; // # of particles //
36 // Int_t isthep[4000]; // Status flag of i'th particle //
37 // Int_t idhep[4000]; // PDG # of particle //
38 // Int_t jmohep[4000][2]; // 1st & 2nd mother particle # //
39 // Int_t jdahep[4000][2]; // 1st & 2nd daughter particle # //
40 // Double_t phep[4000][5]; // 4-momentum and 1 word //
41 // Double_t vhep[4000][4]; // 4-position of production //
42 // } HEPEVT_DEF; //
43 // //
44 // //
45 // C Fortran //
46 // COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(4000),IDHEP(4000), //
47 // + JMOHEP(2,4000),JDAHEP(2,4000),PHEP(5,4000),VHEP(4,4000) //
48 // INTEGER NEVHEP,NHEP,ISTHEP,IDHEP,JMOHEP,JDAHEP //
49 // DOUBLE PRECISION PHEP,VHEP //
50 // //
51 // The generic member functions SetParameter and GetParameter can be //
52 // overloaded to set and get parameters of the event generator. //
53 // //
54 // Note, if the derived class interfaces a (set of) Fortran common //
55 // blocks (like TPythia, TVenus does), one better make the derived //
56 // class a singleton. That is, something like //
57 // //
58 // class MyGenerator : public TGenerator //
59 // { //
60 // public: //
61 // static MyGenerator* Instance() //
62 // { //
63 // if (!fgInstance) fgInstance = new MyGenerator; //
64 // return fgInstance; //
65 // } //
66 // void GenerateEvent() { ... } //
67 // void ImportParticles(TClonesArray* a, Option_t opt="") {...} //
68 // Int_t ImportParticles(Option_t opt="") { ... } //
69 // Int_t SetParameter(const char* name, Double_t val) { ... } //
70 // Double_t GetParameter(const char* name) { ... } //
71 // virtual ~MyGenerator() { ... } //
72 // protected: //
73 // MyGenerator() { ... } //
74 // MyGenerator(const MyGenerator& o) { ... } //
75 // MyGenerator& operator=(const MyGenerator& o) { ... } //
76 // static MyGenerator* fgInstance; //
77 // ClassDef(MyGenerator,0); //
78 // }; //
79 // //
80 // Having multiple objects accessing the same common blocks is not //
81 // safe. //
82 // //
83 // concrete TGenerator classes can be loaded in scripts and subseqent- //
84 // ly used in compiled code: //
85 // //
86 // // MyRun.h //
87 // class MyRun : public TObject //
88 // { //
89 // public: //
90 // static MyRun* Instance() { ... } //
91 // void SetGenerator(TGenerator* g) { fGenerator = g; } //
92 // void Run(Int_t n, Option_t* option="") //
93 // { //
94 // TFile* file = TFile::Open("file.root","RECREATE"); //
95 // TTree* tree = new TTree("T","T"); //
96 // TClonesArray* p = new TClonesArray("TParticles"); //
97 // tree->Branch("particles", &p); //
98 // for (Int_t event = 0; event < n; event++) { //
99 // fGenerator->GenerateEvent(); //
100 // fGenerator->ImportParticles(p,option); //
101 // tree->Fill(); //
102 // } //
103 // file->Write(); //
104 // file->Close(); //
105 // } //
106 // ... //
107 // protected: //
108 // TGenerator* fGenerator; //
109 // ClassDef(MyRun,0); //
110 // }; //
111 // //
112 // // Config.C //
113 // void Config() //
114 // { //
115 // MyRun* run = MyRun::Instance(); //
116 // run->SetGenerator(MyGenerator::Instance()); //
117 // } //
118 // //
119 // // main.cxx //
120 // int //
121 // main(int argc, char** argv) //
122 // { //
123 // TApplication app("", 0, 0); //
124 // gSystem->ProcessLine(".x Config.C"); //
125 // MyRun::Instance()->Run(10); //
126 // return 0; //
127 // } //
128 // //
129 // This is especially useful for example with TVirtualMC or similar. //
130 // //
131 //////////////////////////////////////////////////////////////////////////
132 
133 #include "TROOT.h"
134 #include "TGenerator.h"
135 #include "TDatabasePDG.h"
136 #include "TParticlePDG.h"
137 #include "TParticle.h"
138 #include "TObjArray.h"
139 #include "Hepevt.h"
140 #include "TVirtualPad.h"
141 #include "TView.h"
142 #include "TText.h"
143 #include "TPaveText.h"
144 #include "TClonesArray.h"
145 #include "Riostream.h"
146 
147 
149 
150 ////////////////////////////////////////////////////////////////////////////////
151 /// Event generator default constructor
152 ///
153 
154 TGenerator::TGenerator(const char *name,const char *title): TNamed(name,title)
155 {
156  // Initialize particles table
158  //TDatabasePDG *pdg = TDatabasePDG::Instance();
159  //if (!pdg->ParticleList()) pdg->Init();
160 
161  fPtCut = 0;
162  fShowNeutrons = kTRUE;
163  fParticles = new TObjArray(10000);
164 }
165 
166 ////////////////////////////////////////////////////////////////////////////////
167 /// Event generator default destructor
168 ///
169 
171 {
172  //do nothing
173  if (fParticles) {
174  fParticles->Delete();
175  delete fParticles;
176  fParticles = 0;
177  }
178 }
179 
180 ////////////////////////////////////////////////////////////////////////////////
181 /// must be implemented in concrete class (see eg TPythia6)
182 
184 {
185 }
186 
187 ////////////////////////////////////////////////////////////////////////////////
188 ///
189 /// It reads the /HEPEVT/ common block which has been filled by the
190 /// GenerateEvent method. If the event generator does not use the
191 /// HEPEVT common block, This routine has to be overloaded by the
192 /// subclasses.
193 ///
194 /// The default action is to store only the stable particles (ISTHEP =
195 /// 1) This can be demanded explicitly by setting the option = "Final"
196 /// If the option = "All", all the particles are stored.
197 ///
198 
200 {
201  fParticles->Clear();
202  Int_t numpart = HEPEVT.nhep;
203  if (!strcmp(option,"") || !strcmp(option,"Final")) {
204  for (Int_t i = 0; i<numpart; i++) {
205  if (HEPEVT.isthep[i] == 1) {
206 //
207 // Use the common block values for the TParticle constructor
208 //
209  TParticle *p = new TParticle(
210  HEPEVT.idhep[i],
211  HEPEVT.isthep[i],
212  HEPEVT.jmohep[i][0]-1,
213  HEPEVT.jmohep[i][1]-1,
214  HEPEVT.jdahep[i][0]-1,
215  HEPEVT.jdahep[i][1]-1,
216  HEPEVT.phep[i][0],
217  HEPEVT.phep[i][1],
218  HEPEVT.phep[i][2],
219  HEPEVT.phep[i][3],
220  HEPEVT.vhep[i][0],
221  HEPEVT.vhep[i][1],
222  HEPEVT.vhep[i][2],
223  HEPEVT.vhep[i][3]);
224  fParticles->Add(p);
225  }
226  }
227  } else if (!strcmp(option,"All")) {
228  for (Int_t i = 0; i<numpart; i++) {
229  TParticle *p = new TParticle(
230  HEPEVT.idhep[i],
231  HEPEVT.isthep[i],
232  HEPEVT.jmohep[i][0]-1,
233  HEPEVT.jmohep[i][1]-1,
234  HEPEVT.jdahep[i][0]-1,
235  HEPEVT.jdahep[i][1]-1,
236  HEPEVT.phep[i][0],
237  HEPEVT.phep[i][1],
238  HEPEVT.phep[i][2],
239  HEPEVT.phep[i][3],
240  HEPEVT.vhep[i][0],
241  HEPEVT.vhep[i][1],
242  HEPEVT.vhep[i][2],
243  HEPEVT.vhep[i][3]);
244  fParticles->Add(p);
245  }
246  }
247  return fParticles;
248 }
249 
250 ////////////////////////////////////////////////////////////////////////////////
251 ///
252 /// It reads the /HEPEVT/ common block which has been filled by the
253 /// GenerateEvent method. If the event generator does not use the
254 /// HEPEVT common block, This routine has to be overloaded by the
255 /// subclasses.
256 ///
257 /// The function loops on the generated particles and store them in
258 /// the TClonesArray pointed by the argument particles. The default
259 /// action is to store only the stable particles (ISTHEP = 1) This can
260 /// be demanded explicitly by setting the option = "Final" If the
261 /// option = "All", all the particles are stored.
262 ///
263 
265 {
266  if (particles == 0) return 0;
267  TClonesArray &clonesParticles = *particles;
268  clonesParticles.Clear();
269  Int_t numpart = HEPEVT.nhep;
270  if (!strcmp(option,"") || !strcmp(option,"Final")) {
271  for (Int_t i = 0; i<numpart; i++) {
272  if (HEPEVT.isthep[i] == 1) {
273 //
274 // Use the common block values for the TParticle constructor
275 //
276  new(clonesParticles[i]) TParticle(
277  HEPEVT.idhep[i],
278  HEPEVT.isthep[i],
279  HEPEVT.jmohep[i][0]-1,
280  HEPEVT.jmohep[i][1]-1,
281  HEPEVT.jdahep[i][0]-1,
282  HEPEVT.jdahep[i][1]-1,
283  HEPEVT.phep[i][0],
284  HEPEVT.phep[i][1],
285  HEPEVT.phep[i][2],
286  HEPEVT.phep[i][3],
287  HEPEVT.vhep[i][0],
288  HEPEVT.vhep[i][1],
289  HEPEVT.vhep[i][2],
290  HEPEVT.vhep[i][3]);
291  }
292  }
293  } else if (!strcmp(option,"All")) {
294  for (Int_t i = 0; i<numpart; i++) {
295  new(clonesParticles[i]) TParticle(
296  HEPEVT.idhep[i],
297  HEPEVT.isthep[i],
298  HEPEVT.jmohep[i][0]-1,
299  HEPEVT.jmohep[i][1]-1,
300  HEPEVT.jdahep[i][0]-1,
301  HEPEVT.jdahep[i][1]-1,
302  HEPEVT.phep[i][0],
303  HEPEVT.phep[i][1],
304  HEPEVT.phep[i][2],
305  HEPEVT.phep[i][3],
306  HEPEVT.vhep[i][0],
307  HEPEVT.vhep[i][1],
308  HEPEVT.vhep[i][2],
309  HEPEVT.vhep[i][3]);
310  }
311  }
312  return numpart;
313 }
314 
315 ////////////////////////////////////////////////////////////////////////////////
316 ///browse generator
317 
319 {
320  Draw();
321  gPad->Update();
322 }
323 
324 ////////////////////////////////////////////////////////////////////////////////
325 ///*-*-*-*-*-*-*-*Compute distance from point px,py to objects in event*-*-*-*
326 ///*-* =====================================================
327 ///*-*
328 
330 {
331  const Int_t big = 9999;
332  const Int_t inview = 0;
333  Int_t dist = big;
334  if (px > 50 && py > 50) dist = inview;
335  return dist;
336 }
337 
338 ////////////////////////////////////////////////////////////////////////////////
339 ///
340 /// Insert one event in the pad list
341 ///
342 
344 {
345  // Create a default canvas if a canvas does not exist
346  if (!gPad) {
347  gROOT->MakeDefCanvas();
348  if (gPad->GetVirtCanvas())
349  gPad->GetVirtCanvas()->SetFillColor(13);
350  }
351 
352  static Float_t rbox = 1000;
353  Float_t rmin[3],rmax[3];
354  TView *view = gPad->GetView();
355  if (!strstr(option,"same")) {
356  if (view) { view->GetRange(rmin,rmax); rbox = rmax[2];}
357  gPad->Clear();
358  }
359 
360  AppendPad(option);
361 
362  view = gPad->GetView();
363  // compute 3D view
364  if (view) {
365  view->GetRange(rmin,rmax);
366  rbox = rmax[2];
367  } else {
368  view = TView::CreateView(1,0,0);
369  if (view) view->SetRange(-rbox,-rbox,-rbox, rbox,rbox,rbox );
370  }
371  const Int_t kColorProton = 4;
372  const Int_t kColorNeutron = 5;
373  const Int_t kColorAntiProton= 3;
374  const Int_t kColorPionPlus = 6;
375  const Int_t kColorPionMinus = 2;
376  const Int_t kColorKaons = 7;
377  const Int_t kColorElectrons = 0;
378  const Int_t kColorGamma = 18;
379 
380  Int_t nProtons = 0;
381  Int_t nNeutrons = 0;
382  Int_t nAntiProtons= 0;
383  Int_t nPionPlus = 0;
384  Int_t nPionMinus = 0;
385  Int_t nKaons = 0;
386  Int_t nElectrons = 0;
387  Int_t nGammas = 0;
388 
389  Int_t ntracks = fParticles->GetEntriesFast();
390  Int_t i,lwidth,color,lstyle;
391  TParticlePDG *ap;
392  TParticle *p;
393  const char *name;
394  Double_t etot,vx,vy,vz;
395  Int_t ninvol = 0;
396  for (i=0;i<ntracks;i++) {
397  p = (TParticle*)fParticles->UncheckedAt(i);
398  if(!p) continue;
399  ap = (TParticlePDG*)p->GetPDG();
400  vx = p->Vx();
401  vy = p->Vy();
402  vz = p->Vz();
403  if (vx*vx+vy*vy+vz*vz > rbox*rbox) continue;
404  Float_t pt = p->Pt();
405  if (pt < fPtCut) continue;
406  etot = p->Energy();
407  if (etot > 0.1) lwidth = Int_t(6*TMath::Log10(etot));
408  else lwidth = 1;
409  if (lwidth < 1) lwidth = 1;
410  lstyle = 1;
411  color = 0;
412  name = ap->GetName();
413  if (!strcmp(name,"n")) { if (!fShowNeutrons) continue;
414  color = kColorNeutron; nNeutrons++;}
415  if (!strcmp(name,"p")) { color = kColorProton; nProtons++;}
416  if (!strcmp(name,"p bar")) { color = kColorAntiProton; nAntiProtons++;}
417  if (!strcmp(name,"pi+")) { color = kColorPionPlus; nPionPlus++;}
418  if (!strcmp(name,"pi-")) { color = kColorPionMinus; nPionMinus++;}
419  if (!strcmp(name,"e+")) { color = kColorElectrons; nElectrons++;}
420  if (!strcmp(name,"e-")) { color = kColorElectrons; nElectrons++;}
421  if (!strcmp(name,"gamma")) { color = kColorGamma; nGammas++; lstyle = 3; }
422  if ( strstr(name,"K")) { color = kColorKaons; nKaons++;}
423  p->SetLineColor(color);
424  p->SetLineStyle(lstyle);
425  p->SetLineWidth(lwidth);
426  p->AppendPad();
427  ninvol++;
428  }
429 
430  // event title
431  TPaveText *pt = new TPaveText(-0.94,0.85,-0.25,0.98,"br");
432  pt->AddText((char*)GetName());
433  pt->AddText((char*)GetTitle());
434  pt->SetFillColor(42);
435  pt->Draw();
436 
437  // Annotate color codes
438  Int_t tcolor = 5;
439  if (gPad->GetFillColor() == 10) tcolor = 4;
440  TText *text = new TText(-0.95,-0.47,"Particles");
441  text->SetTextAlign(12);
442  text->SetTextSize(0.025);
443  text->SetTextColor(tcolor);
444  text->Draw();
445  text->SetTextColor(kColorGamma); text->DrawText(-0.95,-0.52,"(on screen)");
446  text->SetTextColor(kColorGamma); text->DrawText(-0.95,-0.57,"Gamma");
447  text->SetTextColor(kColorProton); text->DrawText(-0.95,-0.62,"Proton");
448  text->SetTextColor(kColorNeutron); text->DrawText(-0.95,-0.67,"Neutron");
449  text->SetTextColor(kColorAntiProton); text->DrawText(-0.95,-0.72,"AntiProton");
450  text->SetTextColor(kColorPionPlus); text->DrawText(-0.95,-0.77,"Pion +");
451  text->SetTextColor(kColorPionMinus); text->DrawText(-0.95,-0.82,"Pion -");
452  text->SetTextColor(kColorKaons); text->DrawText(-0.95,-0.87,"Kaons");
453  text->SetTextColor(kColorElectrons); text->DrawText(-0.95,-0.92,"Electrons,etc.");
454 
455  text->SetTextColor(tcolor);
456  text->SetTextAlign(32);
457  char tcount[32];
458  snprintf(tcount,12,"%d",ntracks); text->DrawText(-0.55,-0.47,tcount);
459  snprintf(tcount,12,"%d",ninvol); text->DrawText(-0.55,-0.52,tcount);
460  snprintf(tcount,12,"%d",nGammas); text->DrawText(-0.55,-0.57,tcount);
461  snprintf(tcount,12,"%d",nProtons); text->DrawText(-0.55,-0.62,tcount);
462  snprintf(tcount,12,"%d",nNeutrons); text->DrawText(-0.55,-0.67,tcount);
463  snprintf(tcount,12,"%d",nAntiProtons); text->DrawText(-0.55,-0.72,tcount);
464  snprintf(tcount,12,"%d",nPionPlus); text->DrawText(-0.55,-0.77,tcount);
465  snprintf(tcount,12,"%d",nPionMinus); text->DrawText(-0.55,-0.82,tcount);
466  snprintf(tcount,12,"%d",nKaons); text->DrawText(-0.55,-0.87,tcount);
467  snprintf(tcount,12,"%d",nElectrons); text->DrawText(-0.55,-0.92,tcount);
468 
469  text->SetTextAlign(12);
470  if (nPionPlus+nPionMinus) {
471  snprintf(tcount,31,"Protons/Pions= %4f",Float_t(nProtons)/Float_t(nPionPlus+nPionMinus));
472  } else {
473  strlcpy(tcount,"Protons/Pions= inf",31);
474  }
475  text->DrawText(-0.45,-0.92,tcount);
476 
477  if (nPionPlus+nPionMinus) {
478  snprintf(tcount,12,"Kaons/Pions= %4f",Float_t(nKaons)/Float_t(nPionPlus+nPionMinus));
479  } else {
480  strlcpy(tcount,"Kaons/Pions= inf",31);
481  }
482  text->DrawText(0.30,-0.92,tcount);
483 }
484 
485 
486 ////////////////////////////////////////////////////////////////////////////////
487 ///*-*-*-*-*-*-*-*-*-*-*Execute action corresponding to one event*-*-*-*
488 ///*-* =========================================
489 
491 {
492  if (gPad->GetView()) {
493  gPad->GetView()->ExecuteRotateView(event, px, py);
494  return;
495  }
496 }
497 
498 ////////////////////////////////////////////////////////////////////////////////
499 /// Return the number of particles in the stack
500 
502 {
503  return fParticles->GetLast()+1;
504 }
505 
506 ////////////////////////////////////////////////////////////////////////////////
507 /// Returns pointer to primary number i;
508 ///
509 
511 {
512  if (!fParticles) return 0;
513  Int_t n = fParticles->GetLast();
514  if (i < 0 || i > n) return 0;
515  return (TParticle*)fParticles->UncheckedAt(i);
516 }
517 
518 ////////////////////////////////////////////////////////////////////////////////
519 ///
520 /// Paint one event
521 ///
522 
524 {
525 }
526 
527 ////////////////////////////////////////////////////////////////////////////////
528 ///
529 /// Set Pt threshold below which primaries are not drawn
530 ///
531 
533 {
534  fPtCut = ptcut;
535  Draw();
536  gPad->Update();
537 }
538 
539 ////////////////////////////////////////////////////////////////////////////////
540 ///
541 /// Set lower and upper values of the view range
542 ///
543 
545 {
546  SetViewRange(-rbox,-rbox,-rbox,rbox,rbox,rbox);
547 }
548 
549 ////////////////////////////////////////////////////////////////////////////////
550 ///
551 /// Set lower and upper values of the view range
552 ///
553 
555 {
556  TView *view = gPad->GetView();
557  if (!view) return;
558  view->SetRange(xmin,ymin,zmin,xmax,ymax,zmax);
559 
560  Draw();
561  gPad->Update();
562 }
563 
564 ////////////////////////////////////////////////////////////////////////////////
565 ///
566 /// Set flag to display or not neutrons
567 ///
568 
570 {
571  fShowNeutrons = show;
572  Draw();
573  gPad->Update();
574 }
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:52
virtual void SetLineWidth(Width_t lwidth)
Definition: TAttLine.h:57
Double_t Pt() const
Definition: TParticle.h:122
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
An array of TObjects.
Definition: TObjArray.h:39
float xmin
Definition: THbookFile.cxx:93
virtual void Draw(Option_t *option="")
Draw this pavetext with its current attributes.
Definition: TPaveText.cxx:211
Float_t fPtCut
Definition: TGenerator.h:149
Int_t GetNumberOfParticles() const
Return the number of particles in the stack.
Definition: TGenerator.cxx:501
ClassImp(TGenerator) TGenerator
Event generator default constructor.
Definition: TGenerator.cxx:148
Double_t Vx() const
Definition: TParticle.h:112
virtual void Clear(Option_t *option="")
Remove all objects from the array.
Definition: TObjArray.cxx:297
float Float_t
Definition: RtypesCore.h:53
Double_t Vy() const
Definition: TParticle.h:113
Int_t GetLast() const
Return index of last object in array.
Definition: TObjArray.cxx:527
Bool_t fShowNeutrons
Pt cut. Do not show primaries below.
Definition: TGenerator.h:150
const char Option_t
Definition: RtypesCore.h:62
virtual void Delete(Option_t *option="")
Remove all objects from the array AND delete all heap based objects.
Definition: TObjArray.cxx:328
float ymin
Definition: THbookFile.cxx:93
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
*-*-*-*-*-*-*-*-*-*-*Execute action corresponding to one event*-*-*-* *-* ===========================...
Definition: TGenerator.cxx:490
See TView3D.
Definition: TView.h:36
virtual TText * AddText(Double_t x1, Double_t y1, const char *label)
Add a new Text line to this pavetext at given coordinates.
Definition: TPaveText.cxx:160
virtual Int_t ImportParticles(TClonesArray *particles, Option_t *option="")
It reads the /HEPEVT/ common block which has been filled by the GenerateEvent method.
Definition: TGenerator.cxx:264
#define gROOT
Definition: TROOT.h:340
virtual void SetRange(const Double_t *min, const Double_t *max)=0
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition: TObject.cxx:254
Int_t GetEntriesFast() const
Definition: TObjArray.h:66
virtual void SetViewRange(Float_t xmin=-10000, Float_t ymin=-10000, Float_t zmin=-10000, Float_t xmax=10000, Float_t ymax=10000, Float_t zmax=10000)
Set lower and upper values of the view range.
Definition: TGenerator.cxx:554
Int_t jmohep[4000][2]
Definition: Hepevt.h:29
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:164
HEPEVT_DEF HEPEVT
Definition: Hepevt.h:39
static TDatabasePDG * Instance()
static function
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
Double_t phep[4000][5]
Definition: Hepevt.h:31
Double_t vhep[4000][4]
Definition: Hepevt.h:32
Double_t Log10(Double_t x)
Definition: TMath.h:529
Base class for several text objects.
Definition: TText.h:42
Int_t isthep[4000]
Definition: Hepevt.h:27
virtual void ShowNeutrons(Bool_t show=1)
Set flag to display or not neutrons.
Definition: TGenerator.cxx:569
virtual void GenerateEvent()
must be implemented in concrete class (see eg TPythia6)
Definition: TGenerator.cxx:183
virtual void SetTextAlign(Short_t align=11)
Definition: TAttText.h:55
virtual void Clear(Option_t *option="")
Clear the clones array.
TObject * UncheckedAt(Int_t i) const
Definition: TObjArray.h:91
virtual void SetLineColor(Color_t lcolor)
Definition: TAttLine.h:54
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:41
float ymax
Definition: THbookFile.cxx:93
TPaveText * pt
virtual void Draw(Option_t *option="")
Insert one event in the pad list.
Definition: TGenerator.cxx:343
virtual void SetFillColor(Color_t fcolor)
Definition: TAttFill.h:50
virtual void Paint(Option_t *option="")
Paint one event.
Definition: TGenerator.cxx:523
TParticlePDG * GetPDG(Int_t mode=0) const
returns a pointer to the TParticlePDG object using the pdgcode if mode == 0 (default) always get a fr...
Definition: TParticle.cxx:266
virtual TParticle * GetParticle(Int_t i) const
Returns pointer to primary number i;.
Definition: TGenerator.cxx:510
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
float xmax
Definition: THbookFile.cxx:93
virtual void SetViewRadius(Float_t rbox=1000)
Set lower and upper values of the view range.
Definition: TGenerator.cxx:544
Int_t nhep
Definition: Hepevt.h:26
Int_t jdahep[4000][2]
Definition: Hepevt.h:30
A Pave (see TPave) with text, lines or/and boxes inside.
Definition: TPaveText.h:35
Int_t idhep[4000]
Definition: Hepevt.h:28
double Double_t
Definition: RtypesCore.h:55
TObjArray * fParticles
display neutrons if true
Definition: TGenerator.h:151
TText * text
virtual void SetPtCut(Float_t ptcut=0)
Set Pt threshold below which primaries are not drawn.
Definition: TGenerator.cxx:532
Double_t Vz() const
Definition: TParticle.h:114
virtual void SetLineStyle(Style_t lstyle)
Definition: TAttLine.h:56
static TView * CreateView(Int_t system=1, const Double_t *rmin=0, const Double_t *rmax=0)
Create a concrete default 3-d view via the plug-in manager.
Definition: TView.cxx:36
#define name(a, b)
Definition: linkTestLib0.cpp:5
Double_t Energy() const
Definition: TParticle.h:123
An array of clone (identical) objects.
Definition: TClonesArray.h:32
virtual void Browse(TBrowser *b)
browse generator
Definition: TGenerator.cxx:318
virtual void GetRange(Float_t *min, Float_t *max)=0
virtual TText * DrawText(Double_t x, Double_t y, const char *text)
Draw this text with new coordinates.
Definition: TText.cxx:171
#define gPad
Definition: TVirtualPad.h:288
virtual void SetTextColor(Color_t tcolor=1)
Definition: TAttText.h:57
void Add(TObject *obj)
Definition: TObjArray.h:75
virtual ~TGenerator()
Event generator default destructor.
Definition: TGenerator.cxx:170
virtual void SetTextSize(Float_t tsize=1)
Definition: TAttText.h:60
const Bool_t kTRUE
Definition: Rtypes.h:91
const Int_t n
Definition: legend1.C:16
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
*-*-*-*-*-*-*-*Compute distance from point px,py to objects in event*-*-*-* *-* =====================...
Definition: TGenerator.cxx:329