PYTHIA8/pythia8175/examples/main24.cc

   1 // main24.cc is a part of the PYTHIA event generator.
   2 // Copyright (C) 2013 Peter Skands, 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 to run SUSY processes in Pythia8.
   8 // All input is specified in the main22.cmnd file.
   9
  10 #include "Pythia.h"
  11
  12 using namespace Pythia8;
  13
  14 int main() {
  15
  16   // Generator. Shorthand for the event.
  17   Pythia pythia;
  18   Event& event = pythia.event;
  19
  20   // Read in commands from external file.
  21   pythia.readFile("main24.cmnd");
  22
  23   // Extract settings to be used in the main program.
  24   int nEvent   = pythia.mode("Main:numberOfEvents");
  25   int nAbort   = pythia.mode("Main:timesAllowErrors");
  26   double eCM   = pythia.parm("Beams:eCM");
  27
  28   // Initialize.
  29   pythia.init();
  30
  31   // Histograms.
  32   double epTol = 1e-6 * eCM;
  33   Hist epCons("deviation from energy-momentum conservation",100,0.,epTol);
  34   Hist nFinal("final particle multiplicity",100,-0.5,799.5);
  35   Hist dnparticledy("dn/dy for particles",100,-10.,10.);
  36
  37   // Begin event loop.
  38   int iAbort = 0;
  39   for (int iEvent = 0; iEvent < nEvent; ++iEvent) {
  40
  41     // Generate events. Quit if failure.
  42     if (!pythia.next()) {
  43       if (++iAbort < nAbort) continue;
  44       cout << " Event generation aborted prematurely, owing to error!\n";
  45       break;
  46     }
  47
  48     // Loop over final particles in the event.
  49     int nFin = 0;
  50     Vec4 pSum;
  51     for (int i = 0; i < event.size(); ++i) if (event[i].isFinal()) {
  52       nFin++;
  53       pSum += event[i].p();
  54       dnparticledy.fill(event[i].y());
  55     }
  56
  57     // Check and print event with too big energy-momentum deviation.
  58     nFinal.fill(nFin);
  59     double epDev = abs(pSum.e() - eCM) + abs(pSum.px()) + abs(pSum.py())
  60       + abs(pSum.pz());
  61     epCons.fill(epDev);
  62     if (epDev > epTol) {
  63       cout << " Warning! Event with epDev = " << scientific
  64            << setprecision(4) << epDev << " now listed:";
  65       event.list();
  66     }
  67
  68   // End of event loop.
  69   }
  70
  71   // Final statistics and histogram output.
  72   pythia.stat();
  73   cout << epCons << nFinal << dnparticledy;
  74
  75   return 0;
  76 }
  77