source: HiSusy/trunk/Pythia8/pythia8170/include/HadronScatter.h @ 1

Last change on this file since 1 was 1, checked in by zerwas, 11 years ago

first import of structure, PYTHIA8 and DELPHES

File size: 5.5 KB
Line 
1// HadronScatter.h 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#ifndef Pythia8_HadronScatter_H
7#define Pythia8_HadronScatter_H
8
9#include "Event.h"
10#include "Info.h"
11#include "ParticleData.h"
12#include "PythiaComplex.h"
13
14#include <algorithm>
15#include <map>
16#include <set>
17using std::map;
18using std::set;
19using std::sort;
20
21namespace Pythia8 {
22
23class SigmaPartialWave {
24public:
25  // Initialisation
26  bool init(int, string, string, Info *, ParticleData *, Rndm *);
27
28  // Read data file
29  bool readFile(string, string);
30
31  // Set the subprocess/incoming particles
32  bool setSubprocess(int);
33  bool setSubprocess(int, int);
34
35  // Return sigma total/elastic, dSigma/dCos(theta)
36  double sigmaEl(double Wcm)  { return sigma(0, Wcm); }
37  double sigmaTot(double Wcm) { return sigma(1, Wcm); }
38  double dSigma(double Wcm, double cTheta) { return sigma(2, Wcm, cTheta); }
39
40  // Return a cos(theta) value
41  double pickCosTheta(double);
42
43  // Return maximum sigma elastic
44  double getSigmaElMax() { return sigElMax; }
45
46private:
47  // Pointers
48  Info         *infoPtr;
49  ParticleData *particleDataPtr;
50  Rndm         *rndmPtr;
51
52  // Constants
53  static const int     LSHIFT, ISHIFT, SUBBIN, ITER;
54  static const double  CONVERT2MB, WCMBIN, CTBIN, MASSSAFETY, GRIDSAFETY;
55  static const complex BINEND;
56
57  // Settings
58  int process, subprocess, subprocessMax, norm;
59
60  // Current incoming and maximum L/I values
61  int idA, idB, Lmax, Imax;
62
63  // Masses of incoming, last bin, maximum sigma elastic
64  double mA, mB, binMax, sigElMax;
65
66  // Map subprocess to incoming and vice versa:
67  //   sp2in[subprocess] = idA, idB
68  //   in2sp[idA, idB]   = subprocess
69  map < int, pair < int, int > > sp2in;
70  map < pair < int, int >, int > in2sp;
71
72  // Isospin coefficients isoCoeff[subprocess][2I]
73  map < int, map < int, double > > isoCoeff;
74
75  // Storage for partial wave data:
76  //   pwData[L * LSHIFT + 2I * ISHIFT + 2J][eNow] = T
77  map < int, map < double, complex > > pwData;
78
79  // Values of Pl and Pl' as computed by legendreP
80  vector < double > PlVec, PlpVec;
81
82  // Integration grid [subprocess][WcmBin][cosThetaBin]
83  vector < vector < vector < double > > > gridMax;
84  vector < vector < double > >            gridNorm;
85
86  // Setup subprocesses (including isospin coefficients)
87  void setupSubprocesses();
88
89  // Setup grids for integration
90  void setupGrid();
91
92  // Routine for calculating sigma total/elastic and dSigma/dCos(theta)
93  double sigma(int, double, double = 0.);
94
95  // Generate Legendre polynomials (and optionally derivatives)
96  void legendreP(double, bool = false);
97
98};
99
100
101//==========================================================================
102
103// HadronScatterPair class
104//   Simple class to hold details of a pair of hadrons which will scatter.
105//   Stores indices in event record and the measure used for ordering
106
107// Store a pair of indices
108typedef pair < int, int > HSIndex;
109
110class HadronScatterPair {
111public:
112  // Constructor
113  HadronScatterPair(const HSIndex &i1in, int yt1in, int pt1in,
114                    const HSIndex &i2in, int yt2in, int pt2in,
115                    double measureIn) :
116      i1(i1in), yt1(yt1in), pt1(pt1in),
117      i2(i2in), yt2(yt2in), pt2(pt2in),
118      measure(measureIn) {}
119
120  // Operator for sorting according to ordering measure
121  bool operator<(const HadronScatterPair& in) const {
122    return this->measure < in.measure;
123  }
124
125  // Indices into event record of hadrons to scatter
126  HSIndex i1;
127  int     yt1, pt1;
128  HSIndex i2;
129  int     yt2, pt2;
130  // Ordering measure
131  double measure;
132};
133
134
135//==========================================================================
136
137// HadronScatter class
138
139class HadronScatter {
140
141public:
142
143  // Constructor.
144  HadronScatter() {}
145
146  // Initialisation
147  bool init(Info* infoPtrIn, Settings& settings, Rndm* rndmPtrIn,
148            ParticleData *particleDataPtr);
149
150  // Perform all hadron scatterings
151  void scatter(Event&);
152 
153private: 
154
155  // Pointer to various information on the generation.
156  Info* infoPtr;
157  Rndm* rndmPtr;
158
159  // Settings
160  bool   doHadronScatter, afterDecay, allowDecayProd,
161         scatterRepeat, doTile;
162  int    hadronSelect, scatterProb;
163  double Npar, kPar, pPar, jPar, rMax, rMax2;
164  double pTsigma, pTsigma2, pT0MPI;
165
166  // Tiling
167  int    ytMax, ptMax;
168  double yMin, yMax, ytSize, ptSize;
169  vector < vector < set < HSIndex > > > tile;
170
171  // Keep track of scattered pairs
172  set < HSIndex > scattered;
173
174  // Partial wave amplitudes
175  SigmaPartialWave sigmaPW[3];
176
177  // Maximum sigma elastic
178  double sigElMax;
179
180  // Decide if a hadron can scatter
181  bool canScatter(Event &, int);
182
183  // Probability for a pair of hadrons to scatter
184  bool doesScatter(Event &, const HSIndex &, const HSIndex &);
185
186  // Calculate measure for ordering of scatterings
187  double measure(Event &, int, int);
188
189  // Perform a single hadron scattering
190  bool hadronScatter(Event &, HadronScatterPair &);
191
192  // Tiling functions
193  bool tileIntProb(vector < HadronScatterPair > &, Event &,
194                   const HSIndex &, int, int, bool);
195  int yTile(Event& event, int idx) {
196    return (doTile) ? int((event[idx].y() - yMin) / ytSize) : 0;
197  }
198  int pTile(Event& event, int idx) {
199    return (doTile) ? int((event[idx].phi() + M_PI) / ptSize) : 0;
200  }
201 
202  // Debug
203  void debugOutput();
204};
205 
206//==========================================================================
207
208
209} // end namespace Pythia8
210
211#endif // Pythia8_HadronScatter_H
212
Note: See TracBrowser for help on using the repository browser.