Logo ROOT   6.18/05
Reference Guide
pythia8.C File Reference

Detailed Description

pythia8 basic example

to run, do:

root > .x pythia8.C
#include "TSystem.h"
#include "TH1F.h"
#include "TClonesArray.h"
#include "TPythia8.h"
#include "TParticle.h"
#include "TDatabasePDG.h"
#include "TCanvas.h"
void pythia8(Int_t nev = 100, Int_t ndeb = 1)
{
// Load libraries
gSystem->Load("libEG");
gSystem->Load("libEGPythia8");
// Histograms
TH1F* etaH = new TH1F("etaH", "Pseudorapidity", 120, -12., 12.);
TH1F* ptH = new TH1F("ptH", "pt", 50, 0., 10.);
// Array of particles
TClonesArray* particles = new TClonesArray("TParticle", 1000);
// Create pythia8 object
TPythia8* pythia8 = new TPythia8();
#if PYTHIA_VERSION_INTEGER == 8235
// Pythia 8.235 is known to cause crashes:
printf("ABORTING PYTHIA8 TUTORIAL!\n");
printf("The version of Pythia you use is known to case crashes due to memory errors.\n");
printf("They have been reported to the authors; the Pythia versions 8.1... are known to work.\n");
return;
#endif
// Configure
pythia8->ReadString("HardQCD:all = on");
pythia8->ReadString("Random:setSeed = on");
// use a reproducible seed: always the same results for the tutorial.
pythia8->ReadString("Random:seed = 42");
// Initialize
pythia8->Initialize(2212 /* p */, 2212 /* p */, 14000. /* TeV */);
// Event loop
for (Int_t iev = 0; iev < nev; iev++) {
pythia8->GenerateEvent();
if (iev < ndeb) pythia8->EventListing();
pythia8->ImportParticles(particles,"All");
Int_t np = particles->GetEntriesFast();
// Particle loop
for (Int_t ip = 0; ip < np; ip++) {
TParticle* part = (TParticle*) particles->At(ip);
Int_t ist = part->GetStatusCode();
// Positive codes are final particles.
if (ist <= 0) continue;
Int_t pdg = part->GetPdgCode();
if (charge == 0.) continue;
Float_t eta = part->Eta();
Float_t pt = part->Pt();
etaH->Fill(eta);
if (pt > 0.) ptH->Fill(pt, 1./(2. * pt));
}
}
pythia8->PrintStatistics();
TCanvas* c1 = new TCanvas("c1","Pythia8 test example",800,800);
c1->Divide(1, 2);
c1->cd(1);
etaH->Scale(5./Float_t(nev));
etaH->Draw();
etaH->SetXTitle("#eta");
etaH->SetYTitle("dN/d#eta");
c1->cd(2);
gPad->SetLogy();
ptH->Scale(5./Float_t(nev));
ptH->Draw();
ptH->SetXTitle("p_{t} [GeV/c]");
ptH->SetYTitle("dN/dp_{t}^{2} [GeV/c]^{-2}");
}
int Int_t
Definition: RtypesCore.h:41
float Float_t
Definition: RtypesCore.h:53
R__EXTERN TSystem * gSystem
Definition: TSystem.h:560
#define gPad
Definition: TVirtualPad.h:286
The Canvas class.
Definition: TCanvas.h:31
An array of clone (identical) objects.
Definition: TClonesArray.h:32
static TDatabasePDG * Instance()
static function
TParticlePDG * GetParticle(Int_t pdgCode) const
Get a pointer to the particle object according to the MC code number.
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:571
virtual void SetXTitle(const char *title)
Definition: TH1.h:409
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3258
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2981
virtual void SetYTitle(const char *title)
Definition: TH1.h:410
virtual void Scale(Double_t c1=1, Option_t *option="")
Multiply this histogram by a constant c1.
Definition: TH1.cxx:6218
Int_t GetEntriesFast() const
Definition: TObjArray.h:64
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
Double_t Charge() const
Definition: TParticlePDG.h:68
Description of the dynamic properties of a particle.
Definition: TParticle.h:26
Int_t GetStatusCode() const
Definition: TParticle.h:82
Int_t GetPdgCode() const
Definition: TParticle.h:83
Double_t Pt() const
Definition: TParticle.h:135
Double_t Eta() const
Definition: TParticle.h:137
TPythia8 is an interface class to C++ version of Pythia 8.1 event generators, written by T....
Definition: TPythia8.h:77
void ReadString(const char *string) const
Configuration.
Definition: TPythia8.cxx:299
Bool_t Initialize(Int_t idAin, Int_t idBin, Double_t ecms)
Initialization.
Definition: TPythia8.cxx:146
void EventListing() const
Event listing.
Definition: TPythia8.cxx:363
void PrintStatistics() const
Print end of run statistics.
Definition: TPythia8.cxx:355
virtual void GenerateEvent()
Generate the next event.
Definition: TPythia8.cxx:184
virtual Int_t ImportParticles(TClonesArray *particles, Option_t *option="")
Import particles from Pythia stack.
Definition: TPythia8.cxx:193
virtual int Load(const char *module, const char *entry="", Bool_t system=kFALSE)
Load a shared library.
Definition: TSystem.cxx:1843
TPaveText * pt
return c1
Definition: legend1.C:41
Author
Andreas Morsch

Definition in file pythia8.C.