ROOT logo
// pythia8 basic example
//Author: Andreas Morsch
//
// to run, do
//  root > .x pythia8
//
// Note that before executing this script, 
//   -the env variable PYTHIA8 must point to the pythia8100 (or newer) directoryDATA
//   -the env variable PYTHIA8DATA must be defined and it must point to $PYTHIA8/xmldoc
//
void pythia8(Int_t nev  = 100, Int_t ndeb = 1) {
    
   char* path = gSystem->ExpandPathName("$PYTHIA8DATA");
   if (gSystem->AccessPathName(path)) {
      Warning("pythia8.C", 
              "Environment variable PYTHIA8DATA must contain path to pythi8100/xmldoc directory !");
      return;
   }
    
// Load libraries
   gSystem->Load("$PYTHIA8/lib/libpythia8");
   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();
    
// Configure    
   pythia8->ReadString("SoftQCD:minBias = on");
   pythia8->ReadString("SoftQCD:singleDiffractive = on");
   pythia8->ReadString("SoftQCD:doubleDiffractive = on");

// 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();
         Int_t pdg = part->GetPdgCode();
         if (ist != 1) continue;
         Float_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
         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}");
 }
 pythia8.C:1
 pythia8.C:2
 pythia8.C:3
 pythia8.C:4
 pythia8.C:5
 pythia8.C:6
 pythia8.C:7
 pythia8.C:8
 pythia8.C:9
 pythia8.C:10
 pythia8.C:11
 pythia8.C:12
 pythia8.C:13
 pythia8.C:14
 pythia8.C:15
 pythia8.C:16
 pythia8.C:17
 pythia8.C:18
 pythia8.C:19
 pythia8.C:20
 pythia8.C:21
 pythia8.C:22
 pythia8.C:23
 pythia8.C:24
 pythia8.C:25
 pythia8.C:26
 pythia8.C:27
 pythia8.C:28
 pythia8.C:29
 pythia8.C:30
 pythia8.C:31
 pythia8.C:32
 pythia8.C:33
 pythia8.C:34
 pythia8.C:35
 pythia8.C:36
 pythia8.C:37
 pythia8.C:38
 pythia8.C:39
 pythia8.C:40
 pythia8.C:41
 pythia8.C:42
 pythia8.C:43
 pythia8.C:44
 pythia8.C:45
 pythia8.C:46
 pythia8.C:47
 pythia8.C:48
 pythia8.C:49
 pythia8.C:50
 pythia8.C:51
 pythia8.C:52
 pythia8.C:53
 pythia8.C:54
 pythia8.C:55
 pythia8.C:56
 pythia8.C:57
 pythia8.C:58
 pythia8.C:59
 pythia8.C:60
 pythia8.C:61
 pythia8.C:62
 pythia8.C:63
 pythia8.C:64
 pythia8.C:65
 pythia8.C:66
 pythia8.C:67
 pythia8.C:68
 pythia8.C:69
 pythia8.C:70
 pythia8.C:71
 pythia8.C:72
 pythia8.C:73
 pythia8.C:74
 pythia8.C:75
 pythia8.C:76
 pythia8.C:77
 pythia8.C:78
 pythia8.C:79
 pythia8.C:80
 pythia8.C:81
 pythia8.C:82
thumb