StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
FragmentationSystems.h
1 // FragmentationSystems.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 // This file contains auxiliary classes in the fragmentation process.
7 // ColSinglet contains info on an individual singlet.
8 // ColConfig describes the colour configuration of the whole event.
9 // StringRegion keeps track on string momenta and directions.
10 // StringSystem contains all the StringRegions of the colour singlet.
11 
12 #ifndef Pythia8_FragmentationSystems_H
13 #define Pythia8_FragmentationSystems_H
14 
15 #include "Basics.h"
16 #include "Event.h"
17 #include "FragmentationFlavZpT.h"
18 #include "Info.h"
19 #include "ParticleData.h"
20 #include "PythiaStdlib.h"
21 #include "Settings.h"
22 
23 namespace Pythia8 {
24 
25 //==========================================================================
26 
27 // The ColSinglet class contains info on an individual singlet.
28 // Only to be used inside ColConfig, so no private members.
29 
30 class ColSinglet {
31 
32 public:
33 
34  // Constructors.
35  ColSinglet() : pSum(0., 0., 0., 0.), mass(0.), massExcess(0.),
36  hasJunction(false), isClosed(false), isCollected(false) {}
37  ColSinglet(vector<int>& iPartonIn, Vec4 pSumIn, double massIn,
38  double massExcessIn, bool hasJunctionIn = false,
39  bool isClosedIn = false, bool isCollectedIn = false)
40  : iParton(iPartonIn), pSum(pSumIn), mass(massIn),
41  massExcess(massExcessIn), hasJunction(hasJunctionIn),
42  isClosed(isClosedIn), isCollected(isCollectedIn) {}
43 
44  // Size of iParton array.
45  int size() const { return iParton.size();}
46 
47  // Stored quantities.
48  vector<int> iParton;
49  Vec4 pSum;
50  double mass, massExcess;
51  bool hasJunction, isClosed, isCollected;
52 
53 };
54 
55 //==========================================================================
56 
57 // The ColConfig class describes the colour configuration of the whole event.
58 
59 class ColConfig {
60 
61 public:
62 
63  // Constructor.
64  ColConfig() {singlets.resize(0);}
65 
66  // Initialize and save pointers.
67  void init(Info* infoPtrIn, Settings& settings, StringFlav* flavSelPtrIn);
68 
69  // Number of colour singlets.
70  int size() const {return singlets.size();}
71 
72  // Overload index operator to access separate colour singlets.
73  ColSinglet& operator[](int iSub) {return singlets[iSub];}
74 
75  // Clear contents.
76  void clear() {singlets.resize(0);}
77 
78  // Insert a new colour singlet system in ascending mass order.
79  // Calculate its properties. Join nearby partons.
80  bool insert( vector<int>& iPartonIn, Event& event);
81 
82  // Erase a colour singlet system. (Rare operation.)
83  void erase(int iSub) {singlets.erase(singlets.begin() + iSub);}
84 
85  // Collect all partons of singlet to be consecutively ordered.
86  void collect(int iSub, Event& event, bool skipTrivial = true);
87 
88  // Find to which singlet system a particle belongs.
89  int findSinglet(int i);
90 
91  // List all currently identified singlets.
92  void list(ostream& os = cout) const;
93 
94 private:
95 
96  // Constants: could only be changed in the code itself.
97  static const double CONSTITUENTMASS;
98 
99  // Pointer to various information on the generation.
100  Info* infoPtr;
101 
102  // Pointer to class for flavour generation.
103  StringFlav* flavSelPtr;
104 
105  // Initialization data, to be read from Settings.
106  double mJoin, mJoinJunction, mStringMin;
107 
108  // List of all separate colour singlets.
109  vector<ColSinglet> singlets;
110 
111  // Join two legs of junction to a diquark for small invariant masses.
112  bool joinJunction( vector<int>& iPartonIn, Event& event,
113  double massExcessIn);
114 
115 };
116 
117 //==========================================================================
118 
119 // The StringRegion class contains the information related to
120 // one string section in the evolution of a multiparton system.
121 // Only to be used inside StringFragmentation and MiniStringFragmentation,
122 // so no private members.
123 
125 
126 public:
127 
128  // Constructor.
129  StringRegion() : isSetUp(false), isEmpty(true) {}
130 
131  // Constants: could only be changed in the code itself.
132  static const double MJOIN, TINY;
133 
134  // Data members.
135  bool isSetUp, isEmpty;
136  Vec4 pPos, pNeg, eX, eY;
137  double w2, xPosProj, xNegProj, pxProj, pyProj;
138 
139  // Set up four-vectors for longitudinal and transverse directions.
140  void setUp(Vec4 p1, Vec4 p2, bool isMassless = false);
141 
142  // Construct a four-momentum from (x+, x-, px, py).
143  Vec4 pHad( double xPosIn, double xNegIn, double pxIn, double pyIn)
144  { return xPosIn * pPos + xNegIn * pNeg + pxIn * eX + pyIn * eY; }
145 
146  // Project a four-momentum onto (x+, x-, px, py). Read out projection.
147  void project(Vec4 pIn);
148  void project( double pxIn, double pyIn, double pzIn, double eIn)
149  { project( Vec4( pxIn, pyIn, pzIn, eIn) ); }
150  double xPos() const {return xPosProj;}
151  double xNeg() const {return xNegProj;}
152  double px() const {return pxProj;}
153  double py() const {return pyProj;}
154 
155 };
156 
157 //==========================================================================
158 
159 // The StringSystem class contains the complete set of all string regions.
160 // Only to be used inside StringFragmentation, so no private members.
161 
163 
164 public:
165 
166  // Constructor.
167  StringSystem() {}
168 
169  // Set up system from parton list.
170  void setUp(vector<int>& iSys, Event& event);
171 
172  // Calculate string region from (iPos, iNeg) pair.
173  int iReg( int iPos, int iNeg) const
174  {return (iPos * (indxReg - iPos)) / 2 + iNeg;}
175 
176  // Reference to string region specified by (iPos, iNeg) pair.
177  StringRegion& region(int iPos, int iNeg) {return system[iReg(iPos, iNeg)];}
178 
179  // Reference to low string region specified either by iPos or iNeg.
180  StringRegion& regionLowPos(int iPos) {
181  return system[iReg(iPos, iMax - iPos)]; }
182  StringRegion& regionLowNeg(int iNeg) {
183  return system[iReg(iMax - iNeg, iNeg)]; }
184 
185  // Main content: a vector with all the string regions of the system.
186  vector<StringRegion> system;
187 
188  // Other data members.
189  int sizePartons, sizeStrings, sizeRegions, indxReg, iMax;
190  double mJoin, m2Join;
191 
192 };
193 
194 //==========================================================================
195 
196 } // end namespace Pythia8
197 
198 #endif // Pythia8_FragmentationSystems_H