Line data Source code
1 : // HadronLevel.cc is a part of the PYTHIA event generator.
2 : // Copyright (C) 2015 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 : // Function definitions (not found in the header) for the HadronLevel class.
7 :
8 : #include "Pythia8/HadronLevel.h"
9 :
10 : namespace Pythia8 {
11 :
12 : //==========================================================================
13 :
14 : // The HadronLevel class.
15 :
16 : //--------------------------------------------------------------------------
17 :
18 : // Find settings. Initialize HadronLevel classes as required.
19 :
20 : bool HadronLevel::init(Info* infoPtrIn, Settings& settings,
21 : ParticleData* particleDataPtrIn, Rndm* rndmPtrIn,
22 : Couplings* couplingsPtrIn, TimeShower* timesDecPtr,
23 : RHadrons* rHadronsPtrIn, DecayHandler* decayHandlePtr,
24 : vector<int> handledParticles) {
25 :
26 : // Save pointers.
27 0 : infoPtr = infoPtrIn;
28 0 : particleDataPtr = particleDataPtrIn;
29 0 : rndmPtr = rndmPtrIn;
30 0 : couplingsPtr = couplingsPtrIn;
31 0 : rHadronsPtr = rHadronsPtrIn;
32 :
33 : // Main flags.
34 0 : doHadronize = settings.flag("HadronLevel:Hadronize");
35 0 : doDecay = settings.flag("HadronLevel:Decay");
36 0 : doBoseEinstein = settings.flag("HadronLevel:BoseEinstein");
37 :
38 : // Boundary mass between string and ministring handling.
39 0 : mStringMin = settings.parm("HadronLevel:mStringMin");
40 :
41 : // For junction processing.
42 0 : eNormJunction = settings.parm("StringFragmentation:eNormJunction");
43 :
44 : // Allow R-hadron formation.
45 0 : allowRH = settings.flag("RHadrons:allow");
46 :
47 : // Particles that should decay or not before Bose-Einstein stage.
48 0 : widthSepBE = settings.parm("BoseEinstein:widthSep");
49 :
50 : // Hadron scattering --rjc
51 0 : doHadronScatter = settings.flag("HadronScatter:scatter");
52 0 : hsAfterDecay = settings.flag("HadronScatter:afterDecay");
53 :
54 : // Initialize auxiliary fragmentation classes.
55 0 : flavSel.init(settings, rndmPtr);
56 0 : pTSel.init(settings, *particleDataPtr, rndmPtr);
57 0 : zSel.init(settings, *particleDataPtr, rndmPtr);
58 :
59 : // Initialize auxiliary administrative class.
60 0 : colConfig.init(infoPtr, settings, &flavSel);
61 :
62 : // Initialize string and ministring fragmentation.
63 0 : stringFrag.init(infoPtr, settings, particleDataPtr, rndmPtr,
64 : &flavSel, &pTSel, &zSel);
65 0 : ministringFrag.init(infoPtr, settings, particleDataPtr, rndmPtr,
66 : &flavSel, &pTSel, &zSel);
67 :
68 : // Initialize particle decays.
69 0 : decays.init(infoPtr, settings, particleDataPtr, rndmPtr, couplingsPtr,
70 0 : timesDecPtr, &flavSel, decayHandlePtr, handledParticles);
71 :
72 : // Initialize BoseEinstein.
73 0 : boseEinstein.init(infoPtr, settings, *particleDataPtr);
74 :
75 : // Initialize HadronScatter --rjc
76 0 : if (doHadronScatter)
77 0 : hadronScatter.init(infoPtr, settings, rndmPtr, particleDataPtr);
78 :
79 : // Initialize Hidden-Valley fragmentation, if necessary.
80 0 : useHiddenValley = hiddenvalleyFrag.init(infoPtr, settings,
81 0 : particleDataPtr, rndmPtr);
82 :
83 : // Send flavour and z selection pointers to R-hadron machinery.
84 0 : rHadronsPtr->fragPtrs( &flavSel, &zSel);
85 :
86 : // Initialize the colour tracing class.
87 0 : colTrace.init(infoPtr);
88 :
89 : // Initialize the junction splitting class.
90 0 : junctionSplitting.init(infoPtr, settings, rndmPtr, particleDataPtr);
91 :
92 : // Done.
93 0 : return true;
94 :
95 0 : }
96 :
97 : //--------------------------------------------------------------------------
98 :
99 : // Hadronize and decay the next parton-level.
100 :
101 : bool HadronLevel::next( Event& event) {
102 :
103 : // Store current event size to mark Parton Level content.
104 0 : event.savePartonLevelSize();
105 :
106 : // Do Hidden-Valley fragmentation, if necessary.
107 0 : if (useHiddenValley) hiddenvalleyFrag.fragment(event);
108 :
109 : // Colour-octet onia states must be decayed to singlet + gluon.
110 0 : if (!decayOctetOnia(event)) return false;
111 :
112 : // remove junction structures.
113 0 : if (!junctionSplitting.checkColours(event)) {
114 0 : infoPtr->errorMsg("Error in HadronLevel::next: "
115 : "failed colour/junction check");
116 0 : return false;
117 : }
118 :
119 : // Possibility of hadronization inside decay, but then no BE second time.
120 : // Hadron scattering, first pass only --rjc
121 : bool moreToDo, firstPass = true;
122 0 : bool doBoseEinsteinNow = doBoseEinstein;
123 0 : do {
124 : moreToDo = false;
125 :
126 : // First part: string fragmentation.
127 0 : if (doHadronize) {
128 :
129 : // Find the complete colour singlet configuration of the event.
130 0 : if (!findSinglets( event)) return false;
131 :
132 : // Fragment off R-hadrons, if necessary.
133 0 : if (allowRH && !rHadronsPtr->produce( colConfig, event))
134 0 : return false;
135 :
136 : // Process all colour singlet (sub)systems.
137 0 : for (int iSub = 0; iSub < colConfig.size(); ++iSub) {
138 :
139 : // Collect sequentially all partons in a colour singlet subsystem.
140 0 : colConfig.collect(iSub, event);
141 :
142 : // String fragmentation of each colour singlet (sub)system.
143 0 : if ( colConfig[iSub].massExcess > mStringMin ) {
144 0 : if (!stringFrag.fragment( iSub, colConfig, event)) return false;
145 :
146 : // Low-mass string treated separately. Tell if diffractive system.
147 : } else {
148 0 : bool isDiff = infoPtr->isDiffractiveA() || infoPtr->isDiffractiveB();
149 0 : if (!ministringFrag.fragment( iSub, colConfig, event, isDiff))
150 0 : return false;
151 0 : }
152 : }
153 : }
154 :
155 : // Hadron scattering --rjc
156 0 : if (doHadronScatter && !hsAfterDecay && firstPass)
157 0 : hadronScatter.scatter(event);
158 :
159 : // Second part: sequential decays of short-lived particles (incl. K0).
160 0 : if (doDecay) {
161 :
162 : // Loop through all entries to find those that should decay.
163 : int iDec = 0;
164 0 : do {
165 0 : Particle& decayer = event[iDec];
166 0 : if ( decayer.isFinal() && decayer.canDecay() && decayer.mayDecay()
167 0 : && (decayer.mWidth() > widthSepBE || decayer.idAbs() == 311) ) {
168 0 : decays.decay( iDec, event);
169 0 : if (decays.moreToDo()) moreToDo = true;
170 : }
171 0 : } while (++iDec < event.size());
172 0 : }
173 :
174 : // Hadron scattering --rjc
175 0 : if (doHadronScatter && hsAfterDecay && firstPass)
176 0 : hadronScatter.scatter(event);
177 :
178 : // Third part: include Bose-Einstein effects among current particles.
179 0 : if (doBoseEinsteinNow) {
180 0 : if (!boseEinstein.shiftEvent(event)) return false;
181 : doBoseEinsteinNow = false;
182 0 : }
183 :
184 : // Fourth part: sequential decays also of long-lived particles.
185 0 : if (doDecay) {
186 :
187 : // Loop through all entries to find those that should decay.
188 : int iDec = 0;
189 0 : do {
190 0 : Particle& decayer = event[iDec];
191 0 : if ( decayer.isFinal() && decayer.canDecay() && decayer.mayDecay() ) {
192 0 : decays.decay( iDec, event);
193 0 : if (decays.moreToDo()) moreToDo = true;
194 : }
195 0 : } while (++iDec < event.size());
196 0 : }
197 :
198 : // Normally done first time around, but sometimes not (e.g. Upsilon).
199 0 : } while (moreToDo);
200 :
201 : // Done.
202 0 : return true;
203 :
204 0 : }
205 :
206 : //--------------------------------------------------------------------------
207 :
208 : // Allow more decays if on/off switches changed.
209 : // Note: does not do sequential hadronization, e.g. for Upsilon.
210 :
211 : bool HadronLevel::moreDecays( Event& event) {
212 :
213 : // Colour-octet onia states must be decayed to singlet + gluon.
214 0 : if (!decayOctetOnia(event)) return false;
215 :
216 : // Loop through all entries to find those that should decay.
217 : int iDec = 0;
218 0 : do {
219 0 : if ( event[iDec].isFinal() && event[iDec].canDecay()
220 0 : && event[iDec].mayDecay() ) decays.decay( iDec, event);
221 0 : } while (++iDec < event.size());
222 :
223 : // Done.
224 : return true;
225 :
226 0 : }
227 :
228 : //--------------------------------------------------------------------------
229 :
230 : // Decay colour-octet onium states.
231 :
232 : bool HadronLevel::decayOctetOnia(Event& event) {
233 :
234 : // Loop over particles and decay any onia encountered.
235 0 : for (int iDec = 0; iDec < event.size(); ++iDec)
236 0 : if (event[iDec].isFinal()
237 0 : && particleDataPtr->isOctetHadron(event[iDec].id())) {
238 0 : if (!decays.decay( iDec, event)) return false;
239 :
240 : // Set colour flow by hand: gluon inherits octet-onium state.
241 0 : int iGlu = event.size() - 1;
242 0 : event[iGlu].cols( event[iDec].col(), event[iDec].acol() );
243 0 : }
244 :
245 : // Done.
246 0 : return true;
247 :
248 0 : }
249 :
250 : //--------------------------------------------------------------------------
251 :
252 : // Trace colour flow in the event to form colour singlet subsystems.
253 :
254 : bool HadronLevel::findSinglets(Event& event) {
255 :
256 : // Clear up storage.
257 0 : colConfig.clear();
258 :
259 : // Find a list of final partons and of all colour ends and gluons.
260 0 : if (colTrace.setupColList(event)) return true;
261 :
262 : // Begin arrange the partons into separate colour singlets.
263 :
264 : // Junctions: loop over them, and identify kind.
265 0 : for (int iJun = 0; iJun < event.sizeJunction(); ++iJun)
266 0 : if (event.remainsJunction(iJun)) {
267 0 : event.remainsJunction(iJun, false);
268 0 : int kindJun = event.kindJunction(iJun);
269 0 : iParton.resize(0);
270 :
271 : // Loop over junction legs.
272 0 : for (int iCol = 0; iCol < 3; ++iCol) {
273 0 : int indxCol = event.colJunction(iJun, iCol);
274 0 : iParton.push_back( -(10 + 10 * iJun + iCol) );
275 : // Junctions: find color ends.
276 0 : if (kindJun % 2 == 1 && !colTrace.traceFromAcol(indxCol, event, iJun,
277 0 : iCol, iParton)) return false;
278 : // Antijunctions: find anticolor ends.
279 0 : if (kindJun % 2 == 0 && !colTrace.traceFromCol(indxCol, event, iJun,
280 0 : iCol, iParton)) return false;
281 0 : }
282 :
283 : // A junction may be eliminated by insert if two quarks are nearby.
284 0 : int nJunOld = event.sizeJunction();
285 0 : if (!colConfig.insert(iParton, event)) return false;
286 0 : if (event.sizeJunction() < nJunOld) --iJun;
287 0 : }
288 :
289 : // Open strings: pick up each colour end and trace to its anticolor end.
290 0 : while (!colTrace.colFinished()) {
291 0 : iParton.resize(0);
292 0 : if (!colTrace.traceFromCol( -1, event, -1, -1, iParton)) return false;
293 :
294 : // Store found open string system. Analyze its properties.
295 0 : if (!colConfig.insert(iParton, event)) return false;
296 : }
297 :
298 : // Closed strings : begin at any gluon and trace until back at it.
299 0 : while (!colTrace.finished()) {
300 0 : iParton.resize(0);
301 0 : if (!colTrace.traceInLoop(event, iParton)) return false;
302 :
303 : // Store found closed string system. Analyze its properties.
304 0 : if (!colConfig.insert(iParton, event)) return false;
305 : }
306 :
307 : // Done.
308 0 : return true;
309 :
310 0 : }
311 :
312 : //==========================================================================
313 :
314 : } // end namespace Pythia8
|