Line data Source code
1 : // SigmaProcess.h 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 : // Header file for hard-process differential cross sections.
7 : // SigmaProcess: base class for cross sections.
8 : // Sigma0Process: base class for unresolved processes, derived from above.
9 : // Sigma1Process: base class for 2 -> 1 processes, derived from above.
10 : // Sigma2Process: base class for 2 -> 2 processes, derived from above.
11 : // Sigma3Process: base class for 2 -> 3 processes, derived from above.
12 : // SigmaLHAProcess: wrapper class for Les Houches Accord external input.
13 : // Actual physics processes are found in separate files:
14 : // SigmaQCD for QCD processes;
15 : // SigmaEW for electroweak processes (including photon production);
16 : // SigmaOnia for charmonium and bottomonium processes;
17 : // SigmaHiggs for Higgs processes;
18 : // SigmaSUSY for supersymmetric production;
19 : // SigmaLeftRightSym for processes in left-right-symmetric scenarios;
20 : // SigmaLeptoquark for leptoquark production.
21 : // SigmaExtraDim for processes in scenarios for extra dimensions;
22 : // SigmaGeneric for generic scalar/fermion/vector pair production.
23 :
24 : #ifndef Pythia8_SigmaProcess_H
25 : #define Pythia8_SigmaProcess_H
26 :
27 : #include "Pythia8/Basics.h"
28 : #include "Pythia8/BeamParticle.h"
29 : #include "Pythia8/Event.h"
30 : #include "Pythia8/Info.h"
31 : #include "Pythia8/LesHouches.h"
32 : #include "Pythia8/ParticleData.h"
33 : #include "Pythia8/PartonDistributions.h"
34 : #include "Pythia8/PythiaComplex.h"
35 : #include "Pythia8/PythiaStdlib.h"
36 : #include "Pythia8/ResonanceWidths.h"
37 : #include "Pythia8/Settings.h"
38 : #include "Pythia8/SigmaTotal.h"
39 : #include "Pythia8/StandardModel.h"
40 : #include "Pythia8/SLHAinterface.h"
41 : #include "Pythia8/SusyLesHouches.h"
42 :
43 : namespace Pythia8 {
44 :
45 : //==========================================================================
46 :
47 : // InBeam is a simple helper class for partons and their flux in a beam.
48 :
49 : class InBeam {
50 :
51 : public:
52 :
53 : // Constructor.
54 0 : InBeam( int idIn = 0) : id(idIn), pdf(0.) {}
55 :
56 : // Values.
57 : int id;
58 : double pdf;
59 :
60 : };
61 :
62 : //==========================================================================
63 :
64 : // InPair is a simple helper class for colliding parton pairs and their flux.
65 :
66 : class InPair {
67 :
68 : public:
69 :
70 : // Constructor.
71 0 : InPair( int idAIn = 0, int idBIn = 0) : idA(idAIn), idB(idBIn),
72 0 : pdfA(0.), pdfB(0.), pdfSigma(0.) {}
73 :
74 : // Values.
75 : int idA, idB;
76 : double pdfA, pdfB, pdfSigma;
77 :
78 : };
79 :
80 : //==========================================================================
81 :
82 : // SigmaProcess is the base class for cross section calculations.
83 :
84 : class SigmaProcess {
85 :
86 : public:
87 :
88 : // Destructor.
89 0 : virtual ~SigmaProcess() {}
90 :
91 : // Perform simple initialization and store pointers.
92 : void init(Info* infoPtrIn, Settings* settingsPtrIn,
93 : ParticleData* particleDataPtrIn, Rndm* rndmPtrIn,
94 : BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn, Couplings* couplings,
95 : SigmaTotal* sigmaTotPtrIn = 0, SLHAinterface* slhaInterfacePtrIn = 0);
96 :
97 : // Store or replace Les Houches pointer.
98 0 : void setLHAPtr( LHAup* lhaUpPtrIn) {lhaUpPtr = lhaUpPtrIn;}
99 :
100 : // Initialize process. Only used for some processes.
101 0 : virtual void initProc() {}
102 :
103 : // Set up allowed flux of incoming partons. Default is no flux.
104 : virtual bool initFlux();
105 :
106 : // Input and complement kinematics for resolved 2 -> 1 process.
107 : // Usage: set1Kin( x1in, x2in, sHin).
108 0 : virtual void set1Kin( double , double , double ) {}
109 :
110 : // Input and complement kinematics for resolved 2 -> 2 process.
111 : // Usage: set2Kin( x1in, x2in, sHin, tHin, m3in, m4in, runBW3in, runBW4in).
112 : virtual void set2Kin( double , double , double , double , double ,
113 0 : double, double, double ) {}
114 :
115 : // Ditto, but for Multiparton Interactions applications, so different input.
116 : // Usage: set2KinMPI( x1in, x2in, sHin, tHin, uHin,
117 : // alpSin, alpEMin, needMasses, m3in, m4in)
118 : virtual void set2KinMPI( double , double , double , double ,
119 0 : double , double , double , bool , double , double ) {}
120 :
121 : // Input and complement kinematics for resolved 2 -> 3 process.
122 : // Usage: set3Kin( x1in, x2in, sHin, p3prel, p4prel, p5prel,
123 : // m3in, m4in, m5in, runBW3in, runBW4in, runBW5in);
124 : virtual void set3Kin( double , double , double , Vec4 , Vec4 , Vec4 ,
125 0 : double , double , double , double , double , double ) {}
126 :
127 : // Calculate flavour-independent parts of cross section.
128 0 : virtual void sigmaKin() {}
129 :
130 : // Evaluate sigma for unresolved, sigmaHat(sHat) for 2 -> 1 processes,
131 : // d(sigmaHat)/d(tHat) for (resolved) 2 -> 2 processes, and |M|^2 for
132 : // 2 -> 3 processes. Answer in "native" units, either mb or GeV^-2.
133 0 : virtual double sigmaHat() {return 0.;}
134 :
135 : // Wrapper to sigmaHat, to (a) store current incoming flavours and
136 : // (b) convert from GeV^-2 to mb where required.
137 : // For 2 -> 1/2 also (c) convert from from |M|^2 to d(sigmaHat)/d(tHat).
138 : virtual double sigmaHatWrap(int id1in = 0, int id2in = 0) {
139 0 : id1 = id1in; id2 = id2in;
140 0 : return ( convert2mb() ? CONVERT2MB * sigmaHat() : sigmaHat() ); }
141 :
142 : // Convolute above with parton flux and K factor. Sum over open channels.
143 : virtual double sigmaPDF();
144 :
145 : // Select incoming parton channel and extract parton densities (resolved).
146 : void pickInState(int id1in = 0, int id2in = 0);
147 :
148 : // Select flavour, colour and anticolour.
149 0 : virtual void setIdColAcol() {}
150 :
151 : // Perform kinematics for a Multiparton Interaction, in its rest frame.
152 : virtual bool final2KinMPI( int = 0, int = 0, Vec4 = 0., Vec4 = 0.,
153 0 : double = 0., double = 0.) {return true;}
154 :
155 : // Evaluate weight for simultaneous flavours (only gamma*/Z0 gamma*/Z0).
156 : // Usage: weightDecayFlav( process).
157 0 : virtual double weightDecayFlav( Event&) {return 1.;}
158 :
159 : // Evaluate weight for decay angular configuration.
160 : // Usage: weightDecay( process, iResBeg, iResEnd), where
161 : // iResBeg <= i < iResEnd is range of sister partons to test decays of.
162 0 : virtual double weightDecay( Event&, int, int) {return 1.;}
163 :
164 : // Set scale, when that is missing for an external LHA process.
165 0 : virtual void setScale() {}
166 :
167 : // Process name and code, and the number of final-state particles.
168 0 : virtual string name() const {return "unnamed process";}
169 0 : virtual int code() const {return 0;}
170 0 : virtual int nFinal() const {return 2;}
171 :
172 : // Need to know which incoming partons to set up interaction for.
173 0 : virtual string inFlux() const {return "unknown";}
174 :
175 : // Need to know whether to convert cross section answer from GeV^-2 to mb.
176 0 : virtual bool convert2mb() const {return true;}
177 :
178 : // For 2 -> 2 process optional conversion from |M|^2 to d(sigmaHat)/d(tHat).
179 0 : virtual bool convertM2() const {return false;}
180 :
181 : // Special treatment needed for Les Houches processes.
182 0 : virtual bool isLHA() const {return false;}
183 :
184 : // Special treatment needed for elastic and diffractive processes.
185 0 : virtual bool isNonDiff() const {return false;}
186 0 : virtual bool isResolved() const {return true;}
187 0 : virtual bool isDiffA() const {return false;}
188 0 : virtual bool isDiffB() const {return false;}
189 0 : virtual bool isDiffC() const {return false;}
190 :
191 : // Special treatment needed for SUSY processes.
192 0 : virtual bool isSUSY() const {return false;}
193 :
194 : // Special treatment needed if negative cross sections allowed.
195 0 : virtual bool allowNegativeSigma() const {return false;}
196 :
197 : // Flavours in 2 -> 2/3 processes where masses needed from beginning.
198 : // (For a light quark masses will be used in the final kinematics,
199 : // but not at the matrix-element level. For a gluon no masses at all.)
200 0 : virtual int id3Mass() const {return 0;}
201 0 : virtual int id4Mass() const {return 0;}
202 0 : virtual int id5Mass() const {return 0;}
203 :
204 : // Special treatment needed if process contains an s-channel resonance.
205 0 : virtual int resonanceA() const {return 0;}
206 0 : virtual int resonanceB() const {return 0;}
207 :
208 : // 2 -> 2 and 2 -> 3 processes only through s-channel exchange.
209 0 : virtual bool isSChannel() const {return false;}
210 :
211 : // NOAM: Insert an intermediate resonance in 2 -> 1 -> 2 (or 3) listings.
212 0 : virtual int idSChannel() const {return 0;}
213 :
214 : // QCD 2 -> 3 processes need special phase space selection machinery.
215 0 : virtual bool isQCD3body() const {return false;}
216 :
217 : // Special treatment in 2 -> 3 with two massive propagators.
218 0 : virtual int idTchan1() const {return 0;}
219 0 : virtual int idTchan2() const {return 0;}
220 0 : virtual double tChanFracPow1() const {return 0.3;}
221 0 : virtual double tChanFracPow2() const {return 0.3;}
222 0 : virtual bool useMirrorWeight() const {return false;}
223 :
224 : // Special process-specific gamma*/Z0 choice if >=0 (e.g. f fbar -> H0 Z0).
225 0 : virtual int gmZmode() const {return -1;}
226 :
227 : // Tell whether tHat and uHat are swapped (= same as swap 3 and 4).
228 0 : bool swappedTU() const {return swapTU;}
229 :
230 : // Give back particle properties: flavours, colours, masses, or all.
231 0 : int id(int i) const {return idSave[i];}
232 0 : int col(int i) const {return colSave[i];}
233 0 : int acol(int i) const {return acolSave[i];}
234 0 : double m(int i) const {return mSave[i];}
235 0 : Particle getParton(int i) const {return parton[i];}
236 :
237 : // Give back couplings and parton densities.
238 : // Not all known for nondiffractive.
239 0 : double Q2Ren() const {return Q2RenSave;}
240 0 : double alphaEMRen() const {return alpEM;}
241 0 : double alphaSRen() const {return alpS;}
242 0 : double Q2Fac() const {return Q2FacSave;}
243 0 : double pdf1() const {return pdf1Save;}
244 0 : double pdf2() const {return pdf2Save;}
245 :
246 : // Give back angles; relevant only for multipe-interactions processes.
247 0 : double thetaMPI() const {return atan2( sinTheta, cosTheta);}
248 0 : double phiMPI() const {return phi;}
249 0 : double sHBetaMPI() const {return sHBeta;}
250 : double pT2MPI() const {return pT2Mass;}
251 0 : double pTMPIFin() const {return pTFin;}
252 :
253 : // Save and load kinematics for trial interactions
254 : void saveKin() {
255 0 : for (int i = 0; i < 12; i++) { partonT[i] = parton[i];
256 0 : mSaveT[i] = mSave[i]; }
257 0 : pTFinT = pTFin; phiT = phi; cosThetaT = cosTheta; sinThetaT = sinTheta; }
258 : void loadKin() {
259 : for (int i = 0; i < 12; i++) { parton[i] = partonT[i];
260 : mSave[i] = mSaveT[i]; }
261 : pTFin = pTFinT; cosTheta = cosThetaT; sinTheta = sinThetaT; phi = phiT;
262 : }
263 : void swapKin() {
264 0 : for (int i = 0; i < 12; i++) { swap(parton[i], partonT[i]);
265 0 : swap(mSave[i], mSaveT[i]); }
266 0 : swap(pTFin, pTFinT); swap(cosTheta, cosThetaT);
267 0 : swap(sinTheta, sinThetaT); swap(phi, phiT); }
268 :
269 : protected:
270 :
271 : // Constructor.
272 0 : SigmaProcess() : infoPtr(0), settingsPtr(0), particleDataPtr(0),
273 0 : rndmPtr(0), beamAPtr(0), beamBPtr(0), couplingsPtr(0), sigmaTotPtr(0),
274 0 : slhaPtr(0), lhaUpPtr(0) {for (int i = 0; i < 12; ++i) mSave[i] = 0.;
275 0 : Q2RenSave = alpEM = alpS = Q2FacSave = pdf1Save = pdf2Save = 0.; }
276 :
277 : // Constants: could only be changed in the code itself.
278 : static const double CONVERT2MB, MASSMARGIN, COMPRELERR;
279 : static const int NCOMPSTEP;
280 :
281 : // Pointer to various information on the generation.
282 : Info* infoPtr;
283 :
284 : // Pointer to the settings database.
285 : Settings* settingsPtr;
286 :
287 : // Pointer to the particle data table.
288 : ParticleData* particleDataPtr;
289 :
290 : // Pointer to the random number generator.
291 : Rndm* rndmPtr;
292 :
293 : // Pointers to incoming beams.
294 : BeamParticle* beamAPtr;
295 : BeamParticle* beamBPtr;
296 :
297 : // Pointer to Standard Model couplings, including alphaS and alphaEM.
298 : Couplings* couplingsPtr;
299 :
300 : // Pointer to the total/elastic/diffractive cross section object.
301 : SigmaTotal* sigmaTotPtr;
302 :
303 : // Pointer to an SLHA object.
304 : SusyLesHouches* slhaPtr;
305 :
306 : // Pointer to LHAup for generating external events.
307 : LHAup* lhaUpPtr;
308 :
309 : // Initialization data, normally only set once.
310 : int nQuarkIn, renormScale1, renormScale2, renormScale3, renormScale3VV,
311 : factorScale1, factorScale2, factorScale3, factorScale3VV;
312 : double Kfactor, mcME, mbME, mmuME, mtauME, renormMultFac, renormFixScale,
313 : factorMultFac, factorFixScale;
314 :
315 : // CP violation parameters for Higgs sector, normally only set once.
316 : int higgsH1parity, higgsH2parity, higgsA3parity;
317 : double higgsH1eta, higgsH2eta, higgsA3eta, higgsH1phi, higgsH2phi,
318 : higgsA3phi;
319 :
320 : // Information on incoming beams.
321 : int idA, idB;
322 : double mA, mB;
323 : bool isLeptonA, isLeptonB, hasLeptonBeams;
324 :
325 : // Partons in beams, with PDF's.
326 : vector<InBeam> inBeamA;
327 : vector<InBeam> inBeamB;
328 0 : void addBeamA(int idIn) {inBeamA.push_back(InBeam(idIn));}
329 0 : void addBeamB(int idIn) {inBeamB.push_back(InBeam(idIn));}
330 0 : int sizeBeamA() const {return inBeamA.size();}
331 0 : int sizeBeamB() const {return inBeamB.size();}
332 :
333 : // Allowed colliding parton pairs, with pdf's.
334 : vector<InPair> inPair;
335 : void addPair(int idAIn, int idBIn) {
336 0 : inPair.push_back(InPair(idAIn, idBIn));}
337 0 : int sizePair() const {return inPair.size();}
338 :
339 : // Store common subprocess kinematics quantities.
340 : double mH, sH, sH2;
341 :
342 : // Store Q2 renormalization and factorization scales, and related values.
343 : double Q2RenSave, alpEM, alpS, Q2FacSave, x1Save, x2Save, pdf1Save,
344 : pdf2Save, sigmaSumSave;
345 :
346 : // Store flavour, colour, anticolour, mass, angles and the whole particle.
347 : int id1, id2, id3, id4, id5;
348 : int idSave[12], colSave[12], acolSave[12];
349 : double mSave[12], cosTheta, sinTheta, phi, sHMass, sHBeta, pT2Mass, pTFin;
350 : Particle parton[12];
351 :
352 : // Minimal set of saved kinematics for trial interactions when
353 : // using the x-dependent matter profile of multiparton interactions.
354 : Particle partonT[12];
355 : double mSaveT[12], pTFinT, cosThetaT, sinThetaT, phiT;
356 :
357 : // Calculate and store all modified masses and four-vectors
358 : // intended for matrix elements. Return false if failed.
359 0 : virtual bool setupForME() {return true;}
360 : bool setupForMEin();
361 : double mME[12];
362 : Vec4 pME[12];
363 :
364 : // Store whether tHat and uHat are swapped (= same as swap 3 and 4).
365 : bool swapTU;
366 :
367 : // Set flavour, colour and anticolour.
368 : void setId( int id1in = 0, int id2in = 0, int id3in = 0, int id4in = 0,
369 0 : int id5in = 0) {idSave[1] = id1in; idSave[2] = id2in; idSave[3] = id3in;
370 0 : idSave[4] = id4in; idSave[5] = id5in;}
371 : void setColAcol( int col1 = 0, int acol1 = 0,
372 : int col2 = 0, int acol2 = 0, int col3 = 0, int acol3 = 0,
373 : int col4 = 0, int acol4 = 0, int col5 = 0, int acol5 = 0) {
374 0 : colSave[1] = col1; acolSave[1] = acol1; colSave[2] = col2;
375 0 : acolSave[2] = acol2; colSave[3] = col3; acolSave[3] = acol3;
376 0 : colSave[4] = col4; acolSave[4] = acol4; colSave[5] = col5;
377 0 : acolSave[5] = acol5; }
378 0 : void swapColAcol() { swap(colSave[1], acolSave[1]);
379 0 : swap(colSave[2], acolSave[2]); swap(colSave[3], acolSave[3]);
380 0 : swap(colSave[4], acolSave[4]); swap(colSave[5], acolSave[5]);}
381 0 : void swapCol1234() { swap(colSave[1], colSave[2]);
382 0 : swap(colSave[3], colSave[4]); swap(acolSave[1], acolSave[2]);
383 0 : swap(acolSave[3], acolSave[4]);}
384 0 : void swapCol12() { swap(colSave[1], colSave[2]);
385 0 : swap(acolSave[1], acolSave[2]);}
386 0 : void swapCol34() { swap(colSave[3], colSave[4]);
387 0 : swap(acolSave[3], acolSave[4]);}
388 :
389 : // Common code for top and Higgs secondary decay angular weights.
390 : double weightTopDecay( Event& process, int iResBeg, int iResEnd);
391 : double weightHiggsDecay( Event& process, int iResBeg, int iResEnd);
392 :
393 : };
394 :
395 : //==========================================================================
396 :
397 : // Sigma0Process is the base class for unresolved and minimum-bias processes.
398 : // It is derived from SigmaProcess.
399 :
400 : class Sigma0Process : public SigmaProcess {
401 :
402 : public:
403 :
404 : // Destructor.
405 0 : virtual ~Sigma0Process() {}
406 :
407 : // Number of final-state particles.
408 0 : virtual int nFinal() const {return 2;};
409 :
410 : // No partonic flux to be set up.
411 0 : virtual bool initFlux() {return true;}
412 :
413 : // Evaluate sigma for unresolved processes.
414 0 : virtual double sigmaHat() {return 0.;}
415 :
416 : // Since no PDF's there is no difference from above.
417 0 : virtual double sigmaPDF() {return sigmaHat();}
418 :
419 : // Answer for these processes already in mb, so do not convert.
420 0 : virtual bool convert2mb() const {return false;}
421 :
422 : protected:
423 :
424 : // Constructor.
425 0 : Sigma0Process() {}
426 :
427 : };
428 :
429 : //==========================================================================
430 :
431 : // Sigma1Process is the base class for 2 -> 1 processes.
432 : // It is derived from SigmaProcess.
433 :
434 : class Sigma1Process : public SigmaProcess {
435 :
436 : public:
437 :
438 : // Destructor.
439 0 : virtual ~Sigma1Process() {}
440 :
441 : // Number of final-state particles.
442 0 : virtual int nFinal() const {return 1;};
443 :
444 : // Input and complement kinematics for resolved 2 -> 1 process.
445 : virtual void set1Kin( double x1in, double x2in, double sHin) {
446 0 : store1Kin( x1in, x2in, sHin); sigmaKin();}
447 :
448 : // Evaluate sigmaHat(sHat) for resolved 2 -> 1 processes.
449 0 : virtual double sigmaHat() {return 0.;}
450 :
451 : // Wrapper to sigmaHat, to (a) store current incoming flavours,
452 : // (b) convert from GeV^-2 to mb where required, and
453 : // (c) convert from |M|^2 to d(sigmaHat)/d(tHat) where required.
454 : virtual double sigmaHatWrap(int id1in = 0, int id2in = 0);
455 :
456 : protected:
457 :
458 : // Constructor.
459 0 : Sigma1Process() {}
460 :
461 : // Store kinematics and set scales for resolved 2 -> 1 process.
462 : virtual void store1Kin( double x1in, double x2in, double sHin);
463 :
464 : // Calculate modified masses and four-vectors for matrix elements.
465 : virtual bool setupForME();
466 :
467 : };
468 :
469 : //==========================================================================
470 :
471 : // Sigma2Process is the base class for 2 -> 2 processes.
472 : // It is derived from SigmaProcess.
473 :
474 : class Sigma2Process : public SigmaProcess {
475 :
476 : public:
477 :
478 : // Destructor.
479 0 : virtual ~Sigma2Process() {}
480 :
481 : // Number of final-state particles.
482 0 : virtual int nFinal() const {return 2;};
483 :
484 : // Input and complement kinematics for resolved 2 -> 2 process.
485 : virtual void set2Kin( double x1in, double x2in, double sHin,
486 : double tHin, double m3in, double m4in, double runBW3in,
487 0 : double runBW4in) { store2Kin( x1in, x2in, sHin, tHin, m3in, m4in,
488 0 : runBW3in, runBW4in); sigmaKin();}
489 :
490 : // Ditto, but for Multiparton Interactions applications, so different input.
491 : virtual void set2KinMPI( double x1in, double x2in, double sHin,
492 : double tHin, double uHin, double alpSin, double alpEMin,
493 : bool needMasses, double m3in, double m4in) {
494 0 : store2KinMPI( x1in, x2in, sHin, tHin, uHin, alpSin, alpEMin,
495 0 : needMasses, m3in, m4in); sigmaKin();}
496 :
497 : // Evaluate d(sigmaHat)/d(tHat) for resolved 2 -> 2 processes.
498 0 : virtual double sigmaHat() {return 0.;}
499 :
500 : // Wrapper to sigmaHat, to (a) store current incoming flavours,
501 : // (b) convert from GeV^-2 to mb where required, and
502 : // (c) convert from |M|^2 to d(sigmaHat)/d(tHat) where required.
503 : virtual double sigmaHatWrap(int id1in = 0, int id2in = 0) {
504 0 : id1 = id1in; id2 = id2in; double sigmaTmp = sigmaHat();
505 0 : if (convertM2()) sigmaTmp /= 16. * M_PI * sH2;
506 0 : if (convert2mb()) sigmaTmp *= CONVERT2MB; return sigmaTmp;}
507 :
508 : // Perform kinematics for a Multiparton Interaction, in its rest frame.
509 : virtual bool final2KinMPI( int i1Res = 0, int i2Res = 0, Vec4 p1Res = 0.,
510 : Vec4 p2Res = 0., double m1Res = 0., double m2Res = 0.);
511 :
512 : protected:
513 :
514 : // Constructor.
515 0 : Sigma2Process() : tH(0.), uH(0.), tH2(0.), uH2(0.), m3(0.), s3(0.),
516 0 : m4(0.), s4(0.), pT2(0.), runBW3(0.), runBW4(0.) {}
517 :
518 : // Store kinematics and set scales for resolved 2 -> 2 process.
519 : virtual void store2Kin( double x1in, double x2in, double sHin,
520 : double tHin, double m3in, double m4in, double runBW3in,
521 : double runBW4in);
522 : virtual void store2KinMPI( double x1in, double x2in, double sHin,
523 : double tHin, double uHin, double alpSin, double alpEMin,
524 : bool needMasses, double m3in, double m4in);
525 :
526 : // Calculate modified masses and four-vectors for matrix elements.
527 : virtual bool setupForME();
528 :
529 : // Store subprocess kinematics quantities.
530 : double tH, uH, tH2, uH2, m3, s3, m4, s4, pT2, runBW3, runBW4;
531 :
532 : };
533 :
534 : //==========================================================================
535 :
536 : // Sigma3Process is the base class for 2 -> 3 processes.
537 : // It is derived from SigmaProcess.
538 :
539 : class Sigma3Process : public SigmaProcess {
540 :
541 : public:
542 :
543 : // Destructor.
544 0 : virtual ~Sigma3Process() {}
545 :
546 : // Number of final-state particles.
547 0 : virtual int nFinal() const {return 3;};
548 :
549 : // Input and complement kinematics for resolved 2 -> 3 process.
550 : virtual void set3Kin( double x1in, double x2in, double sHin,
551 : Vec4 p3cmIn, Vec4 p4cmIn, Vec4 p5cmIn, double m3in, double m4in,
552 : double m5in, double runBW3in, double runBW4in, double runBW5in) {
553 0 : store3Kin( x1in, x2in, sHin, p3cmIn, p4cmIn, p5cmIn, m3in, m4in, m5in,
554 0 : runBW3in, runBW4in, runBW5in); sigmaKin();}
555 :
556 : // Evaluate d(sigmaHat)/d(tHat) for resolved 2 -> 3 processes.
557 0 : virtual double sigmaHat() {return 0.;}
558 :
559 : protected:
560 :
561 : // Constructor.
562 0 : Sigma3Process() {}
563 :
564 : // Store kinematics and set scales for resolved 2 -> 3 process.
565 : virtual void store3Kin( double x1in, double x2in, double sHin,
566 : Vec4 p3cmIn, Vec4 p4cmIn, Vec4 p5cmIn, double m3in, double m4in,
567 : double m5in, double runBW3in, double runBW4in, double runBW5in);
568 :
569 : // Calculate modified masses and four-vectors for matrix elements.
570 : virtual bool setupForME();
571 :
572 : // Store subprocess kinematics quantities.
573 : double m3, s3, m4, s4, m5, s5, runBW3, runBW4, runBW5;
574 : Vec4 p3cm, p4cm, p5cm;
575 :
576 : };
577 :
578 : //==========================================================================
579 :
580 : // SigmaLHAProcess is a wrapper class for Les Houches Accord external input.
581 : // It is derived from SigmaProcess.
582 :
583 : class SigmaLHAProcess : public SigmaProcess {
584 :
585 : public:
586 :
587 : // Constructor.
588 0 : SigmaLHAProcess() {}
589 :
590 : // Destructor.
591 0 : virtual ~SigmaLHAProcess() {}
592 :
593 : // No partonic flux to be set up.
594 0 : virtual bool initFlux() {return true;}
595 :
596 : // Dummy function: action is put in PhaseSpaceLHA.
597 0 : virtual double sigmaPDF() {return 1.;}
598 :
599 : // Evaluate weight for decay angular configuration, where relevant.
600 : virtual double weightDecay( Event& process, int iResBeg, int iResEnd);
601 :
602 : // Set scale, when that is missing for an external LHA process.
603 : virtual void setScale();
604 :
605 : // Info on the subprocess.
606 0 : virtual string name() const {return "Les Houches User Process(es)";}
607 0 : virtual int code() const {return 9999;}
608 :
609 : // Number of final-state particles depends on current process choice.
610 : virtual int nFinal() const;
611 :
612 : // Answer for these processes not in GeV^-2, so do not do this conversion.
613 0 : virtual bool convert2mb() const {return false;}
614 :
615 : // Ensure special treatment of Les Houches processes.
616 0 : virtual bool isLHA() const {return true;}
617 :
618 : // Special treatment needed if negative cross sections allowed.
619 : virtual bool allowNegativeSigma() const {
620 0 : return (lhaUpPtr->strategy() < 0);}
621 :
622 : private:
623 :
624 : };
625 :
626 : //==========================================================================
627 :
628 : } // end namespace Pythia8
629 :
630 : #endif // Pythia8_SigmaProcess_H
|
