StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
main27.cc
1 // main27.cc is a part of the PYTHIA event generator.
2 // Copyright (C) 2014 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 // Kaluza-Klein gamma*/Z resonances in TeV-sized extra dimensions.
7 
8 #include <assert.h>
9 #include <time.h>
10 #include <sstream>
11 
12 #include "Pythia8/Pythia.h"
13 using namespace Pythia8;
14 
15 int main() {
16 
17  // Generator.
18  Pythia pythia;
19  ParticleData& pdt = pythia.particleData;
20 
21  // Pick new random number seed for each run, based on clock.
22  pythia.readString("Random:setSeed = on");
23  pythia.readString("Random:seed = 0");
24 
25  // Process selection.
26  // ANY COMBINATION OF THE PROCESSES FLAGS BELOW IS ALLOWED
27  // HERE WE SWITCH ON ONLY THE MU+MU- FINAL STATE.
28  // TO SWITCH ALL POSSIBLE FINAL STATES ON, UNCOMMENT ALL
29  // THE RELEVANT LINES BELOW:
30  //pythia.readString("ExtraDimensionsTEV:ffbar2e+e- = on");
31  pythia.readString("ExtraDimensionsTEV:ffbar2mu+mu- = on");
32  //pythia.readString("ExtraDimensionsTEV:ffbar2tau+tau- = on");
33  //pythia.readString("ExtraDimensionsTEV:ffbar2uubar = on");
34  //pythia.readString("ExtraDimensionsTEV:ffbar2ddbar = on");
35  //pythia.readString("ExtraDimensionsTEV:ffbar2ccbar = on");
36  //pythia.readString("ExtraDimensionsTEV:ffbar2ssbar = on");
37  //pythia.readString("ExtraDimensionsTEV:ffbar2bbbar = on");
38  //pythia.readString("ExtraDimensionsTEV:ffbar2ttbar = on");
39  //pythia.readString("ExtraDimensionsTEV:ffbar2nuenuebar = on");
40  //pythia.readString("ExtraDimensionsTEV:ffbar2numunumubar = on");
41  //pythia.readString("ExtraDimensionsTEV:ffbar2nutaunutaubar = on");
42 
43  // Pick KK mass.
44  double newMass = 4000.; // GeV
45  cout << "|-------------------------" << endl;
46  cout << "| KK mass is: " << newMass << endl;
47  cout << "|-------------------------" << endl;
48  stringstream strm;
49  string sNewMass, sNewWidth, sNewLowBound, sNewHighBound;
50 
51  // Manually set the mass and therefore the width
52  // and the phase space for the sampling
53  strm.clear();
54  strm << newMass;
55  strm >> sNewMass;
56  strm.clear();
57  strm << newMass / pdt.m0(5000023) * pdt.mWidth(5000023);
58  strm >> sNewWidth;
59  strm.clear();
60  strm << newMass/4.;
61  strm >> sNewLowBound;
62  strm.clear();
63  strm << newMass*2.;
64  strm >> sNewHighBound;
65 
66  // Feed in KK state information and other generation specifics.
67  pythia.readString("5000023:m0 = " + sNewMass);
68  pythia.readString("5000023:mWidth = " + sNewWidth);
69  pythia.readString("5000023:mMin = " + sNewLowBound);
70  pythia.readString("5000023:mMax = " + sNewHighBound);
72  pythia.readString("5000023:isResonance = false"); // THIS IS MANDATORY //
74  // 0=(gm+Z), 1=(gm), 2=(Z), 3=(gm+Z+gmKK+ZKK), 4=(m+Z+gmKK), 5=(m+Z+ZKK)
75  pythia.readString("ExtraDimensionsTEV:gmZmode = 3");
76  // min=0, max=100, default=10
77  pythia.readString("ExtraDimensionsTEV:nMax = 100");
78  pythia.readString("ExtraDimensionsTEV:mStar = " + sNewMass);
79  pythia.readString("PhaseSpace:mHatMin = " + sNewLowBound);
80  pythia.readString("PhaseSpace:mHatMax = " + sNewHighBound);
81 
82  // Initialize for LHC.
83  pythia.readString("Beams:eCM = 14000.");
84  pythia.init();
85 
86  // Histograms.
87  Hist mHatHisto("dN/dmHat", 50, newMass/4., newMass*2.);
88  Hist pTmuHisto("(dN/dpT)_mu^-", 50, 1., 2501.);
89 
90  vector<int> moms;
91 
92  // Measure the cpu runtime.
93  clock_t start, stop;
94  double t = 0.0;
95  // Depending on operating system, either of lines below gives warning.
96  //assert((start = clock()) != -1); // Start timer; clock_t signed.
97  //assert((start = clock()) != -1u); // Start timer; clock_t unsigned.
98  // Simpler option, not using assert.
99  start = clock();
100 
101  // Begin event loop. Generate event. Skip if error. List first one.
102  for (int iEvent = 0 ; iEvent < 500 ; ++iEvent) {
103  if (!pythia.next()) continue;
104 
105  // Begin event analysis.
106  bool isZ = false;
107  bool ismu = false;
108  int iZ = 0;
109  int imu = 0;
110  for (int i = 0 ; i < pythia.event.size() ; ++i) {
111  // find the most recent Z
112  if (pythia.event[i].id() == 5000023) {
113  iZ = i;
114  isZ = true;
115  }
116  // find the final muon who's first mother is the Z
117  if (pythia.event[i].id() == 13 && pythia.event[i].isFinal()) {
118  moms.clear();
119  moms = pythia.event.motherList(i);
120  for (int m = 0 ; m < int(moms.size()) ; m++) {
121  if( pythia.event[ moms[m] ].id() == 5000023 ) {
122  imu = i;
123  ismu = true;
124  break;
125  } // end if 5000023
126  } // end for moms.size()
127  } // end if final muon
128  } // end for event.size()
129  if(isZ && ismu) {
130  mHatHisto.fill( pythia.event[iZ].m() );
131  pTmuHisto.fill( pythia.event[imu].pT() );
132  }
133  if(iEvent%10 == 0) cout << "Event: " << iEvent << endl << flush;
134  } // end for iEvent<500
135 
136  // Done. Print results.
137  stop = clock(); // Stop timer
138  t = (double) (stop-start)/CLOCKS_PER_SEC;
139 
140  pythia.stat();
141  cout << mHatHisto;
142  cout << pTmuHisto;
143 
144  cout << "\n" << "|----------------------------------------|" << endl;
145  cout << "| CPU Runtime = " << t << " sec" << endl;
146  cout << "|----------------------------------------|" << "\n" << endl;
147 
148  return 0;
149 }