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