1afbab1cc91403fb3da195efc1f8aa7e4b873de2
[u/mrichter/AliRoot.git] / test / vmctest / ppbench / gen.C
1 // $Id$
2 //
3 // Generation of kinematics tree which can be then used as an extrenal
4 // event generator.
5 // According to: $ALICE_ROOT/test/genkine/gen/fastGen.C
6 //
7 // By I. Hrivnacova, IPN Orsay
8
9 #if !defined(__CINT__) || defined(__MAKECINT__)
10 #include <Riostream.h>
11 #include <TH1F.h>
12 #include <TStopwatch.h>
13 #include <TDatime.h>
14 #include <TRandom.h>
15 #include <TDatabasePDG.h>
16 #include <TParticle.h>
17 #include <TArrayI.h>
18
19 #include "AliGenerator.h"
20 #include "AliPDG.h"
21 #include "AliRunLoader.h"
22 #include "AliRun.h"
23 #include "AliStack.h"
24 #include "AliHeader.h"
25 #include "PYTHIA6/AliGenPythia.h"
26 #include "PYTHIA6/AliPythia.h"
27 #endif
28
29 void gen(Int_t nev = 1, char* filename = "galice.root")
30 {
31   // Load libraries
32   // gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include -I$ALICE_ROOT");
33   gSystem->Load("liblhapdf.so");      // Parton density functions
34   gSystem->Load("libEGPythia6.so");   // TGenerator interface
35   gSystem->Load("libpythia6.so");     // Pythia
36   gSystem->Load("libAliPythia6.so");  // ALICE specific implementations
37
38   AliPDG::AddParticlesToPdgDataBase();
39   TDatabasePDG::Instance();
40
41   // Run loader
42   AliRunLoader* rl = AliRunLoader::Open("galice.root","FASTRUN","recreate");
43   
44   rl->SetCompressionLevel(2);
45   rl->SetNumberOfEventsPerFile(nev);
46   rl->LoadKinematics("RECREATE");
47   rl->MakeTree("E");
48   gAlice->SetRunLoader(rl);
49   
50   //  Create stack
51   rl->MakeStack();
52   AliStack* stack = rl->Stack();
53   
54   //  Header
55   AliHeader* header = rl->GetHeader();
56   
57   //  Create and Initialize Generator
58   gROOT->LoadMacro("$ALICE_ROOT/test/vmctest/ppbench/genPPbenchConfig.C");
59   AliGenerator* gener = genPPbenchConfig();
60
61   // Go to galice.root
62   rl->CdGAFile();
63
64   // Forbid some decays. Do it after gener->Init(0, because
65   // the initialization of the generator includes reading of the decay table.
66   // ...
67
68   //
69   // Event Loop
70   //
71   
72   TStopwatch timer;
73   timer.Start();
74   for (Int_t iev = 0; iev < nev; iev++) {
75     
76     cout <<"Event number "<< iev << endl;
77     
78     // Initialize event
79     header->Reset(0,iev);
80     rl->SetEventNumber(iev);
81     stack->Reset();
82     rl->MakeTree("K");
83     
84     // Generate event
85     stack->Reset();
86     stack->ConnectTree(rl->TreeK());
87     gener->Generate();
88     cout << "Number of particles " << stack->GetNprimary() << endl;
89     
90     // Finish event
91     header->SetNprimary(stack->GetNprimary());
92     header->SetNtrack(stack->GetNtrack());  
93     
94     // I/O
95     stack->FinishEvent();
96     header->SetStack(stack);
97     rl->TreeE()->Fill();
98     rl->WriteKinematics("OVERWRITE");
99     
100   } // event loop
101   timer.Stop();
102   timer.Print();
103   
104   //                         Termination
105   //  Generator
106   gener->FinishRun();
107   //  Write file
108   rl->WriteHeader("OVERWRITE");
109   gener->Write();
110   rl->Write();
111 }