Update to pythi8.170
[u/mrichter/AliRoot.git] / PYTHIA8 / pythia8170 / examples / main03.cc
1 // main03.cc is a part of the PYTHIA event generator.
2 // Copyright (C) 2012 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5
6 // This is a simple test program. 
7 // It illustrates how different processes can be selected and studied.
8 // All input is specified in the main03.cmnd file.
9
10 #include "Pythia.h"
11
12 using namespace Pythia8; 
13
14 int main() {
15
16   // Generator.
17   Pythia pythia;
18
19   // Shorthand for the event record in pythia.
20   Event& event = pythia.event;
21
22   // Read in commands from external file.
23   pythia.readFile("main03.cmnd");    
24
25   // Extract settings to be used in the main program.
26   int nEvent = pythia.mode("Main:numberOfEvents");
27   int nAbort = pythia.mode("Main:timesAllowErrors");
28
29   // Initialize. 
30   pythia.init();
31
32   // Book histograms.
33   Hist pThard("process pT scale", 100, 0., 500.);
34   Hist nFinal("final particle multiplicity", 100, -0.5, 1599.5);
35   Hist nCharged("charged particle multiplicity", 100, -0.5, 799.5);
36   Hist dndy("dn/dy for charged particles", 100, -10., 10.);
37   Hist dndeta("dn/d(eta) for charged particles", 100, -10., 10.);
38   Hist dndpT("dn/dpT for charged particles", 100, 0., 10.);
39   Hist epCons("deviation from energy-momentum conservation", 100, 0., 1e-6);
40
41   // Begin event loop.
42   int iAbort = 0;
43   for (int iEvent = 0; iEvent < nEvent; ++iEvent) {
44
45     // Generate events. Quit if many failures.
46     if (!pythia.next()) {
47       if (++iAbort < nAbort) continue;
48       cout << " Event generation aborted prematurely, owing to error!\n"; 
49       break;
50     }
51
52     // Fill hard scale of event.
53     pThard.fill( pythia.info.pTHat() );
54
55     // Loop over final particles in the event. 
56     int  nFin = 0;
57     int  nChg = 0;
58     Vec4 pSum;
59     for (int i = 0; i < event.size(); ++i) if (event[i].isFinal()) {
60
61        // Analyze all particles.
62       nFin++;
63       pSum += event[i].p();
64        
65       // Analyze charged particles and fill histograms.
66       if (event[i].isCharged()) {
67         ++nChg;
68         dndy.fill( event[i].y() );
69         dndeta.fill( event[i].eta() );
70         dndpT.fill( event[i].pT() );
71       }
72
73     // End of particle loop. Fill global properties.
74     }
75     nFinal.fill( nFin );
76     nCharged.fill( nChg );
77     pSum /= event[0].e();
78     double epDev = abs(pSum.e() - 1.) + abs(pSum.px()) + abs(pSum.py())
79       + abs(pSum.pz());
80     epCons.fill(epDev);
81
82   // End of event loop.
83   }
84    
85   // Final statistics. Normalize and output histograms.
86   pythia.stat();
87   dndy   *=  5. / nEvent;
88   dndeta *=  5. / nEvent;
89   dndpT  *= 10. / nEvent;
90   cout << pThard << nFinal << nCharged << dndy << dndeta << dndpT << epCons; 
91
92   // Done.
93   return 0;
94 }