1 // ProcessContainer.cc 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.
6 // Function definitions (not found in the header) for the
7 // ProcessContainer and SetupContainers classes.
9 #include "ProcessContainer.h"
11 // Internal headers for special processes.
12 #include "SigmaCompositeness.h"
14 #include "SigmaExtraDim.h"
15 #include "SigmaGeneric.h"
16 #include "SigmaHiggs.h"
17 #include "SigmaLeftRightSym.h"
18 #include "SigmaLeptoquark.h"
19 #include "SigmaNewGaugeBosons.h"
20 #include "SigmaOnia.h"
22 #include "SigmaSUSY.h"
26 //==========================================================================
28 // ProcessContainer class.
29 // Information allowing the generation of a specific process.
31 //--------------------------------------------------------------------------
33 // Constants: could be changed here if desired, but normally should not.
34 // These are of technical nature, as described for each.
36 // Number of event tries to check maximization finding reliability.
37 const int ProcessContainer::N12SAMPLE = 100;
39 // Ditto, but increased for 2 -> 3 processes.
40 const int ProcessContainer::N3SAMPLE = 1000;
42 //--------------------------------------------------------------------------
44 // Initialize phase space and counters.
45 // Argument isFirst distinguishes two hard processes in same event.
47 bool ProcessContainer::init(bool isFirst, Info* infoPtrIn,
48 Settings& settings, ParticleData* particleDataPtrIn, Rndm* rndmPtrIn,
49 BeamParticle* beamAPtr, BeamParticle* beamBPtr, Couplings* couplingsPtr,
50 SigmaTotal* sigmaTotPtr, ResonanceDecays* resDecaysPtrIn,
51 SusyLesHouches* slhaPtr, UserHooks* userHooksPtrIn) {
53 // Extract info about current process from SigmaProcess object.
54 isLHA = sigmaProcessPtr->isLHA();
55 isMinBias = sigmaProcessPtr->isMinBias();
56 isResolved = sigmaProcessPtr->isResolved();
57 isDiffA = sigmaProcessPtr->isDiffA();
58 isDiffB = sigmaProcessPtr->isDiffB();
59 isDiffC = sigmaProcessPtr->isDiffC();
60 isQCD3body = sigmaProcessPtr->isQCD3body();
61 int nFin = sigmaProcessPtr->nFinal();
62 lhaStrat = (isLHA) ? lhaUpPtr->strategy() : 0;
63 lhaStratAbs = abs(lhaStrat);
64 allowNegSig = sigmaProcessPtr->allowNegativeSigma();
66 // Flag for maximum violation handling.
67 increaseMaximum = settings.flag("PhaseSpace:increaseMaximum");
69 // Pick and create phase space generator. Send pointers where required.
70 if (isLHA) phaseSpacePtr = new PhaseSpaceLHA();
71 else if (isMinBias) phaseSpacePtr = new PhaseSpace2to2minbias();
72 else if (!isResolved && !isDiffA && !isDiffB && !isDiffC )
73 phaseSpacePtr = new PhaseSpace2to2elastic();
74 else if (!isResolved && !isDiffA && !isDiffB && isDiffC)
75 phaseSpacePtr = new PhaseSpace2to3diffractive();
76 else if (!isResolved) phaseSpacePtr = new PhaseSpace2to2diffractive(
78 else if (nFin == 1) phaseSpacePtr = new PhaseSpace2to1tauy();
79 else if (nFin == 2) phaseSpacePtr = new PhaseSpace2to2tauyz();
80 else if (isQCD3body) phaseSpacePtr = new PhaseSpace2to3yyycyl();
81 else phaseSpacePtr = new PhaseSpace2to3tauycyl();
83 // Store pointers and perform simple initialization.
85 particleDataPtr = particleDataPtrIn;
87 resDecaysPtr = resDecaysPtrIn;
88 userHooksPtr = userHooksPtrIn;
89 canVetoResDecay = (userHooksPtr != 0)
90 ? userHooksPtr->canVetoResonanceDecays() : false;
92 sigmaProcessPtr->setLHAPtr(lhaUpPtr);
93 phaseSpacePtr->setLHAPtr(lhaUpPtr);
95 sigmaProcessPtr->init(infoPtr, &settings, particleDataPtr, rndmPtr,
96 beamAPtr, beamBPtr, couplingsPtr, sigmaTotPtr, slhaPtr);
97 phaseSpacePtr->init( isFirst, sigmaProcessPtr, infoPtr, &settings,
98 particleDataPtr, rndmPtr, beamAPtr, beamBPtr, couplingsPtr, sigmaTotPtr,
101 // Reset cross section statistics.
115 // Initialize process and allowed incoming partons.
116 sigmaProcessPtr->initProc();
117 if (!sigmaProcessPtr->initFlux()) return false;
119 // Find maximum of differential cross section * phasespace.
120 bool physical = phaseSpacePtr->setupSampling();
121 sigmaMx = phaseSpacePtr->sigmaMax();
122 double sigmaHalfWay = sigmaMx;
124 // Separate signed maximum needed for LHA with negative weight.
125 sigmaSgn = phaseSpacePtr->sigmaSumSigned();
127 // Check maximum by a few events, and extrapolate a further increase.
128 if (physical & !isLHA) {
129 int nSample = (nFin < 3) ? N12SAMPLE : N3SAMPLE;
130 for (int iSample = 0; iSample < nSample; ++iSample) {
132 while (!test) test = phaseSpacePtr->trialKin(false);
133 if (iSample == nSample/2) sigmaHalfWay = phaseSpacePtr->sigmaMax();
135 double sigmaFullWay = phaseSpacePtr->sigmaMax();
136 sigmaMx = (sigmaHalfWay > 0.) ? pow2(sigmaFullWay) / sigmaHalfWay
138 phaseSpacePtr->setSigmaMax(sigmaMx);
145 //--------------------------------------------------------------------------
147 // Generate a trial event; selected or not.
149 bool ProcessContainer::trialProcess() {
151 // Loop over tries only occurs for Les Houches strategy = +-2.
152 for (int iTry = 0; ; ++iTry) {
154 // Generate a trial phase space point, if meaningful.
155 if (sigmaMx == 0.) return false;
156 infoPtr->setEndOfFile(false);
157 bool repeatSame = (iTry > 0);
158 bool physical = phaseSpacePtr->trialKin(true, repeatSame);
160 // Possibly fail, e.g. if at end of Les Houches file, else cross section.
161 if (isLHA && !physical) infoPtr->setEndOfFile(true);
163 if (!physical) return false;
164 double sigmaNow = phaseSpacePtr->sigmaNow();
166 // Tell if this event comes with weight from cross section.
167 double sigmaWeight = 1.;
168 if (!isLHA && !increaseMaximum && sigmaNow > sigmaMx)
169 sigmaWeight = sigmaNow / sigmaMx;
170 if ( lhaStrat < 0 && sigmaNow < 0.) sigmaWeight = -1.;
171 if ( lhaStratAbs == 4) sigmaWeight = sigmaNow;
173 // Also compensating weight from biased phase-space selection.
174 double biasWeight = phaseSpacePtr->biasSelectionWeight();
175 weightNow = sigmaWeight * biasWeight;
176 infoPtr->setWeight( weightNow, lhaStrat);
178 // Check that not negative cross section when not allowed.
180 if (sigmaNow < sigmaNeg) {
181 infoPtr->errorMsg("Warning in ProcessContainer::trialProcess: neg"
182 "ative cross section set 0", "for " + sigmaProcessPtr->name() );
185 if (sigmaNow < 0.) sigmaNow = 0.;
188 // Cross section of process may come with a weight. Update sum.
189 double sigmaAdd = sigmaNow * biasWeight;
190 if (lhaStratAbs == 2 || lhaStratAbs == 3) sigmaAdd = sigmaSgn;
191 sigmaSum += sigmaAdd;
192 sigma2Sum += pow2(sigmaAdd);
194 // Check if maximum violated.
195 newSigmaMx = phaseSpacePtr->newSigmaMax();
196 if (newSigmaMx) sigmaMx = phaseSpacePtr->sigmaMax();
198 // Select or reject trial point.
200 if (lhaStratAbs < 3) select
201 = (newSigmaMx || rndmPtr->flat() * abs(sigmaMx) < abs(sigmaNow));
203 if (select || lhaStratAbs != 2) return select;
208 //--------------------------------------------------------------------------
210 // Pick flavours and colour flow of process.
212 bool ProcessContainer::constructState() {
214 // Construct flavour and colours for selected event.
215 if (isResolved && !isMinBias) sigmaProcessPtr->pickInState();
216 sigmaProcessPtr->setIdColAcol();
223 //--------------------------------------------------------------------------
225 // Give the hard subprocess with kinematics.
227 bool ProcessContainer::constructProcess( Event& process, bool isHardest) {
229 // Construct kinematics from selected phase space point.
230 if (!phaseSpacePtr->finalKin()) return false;
231 int nFin = sigmaProcessPtr->nFinal();
233 // Basic info on process.
234 if (isHardest) infoPtr->setType( name(), code(), nFin, isMinBias,
235 isResolved, isDiffA, isDiffB, isDiffC, isLHA);
237 // Let hard process record begin with the event as a whole and
238 // the two incoming beam particles.
239 process.append( 90, -11, 0, 0, 0, 0, 0, 0,
240 Vec4(0., 0., 0., infoPtr->eCM()), infoPtr->eCM(), 0. );
241 process.append( infoPtr->idA(), -12, 0, 0, 0, 0, 0, 0,
242 Vec4(0., 0., infoPtr->pzA(), infoPtr->eA()), infoPtr->mA(), 0. );
243 process.append( infoPtr->idB(), -12, 0, 0, 0, 0, 0, 0,
244 Vec4(0., 0., infoPtr->pzB(), infoPtr->eB()), infoPtr->mB(), 0. );
246 // For minbias process no interaction selected so far, so done.
247 if (isMinBias) return true;
249 // Entries 3 and 4, now to be added, come from 1 and 2.
250 process[1].daughter1(3);
251 process[2].daughter1(4);
254 // For DiffC entries 3 - 5 come jointly from 1 and 2 (to keep HepMC happy).
256 process[1].daughters(3, 5);
257 process[2].daughters(3, 5);
260 // Insert the subprocess partons - resolved processes.
261 int idRes = sigmaProcessPtr->idSChannel();
262 if (isResolved && !isLHA) {
264 // NOAM: Mothers and daughters without/with intermediate state.
268 int m_D2 = (nFin == 1) ? 0 : 4 + nFin;
276 // Find scale from which to begin MPI/ISR/FSR evolution.
277 scale = sqrt(Q2Fac());
278 process.scale( scale );
280 // Loop over incoming and outgoing partons.
281 int colOffset = process.lastColTag();
282 for (int i = 1; i <= 2 + nFin; ++i) {
284 // Read out particle info from SigmaProcess object.
285 int id = sigmaProcessPtr->id(i);
286 int status = (i <= 2) ? -21 : 23;
287 int mother1 = (i <= 2) ? i : m_M1 ;
288 int mother2 = (i <= 2) ? 0 : m_M2 ;
289 int daughter1 = (i <= 2) ? m_D1 : 0;
290 int daughter2 = (i <= 2) ? m_D2 : 0;
291 int col = sigmaProcessPtr->col(i);
292 if (col > 0) col += colOffset;
293 else if (col < 0) col -= colOffset;
294 int acol = sigmaProcessPtr->acol(i);
295 if (acol > 0) acol += colOffset;
296 else if (acol < 0) acol -= colOffset;
298 // Append to process record.
299 int iNow = process.append( id, status, mother1, mother2,
300 daughter1, daughter2, col, acol, phaseSpacePtr->p(i),
301 phaseSpacePtr->m(i), scale);
303 // NOAM: If there is an intermediate state, insert the it in
304 // the process record after the two incoming particles.
305 if (i == 2 && idRes != 0) {
307 // Sign of intermediate state: go by charge.
308 if (particleDataPtr->hasAnti(idRes)
309 && process[3].chargeType() + process[4].chargeType() < 0)
312 // The colour configuration of the intermediate state has to be
313 // resolved separately.
316 int m_col1 = sigmaProcessPtr->col(1);
317 int m_acol1 = sigmaProcessPtr->acol(1);
318 int m_col2 = sigmaProcessPtr->col(2);
319 int m_acol2 = sigmaProcessPtr->acol(2);
320 if (m_col1 == m_acol2 && m_col2 != m_acol1) {
323 } else if (m_col2 == m_acol1 && m_col1 != m_acol2) {
327 if ( col > 0) col += colOffset;
328 else if ( col < 0) col -= colOffset;
329 if (acol > 0) acol += colOffset;
330 else if (acol < 0) acol -= colOffset;
332 // Insert the intermediate state into the event record.
333 Vec4 pIntMed = phaseSpacePtr->p(1) + phaseSpacePtr->p(2);
334 process.append( idRes, -22, 3, 4, 6, 5 + nFin, col, acol,
335 pIntMed, pIntMed.mCalc(), scale);
338 // Pick lifetime where relevant, else not.
339 if (process[iNow].tau0() > 0.) process[iNow].tau(
340 process[iNow].tau0() * rndmPtr->exp() );
344 // Insert the outgoing particles - unresolved processes.
346 int id3 = sigmaProcessPtr->id(3);
347 int status3 = (id3 == process[1].id()) ? 14 : 15;
348 process.append( id3, status3, 1, 0, 0, 0, 0, 0,
349 phaseSpacePtr->p(3), phaseSpacePtr->m(3));
350 int id4 = sigmaProcessPtr->id(4);
351 int status4 = (id4 == process[2].id()) ? 14 : 15;
352 process.append( id4, status4, 2, 0, 0, 0, 0, 0,
353 phaseSpacePtr->p(4), phaseSpacePtr->m(4));
355 // For central diffraction: two scattered protons inserted so far,
356 // but modify mothers, add also centrally-produced hadronic system.
358 process[3].mothers( 1, 2);
359 process[4].mothers( 1, 2);
360 int id5 = sigmaProcessPtr->id(5);
362 process.append( id5, status5, 1, 2, 0, 0, 0, 0,
363 phaseSpacePtr->p(5), phaseSpacePtr->m(5));
367 // Insert the outgoing particles - Les Houches Accord processes.
370 // Since LHA partons may be out of order, determine correct one.
371 // (Recall that zeroth particle is empty.)
373 newPos.reserve(lhaUpPtr->sizePart());
375 for (int iNew = 0; iNew < lhaUpPtr->sizePart(); ++iNew) {
376 // For iNew == 0 look for the two incoming partons, then for
377 // partons having them as mothers, and so on layer by layer.
378 for (int i = 1; i < lhaUpPtr->sizePart(); ++i)
379 if (lhaUpPtr->mother1(i) == newPos[iNew]) newPos.push_back(i);
380 if (int(newPos.size()) <= iNew) break;
383 // Find scale from which to begin MPI/ISR/FSR evolution.
384 scale = lhaUpPtr->scale();
385 double scalePr = (scale < 0.) ? sqrt(Q2Fac()) : scale;
386 process.scale( scalePr);
388 // Copy over info from LHA event to process, in proper order.
389 for (int i = 1; i < lhaUpPtr->sizePart(); ++i) {
390 int iOld = newPos[i];
391 int id = lhaUpPtr->id(iOld);
393 // Translate from LHA status codes.
394 int lhaStatus = lhaUpPtr->status(iOld);
396 if (lhaStatus == 2 || lhaStatus == 3) status = -22;
397 if (lhaStatus == 1) status = 23;
399 // Find where mothers have been moved by reordering.
400 int mother1Old = lhaUpPtr->mother1(iOld);
401 int mother2Old = lhaUpPtr->mother2(iOld);
404 for (int im = 1; im < i; ++im) {
405 if (mother1Old == newPos[im]) mother1 = im + 2;
406 if (mother2Old == newPos[im]) mother2 = im + 2;
408 if (i <= 2) mother1 = i;
410 // Ensure that second mother = 0 except for bona fide carbon copies.
411 if (mother1 > 0 && mother2 == mother1) {
412 int sister1 = process[mother1].daughter1();
413 int sister2 = process[mother1].daughter2();
414 if (sister2 != sister1 && sister2 != 0) mother2 = 0;
417 // Find daughters and where they have been moved by reordering.
418 // (Values shifted two steps to account for inserted beams.)
421 for (int im = i + 1; im < lhaUpPtr->sizePart(); ++im) {
422 if (lhaUpPtr->mother1(newPos[im]) == iOld
423 || lhaUpPtr->mother2(newPos[im]) == iOld) {
424 if (daughter1 == 0 || im + 2 < daughter1) daughter1 = im + 2;
425 if (daughter2 == 0 || im + 2 > daughter2) daughter2 = im + 2;
428 // For 2 -> 1 hard scatterings reset second daughter to 0.
429 if (daughter2 == daughter1) daughter2 = 0;
431 // Colour trivial, except reset irrelevant colour indices.
432 int colType = particleDataPtr->colType(id);
433 int col1 = (colType == 1 || colType == 2 || abs(colType) == 3)
434 ? lhaUpPtr->col1(iOld) : 0;
435 int col2 = (colType == -1 || colType == 2 || abs(colType) == 3)
436 ? lhaUpPtr->col2(iOld) : 0;
439 double px = lhaUpPtr->px(iOld);
440 double py = lhaUpPtr->py(iOld);
441 double pz = lhaUpPtr->pz(iOld);
442 double e = lhaUpPtr->e(iOld);
443 double m = lhaUpPtr->m(iOld);
446 double pol = lhaUpPtr->spin(iOld);
448 // For resonance decay products use resonance mass as scale.
449 double scaleNow = scalePr;
450 if (mother1 > 4) scaleNow = process[mother1].m();
452 // Store Les Houches Accord partons.
453 int iNow = process.append( id, status, mother1, mother2, daughter1,
454 daughter2, col1, col2, Vec4(px, py, pz, e), m, scaleNow, pol);
456 // Check if need to store lifetime.
457 double tau = lhaUpPtr->tau(iOld);
458 if (tau > 0.) process[iNow].tau(tau);
462 // Loop through decay chains and set secondary vertices when needed.
463 for (int i = 3; i < process.size(); ++i) {
464 int iMother = process[i].mother1();
466 // If sister to already assigned vertex then assign same.
467 if ( process[i - 1].mother1() == iMother && process[i - 1].hasVertex() )
468 process[i].vProd( process[i - 1].vProd() );
470 // Else if mother already has vertex and/or lifetime then assign.
471 else if ( process[iMother].hasVertex() || process[iMother].tau() > 0.)
472 process[i].vProd( process[iMother].vDec() );
475 // Further info on process. Reset quantities that may or may not be known.
476 int id1Now = process[3].id();
477 int id2Now = process[4].id();
491 double x1Now, x2Now, Q2FacNow, alphaEM, alphaS, Q2Ren, sHat;
493 // Internally generated and stored information.
497 x1Now = phaseSpacePtr->x1();
498 x2Now = phaseSpacePtr->x2();
501 pdf1 = sigmaProcessPtr->pdf1();
502 pdf2 = sigmaProcessPtr->pdf2();
503 Q2FacNow = sigmaProcessPtr->Q2Fac();
504 alphaEM = sigmaProcessPtr->alphaEMRen();
505 alphaS = sigmaProcessPtr->alphaSRen();
506 Q2Ren = sigmaProcessPtr->Q2Ren();
507 sHat = phaseSpacePtr->sHat();
508 tHat = phaseSpacePtr->tHat();
509 uHat = phaseSpacePtr->uHat();
510 pTHatL = phaseSpacePtr->pTHat();
511 m3 = phaseSpacePtr->m(3);
512 m4 = phaseSpacePtr->m(4);
513 theta = phaseSpacePtr->thetaHat();
514 phi = phaseSpacePtr->phiHat();
517 // Les Houches Accord process partly available, partly to be constructed.
519 x1Now = 2. * process[3].e() / infoPtr->eCM();
520 x2Now = 2. * process[4].e() / infoPtr->eCM();
521 Q2FacNow = (scale < 0.) ? sigmaProcessPtr->Q2Fac() : pow2(scale);
522 alphaEM = lhaUpPtr->alphaQED();
523 if (alphaEM < 0.001) alphaEM = sigmaProcessPtr->alphaEMRen();
524 alphaS = lhaUpPtr->alphaQCD();
525 if (alphaS < 0.001) alphaS = sigmaProcessPtr->alphaSRen();
526 Q2Ren = (scale < 0.) ? sigmaProcessPtr->Q2Ren() : pow2(scale);
527 Vec4 pSum = process[3].p() + process[4].p();
531 for (int i = 5; i < process.size(); ++i)
532 if (process[i].mother1() == 3 && process[i].mother2() == 4) {
534 pTLH += process[i].pT();
536 if (nFinLH > 0) pTHatL = pTLH / nFinLH;
538 // Read info on parton densities if provided.
539 id1pdf = lhaUpPtr->id1pdf();
540 id2pdf = lhaUpPtr->id2pdf();
541 x1pdf = lhaUpPtr->x1pdf();
542 x2pdf = lhaUpPtr->x2pdf();
543 if (lhaUpPtr->pdfIsSet()) {
544 pdf1 = lhaUpPtr->pdf1();
545 pdf2 = lhaUpPtr->pdf2();
546 Q2FacNow = pow2(lhaUpPtr->scalePDF());
549 // Reconstruct kinematics of 2 -> 2 processes from momenta.
551 Vec4 pDifT = process[3].p() - process[5].p();
552 tHat = pDifT * pDifT;
553 Vec4 pDifU = process[3].p() - process[6].p();
554 uHat = pDifU * pDifU;
555 pTHatL = process[5].pT();
558 Vec4 p5 = process[5].p();
561 phi = process[5].phi();
565 // Store information.
567 infoPtr->setPDFalpha( 0, id1pdf, id2pdf, x1pdf, x2pdf, pdf1, pdf2,
568 Q2FacNow, alphaEM, alphaS, Q2Ren);
569 infoPtr->setKin( 0, id1Now, id2Now, x1Now, x2Now, sHat, tHat, uHat,
570 pTHatL, m3, m4, theta, phi);
572 infoPtr->setTypeMPI( code(), pTHatL);
574 // For Les Houches event store subprocess classification.
576 int codeSub = lhaUpPtr->idProcess();
577 ostringstream nameSub;
578 nameSub << "user process " << codeSub;
579 infoPtr->setSubType( 0, nameSub.str(), codeSub, nFin);
587 //--------------------------------------------------------------------------
589 // Give the hard resonance decay chain from Les Houches input.
590 // Note: mildly modified copy of some constructProcess code; to streamline.
592 bool ProcessContainer::constructDecays( Event& process) {
594 // Let hard process record begin with the event as a whole.
596 process.append( 90, -11, 0, 0, 0, 0, 0, 0, Vec4(0., 0., 0., 0.), 0., 0. );
598 // Since LHA partons may be out of order, determine correct one.
599 // (Recall that zeroth particle is empty.)
601 newPos.reserve(lhaUpPtr->sizePart());
603 for (int iNew = 0; iNew < lhaUpPtr->sizePart(); ++iNew) {
604 // For iNew == 0 look for the two incoming partons, then for
605 // partons having them as mothers, and so on layer by layer.
606 for (int i = 1; i < lhaUpPtr->sizePart(); ++i)
607 if (lhaUpPtr->mother1(i) == newPos[iNew]) newPos.push_back(i);
608 if (int(newPos.size()) <= iNew) break;
611 // Find scale from which to begin MPI/ISR/FSR evolution.
612 double scale = lhaUpPtr->scale();
613 process.scale( scale);
616 // Copy over info from LHA event to process, in proper order.
617 for (int i = 1; i < lhaUpPtr->sizePart(); ++i) {
618 int iOld = newPos[i];
619 int id = lhaUpPtr->id(iOld);
621 // Translate from LHA status codes.
622 int lhaStatus = lhaUpPtr->status(iOld);
624 if (lhaStatus == 2 || lhaStatus == 3) status = -22;
625 if (lhaStatus == 1) status = 23;
627 // Find where mothers have been moved by reordering.
628 int mother1Old = lhaUpPtr->mother1(iOld);
629 int mother2Old = lhaUpPtr->mother2(iOld);
632 for (int im = 1; im < i; ++im) {
633 if (mother1Old == newPos[im]) mother1 = im;
634 if (mother2Old == newPos[im]) mother2 = im;
637 // Ensure that second mother = 0 except for bona fide carbon copies.
638 if (mother1 > 0 && mother2 == mother1) {
639 int sister1 = process[mother1].daughter1();
640 int sister2 = process[mother1].daughter2();
641 if (sister2 != sister1 && sister2 != 0) mother2 = 0;
644 // Find daughters and where they have been moved by reordering.
647 for (int im = i + 1; im < lhaUpPtr->sizePart(); ++im) {
648 if (lhaUpPtr->mother1(newPos[im]) == iOld
649 || lhaUpPtr->mother2(newPos[im]) == iOld) {
650 if (daughter1 == 0 || im < daughter1) daughter1 = im;
651 if (daughter2 == 0 || im > daughter2) daughter2 = im;
654 // For 2 -> 1 hard scatterings reset second daughter to 0.
655 if (daughter2 == daughter1) daughter2 = 0;
657 // Colour trivial, except reset irrelevant colour indices.
658 int colType = particleDataPtr->colType(id);
659 int col1 = (colType == 1 || colType == 2 || abs(colType) == 3)
660 ? lhaUpPtr->col1(iOld) : 0;
661 int col2 = (colType == -1 || colType == 2 || abs(colType) == 3)
662 ? lhaUpPtr->col2(iOld) : 0;
665 double px = lhaUpPtr->px(iOld);
666 double py = lhaUpPtr->py(iOld);
667 double pz = lhaUpPtr->pz(iOld);
668 double e = lhaUpPtr->e(iOld);
669 double m = lhaUpPtr->m(iOld);
670 if (status > 0) pSum += Vec4( px, py, pz, e);
673 double pol = lhaUpPtr->spin(iOld);
675 // For resonance decay products use resonance mass as scale.
676 double scaleNow = scale;
677 if (mother1 > 0) scaleNow = process[mother1].m();
679 // Store Les Houches Accord partons.
680 int iNow = process.append( id, status, mother1, mother2, daughter1,
681 daughter2, col1, col2, Vec4(px, py, pz, e), m, scaleNow, pol);
683 // Check if need to store lifetime.
684 double tau = lhaUpPtr->tau(iOld);
685 if (tau > 0.) process[iNow].tau(tau);
688 // Update four-momentum of system as a whole.
690 process[0].m( pSum.mCalc());
692 // Loop through decay chains and set secondary vertices when needed.
693 for (int i = 1; i < process.size(); ++i) {
694 int iMother = process[i].mother1();
696 // If sister to already assigned vertex then assign same.
697 if ( process[i - 1].mother1() == iMother && process[i - 1].hasVertex()
698 && i > 1) process[i].vProd( process[i - 1].vProd() );
700 // Else if mother already has vertex and/or lifetime then assign.
701 else if ( process[iMother].hasVertex() || process[iMother].tau() > 0.)
702 process[i].vProd( process[iMother].vDec() );
710 //--------------------------------------------------------------------------
712 // Handle resonance decays.
714 bool ProcessContainer::decayResonances( Event& process) {
716 // Save current event-record size and status codes.
718 vector<int> statusSave( process.size());
719 for (int i = 0; i < process.size(); ++i)
720 statusSave[i] = process[i].status();
721 bool physical = true;
722 bool newChain = false;
723 bool newFlavours = false;
725 // Do loop over user veto.
728 // Do sequential chain of uncorrelated isotropic decays.
730 physical = resDecaysPtr->next( process);
731 if (!physical) return false;
733 // Check whether flavours should be correlated.
734 // (Currently only relevant for f fbar -> gamma*/Z0 gamma*/Z0.)
735 newFlavours = ( sigmaProcessPtr->weightDecayFlav( process)
738 // Reset the decay chains if have to redo.
740 process.restoreSize();
741 for (int i = 0; i < process.size(); ++i)
742 process[i].status( statusSave[i]);
745 // Loop back where required to generate new decays with new flavours.
746 } while (newFlavours);
748 // Correct to nonisotropic decays.
749 phaseSpacePtr->decayKinematics( process);
751 // Optionally user hooks check/veto on decay chain.
753 newChain = userHooksPtr->doVetoResonanceDecays( process);
755 // Reset the decay chains if have to redo.
757 process.restoreSize();
758 for (int i = 0; i < process.size(); ++i)
759 process[i].status( statusSave[i]);
762 // Loop back where required to generate new decay chain.
770 //--------------------------------------------------------------------------
772 // Reset event generation statistics; but NOT maximum of cross section.
774 void ProcessContainer::reset() {
790 //--------------------------------------------------------------------------
792 // Estimate integrated cross section and its uncertainty.
794 void ProcessContainer::sigmaDelta() {
796 // Initial values. No analysis meaningful unless accepted events.
801 if (nAcc == 0) return;
803 // Average value. No error analysis unless at least two events.
804 double nTryInv = 1. / nTry;
805 double nSelInv = 1. / nSel;
806 double nAccInv = 1. / nAcc;
807 sigmaAvg = sigmaSum * nTryInv;
808 double fracAcc = nAcc * nSelInv;
809 sigmaFin = sigmaAvg * fracAcc;
811 if (nAcc == 1) return;
813 // Estimated error. Quadratic sum of cross section term and
814 // binomial from accept/reject step.
815 double delta2Sig = (lhaStratAbs != 3)
816 ? (sigma2Sum * nTryInv - pow2(sigmaAvg)) * nTryInv / pow2(sigmaAvg)
817 : pow2( lhaUpPtr->xErrSum() / lhaUpPtr->xSecSum());
818 double delta2Veto = (nSel - nAcc) * nAccInv * nSelInv;
819 double delta2Sum = delta2Sig + delta2Veto;
820 deltaFin = sqrtpos(delta2Sum) * sigmaFin;
824 //==========================================================================
826 // SetupContainer class.
827 // Turns list of user-desired processes into a vector of containers.
829 //--------------------------------------------------------------------------
831 // Main routine to initialize list of processes.
833 bool SetupContainers::init(vector<ProcessContainer*>& containerPtrs,
834 Settings& settings, ParticleData* particleDataPtr, Couplings* couplings) {
836 // Reset process list, if filled in previous subrun.
837 if (containerPtrs.size() > 0) {
838 for (int i = 0; i < int(containerPtrs.size()); ++i)
839 delete containerPtrs[i];
840 containerPtrs.clear();
842 SigmaProcess* sigmaPtr;
844 // Set up requested objects for soft QCD processes.
845 bool softQCD = settings.flag("SoftQCD:all");
846 if (softQCD || settings.flag("SoftQCD:minBias")) {
847 sigmaPtr = new Sigma0minBias;
848 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
850 if (softQCD || settings.flag("SoftQCD:elastic")) {
851 sigmaPtr = new Sigma0AB2AB;
852 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
854 if (softQCD || settings.flag("SoftQCD:singleDiffractive")) {
855 sigmaPtr = new Sigma0AB2XB;
856 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
857 sigmaPtr = new Sigma0AB2AX;
858 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
860 if (softQCD || settings.flag("SoftQCD:doubleDiffractive")) {
861 sigmaPtr = new Sigma0AB2XX;
862 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
864 if (softQCD || settings.flag("SoftQCD:centralDiffractive")) {
865 sigmaPtr = new Sigma0AB2AXB;
866 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
869 // Set up requested objects for hard QCD processes.
870 bool hardQCD = settings.flag("HardQCD:all");
871 if (hardQCD || settings.flag("HardQCD:gg2gg")) {
872 sigmaPtr = new Sigma2gg2gg;
873 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
875 if (hardQCD || settings.flag("HardQCD:gg2qqbar")) {
876 sigmaPtr = new Sigma2gg2qqbar;
877 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
879 if (hardQCD || settings.flag("HardQCD:qg2qg")) {
880 sigmaPtr = new Sigma2qg2qg;
881 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
883 if (hardQCD || settings.flag("HardQCD:qq2qq")) {
884 sigmaPtr = new Sigma2qq2qq;
885 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
887 if (hardQCD || settings.flag("HardQCD:qqbar2gg")) {
888 sigmaPtr = new Sigma2qqbar2gg;
889 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
891 if (hardQCD || settings.flag("HardQCD:qqbar2qqbarNew")) {
892 sigmaPtr = new Sigma2qqbar2qqbarNew;
893 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
896 // Set up requested objects for c cbar and b bbar, also hard QCD.
897 if (hardQCD || settings.flag("HardQCD:gg2ccbar")) {
898 sigmaPtr = new Sigma2gg2QQbar(4, 121);
899 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
901 if (hardQCD || settings.flag("HardQCD:qqbar2ccbar")) {
902 sigmaPtr = new Sigma2qqbar2QQbar(4, 122);
903 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
905 if (hardQCD || settings.flag("HardQCD:gg2bbbar")) {
906 sigmaPtr = new Sigma2gg2QQbar(5, 123);
907 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
909 if (hardQCD || settings.flag("HardQCD:qqbar2bbbar")) {
910 sigmaPtr = new Sigma2qqbar2QQbar(5, 124);
911 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
914 // Set up requested objects for hard QCD 2 -> 3 processes.
915 bool hardQCD3parton = settings.flag("HardQCD:3parton");
916 if (hardQCD3parton || settings.flag("HardQCD:gg2ggg")) {
917 sigmaPtr = new Sigma3gg2ggg;
918 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
920 if (hardQCD3parton || settings.flag("HardQCD:qqbar2ggg")) {
921 sigmaPtr = new Sigma3qqbar2ggg;
922 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
924 if (hardQCD3parton || settings.flag("HardQCD:qg2qgg")) {
925 sigmaPtr = new Sigma3qg2qgg;
926 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
928 if (hardQCD3parton || settings.flag("HardQCD:qq2qqgDiff")) {
929 sigmaPtr = new Sigma3qq2qqgDiff;
930 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
932 if (hardQCD3parton || settings.flag("HardQCD:qq2qqgSame")) {
933 sigmaPtr = new Sigma3qq2qqgSame;
934 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
936 if (hardQCD3parton || settings.flag("HardQCD:qqbar2qqbargDiff")) {
937 sigmaPtr = new Sigma3qqbar2qqbargDiff;
938 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
940 if (hardQCD3parton || settings.flag("HardQCD:qqbar2qqbargSame")) {
941 sigmaPtr = new Sigma3qqbar2qqbargSame;
942 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
944 if (hardQCD3parton || settings.flag("HardQCD:gg2qqbarg")) {
945 sigmaPtr = new Sigma3gg2qqbarg;
946 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
948 if (hardQCD3parton || settings.flag("HardQCD:qg2qqqbarDiff")) {
949 sigmaPtr = new Sigma3qg2qqqbarDiff;
950 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
952 if (hardQCD3parton || settings.flag("HardQCD:qg2qqqbarSame")) {
953 sigmaPtr = new Sigma3qg2qqqbarSame;
954 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
957 // Set up requested objects for prompt photon processes.
958 bool promptPhotons = settings.flag("PromptPhoton:all");
960 || settings.flag("PromptPhoton:qg2qgamma")) {
961 sigmaPtr = new Sigma2qg2qgamma;
962 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
965 || settings.flag("PromptPhoton:qqbar2ggamma")) {
966 sigmaPtr = new Sigma2qqbar2ggamma;
967 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
970 || settings.flag("PromptPhoton:gg2ggamma")) {
971 sigmaPtr = new Sigma2gg2ggamma;
972 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
975 || settings.flag("PromptPhoton:ffbar2gammagamma")) {
976 sigmaPtr = new Sigma2ffbar2gammagamma;
977 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
980 || settings.flag("PromptPhoton:gg2gammagamma")) {
981 sigmaPtr = new Sigma2gg2gammagamma;
982 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
985 // Set up requested objects for weak gauge boson t-channel exchange.
986 bool weakBosonExchanges = settings.flag("WeakBosonExchange:all");
987 if (weakBosonExchanges
988 || settings.flag("WeakBosonExchange:ff2ff(t:gmZ)")) {
989 sigmaPtr = new Sigma2ff2fftgmZ;
990 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
992 if (weakBosonExchanges
993 || settings.flag("WeakBosonExchange:ff2ff(t:W)")) {
994 sigmaPtr = new Sigma2ff2fftW;
995 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
998 // Set up requested objects for weak gauge boson processes.
999 bool weakSingleBosons = settings.flag("WeakSingleBoson:all");
1000 if (weakSingleBosons
1001 || settings.flag("WeakSingleBoson:ffbar2gmZ")) {
1002 sigmaPtr = new Sigma1ffbar2gmZ;
1003 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1005 if (weakSingleBosons
1006 || settings.flag("WeakSingleBoson:ffbar2W")) {
1007 sigmaPtr = new Sigma1ffbar2W;
1008 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1011 // Set up requested object for s-channel gamma exchange.
1012 // Subset of gamma*/Z0 above, intended for multiparton interactions.
1013 if (settings.flag("WeakSingleBoson:ffbar2ffbar(s:gm)")) {
1014 sigmaPtr = new Sigma2ffbar2ffbarsgm;
1015 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1018 // Set up requested objects for weak gauge boson pair processes.
1019 bool weakDoubleBosons = settings.flag("WeakDoubleBoson:all");
1020 if (weakDoubleBosons
1021 || settings.flag("WeakDoubleBoson:ffbar2gmZgmZ")) {
1022 sigmaPtr = new Sigma2ffbar2gmZgmZ;
1023 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1025 if (weakDoubleBosons
1026 || settings.flag("WeakDoubleBoson:ffbar2ZW")) {
1027 sigmaPtr = new Sigma2ffbar2ZW;
1028 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1030 if (weakDoubleBosons
1031 || settings.flag("WeakDoubleBoson:ffbar2WW")) {
1032 sigmaPtr = new Sigma2ffbar2WW;
1033 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1036 // Set up requested objects for weak gauge boson + parton processes.
1037 bool weakBosonAndPartons = settings.flag("WeakBosonAndParton:all");
1038 if (weakBosonAndPartons
1039 || settings.flag("WeakBosonAndParton:qqbar2gmZg")) {
1040 sigmaPtr = new Sigma2qqbar2gmZg;
1041 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1043 if (weakBosonAndPartons
1044 || settings.flag("WeakBosonAndParton:qg2gmZq")) {
1045 sigmaPtr = new Sigma2qg2gmZq;
1046 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1048 if (weakBosonAndPartons
1049 || settings.flag("WeakBosonAndParton:ffbar2gmZgm")) {
1050 sigmaPtr = new Sigma2ffbar2gmZgm;
1051 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1053 if (weakBosonAndPartons
1054 || settings.flag("WeakBosonAndParton:fgm2gmZf")) {
1055 sigmaPtr = new Sigma2fgm2gmZf;
1056 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1058 if (weakBosonAndPartons
1059 || settings.flag("WeakBosonAndParton:qqbar2Wg")) {
1060 sigmaPtr = new Sigma2qqbar2Wg;
1061 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1063 if (weakBosonAndPartons
1064 || settings.flag("WeakBosonAndParton:qg2Wq")) {
1065 sigmaPtr = new Sigma2qg2Wq;
1066 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1068 if (weakBosonAndPartons
1069 || settings.flag("WeakBosonAndParton:ffbar2Wgm")) {
1070 sigmaPtr = new Sigma2ffbar2Wgm;
1071 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1073 if (weakBosonAndPartons
1074 || settings.flag("WeakBosonAndParton:fgm2Wf")) {
1075 sigmaPtr = new Sigma2fgm2Wf;
1076 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1079 // Set up requested objects for photon collision processes.
1080 bool photonCollisions = settings.flag("PhotonCollision:all");
1081 if (photonCollisions || settings.flag("PhotonCollision:gmgm2qqbar")) {
1082 sigmaPtr = new Sigma2gmgm2ffbar(1, 261);
1083 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1085 if (photonCollisions || settings.flag("PhotonCollision:gmgm2ccbar")) {
1086 sigmaPtr = new Sigma2gmgm2ffbar(4, 262);
1087 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1089 if (photonCollisions || settings.flag("PhotonCollision:gmgm2bbbar")) {
1090 sigmaPtr = new Sigma2gmgm2ffbar(5, 263);
1091 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1093 if (photonCollisions || settings.flag("PhotonCollision:gmgm2ee")) {
1094 sigmaPtr = new Sigma2gmgm2ffbar(11, 264);
1095 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1097 if (photonCollisions || settings.flag("PhotonCollision:gmgm2mumu")) {
1098 sigmaPtr = new Sigma2gmgm2ffbar(13, 265);
1099 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1101 if (photonCollisions || settings.flag("PhotonCollision:gmgm2tautau")) {
1102 sigmaPtr = new Sigma2gmgm2ffbar(15, 266);
1103 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1106 // Set up requested objects for charmonium production
1107 bool charmoniums = settings.flag("Charmonium:all");
1108 if (charmoniums || settings.flag("Charmonium:gg2QQbar[3S1(1)]g")) {
1109 sigmaPtr = new Sigma2gg2QQbar3S11g(4, 401);
1110 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1112 if (charmoniums || settings.flag("Charmonium:gg2QQbar[3P0(1)]g")) {
1113 sigmaPtr = new Sigma2gg2QQbar3PJ1g(4, 0, 402);
1114 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1116 if (charmoniums || settings.flag("Charmonium:gg2QQbar[3P1(1)]g")) {
1117 sigmaPtr = new Sigma2gg2QQbar3PJ1g(4, 1, 403);
1118 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1120 if (charmoniums || settings.flag("Charmonium:gg2QQbar[3P2(1)]g")) {
1121 sigmaPtr = new Sigma2gg2QQbar3PJ1g(4, 2, 404);
1122 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1124 if (charmoniums || settings.flag("Charmonium:qg2QQbar[3P0(1)]q")) {
1125 sigmaPtr = new Sigma2qg2QQbar3PJ1q(4, 0, 405);
1126 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1128 if (charmoniums || settings.flag("Charmonium:qg2QQbar[3P1(1)]q")) {
1129 sigmaPtr = new Sigma2qg2QQbar3PJ1q(4, 1, 406);
1130 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1132 if (charmoniums || settings.flag("Charmonium:qg2QQbar[3P2(1)]q")) {
1133 sigmaPtr = new Sigma2qg2QQbar3PJ1q(4, 2, 407);
1134 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1136 if (charmoniums || settings.flag("Charmonium:qqbar2QQbar[3P0(1)]g")) {
1137 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(4, 0, 408);
1138 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1140 if (charmoniums || settings.flag("Charmonium:qqbar2QQbar[3P1(1)]g")) {
1141 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(4, 1, 409);
1142 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1144 if (charmoniums || settings.flag("Charmonium:qqbar2QQbar[3P2(1)]g")) {
1145 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(4, 2, 410);
1146 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1148 if (charmoniums || settings.flag("Charmonium:gg2QQbar[3S1(8)]g")) {
1149 sigmaPtr = new Sigma2gg2QQbarX8g(4, 0, 411);
1150 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1152 if (charmoniums || settings.flag("Charmonium:gg2QQbar[1S0(8)]g")) {
1153 sigmaPtr = new Sigma2gg2QQbarX8g(4, 1, 412);
1154 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1156 if (charmoniums || settings.flag("Charmonium:gg2QQbar[3PJ(8)]g")) {
1157 sigmaPtr = new Sigma2gg2QQbarX8g(4, 2, 413);
1158 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1160 if (charmoniums || settings.flag("Charmonium:qg2QQbar[3S1(8)]q")) {
1161 sigmaPtr = new Sigma2qg2QQbarX8q(4, 0, 414);
1162 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1164 if (charmoniums || settings.flag("Charmonium:qg2QQbar[1S0(8)]q")) {
1165 sigmaPtr = new Sigma2qg2QQbarX8q(4, 1, 415);
1166 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1168 if (charmoniums || settings.flag("Charmonium:qg2QQbar[3PJ(8)]q")) {
1169 sigmaPtr = new Sigma2qg2QQbarX8q(4, 2, 416);
1170 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1172 if (charmoniums || settings.flag("Charmonium:qqbar2QQbar[3S1(8)]g")) {
1173 sigmaPtr = new Sigma2qqbar2QQbarX8g(4, 0, 417);
1174 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1176 if (charmoniums || settings.flag("Charmonium:qqbar2QQbar[1S0(8)]g")) {
1177 sigmaPtr = new Sigma2qqbar2QQbarX8g(4, 1, 418);
1178 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1180 if (charmoniums || settings.flag("Charmonium:qqbar2QQbar[3PJ(8)]g")) {
1181 sigmaPtr = new Sigma2qqbar2QQbarX8g(4, 2, 419);
1182 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1185 // Set up requested objects for bottomonium production
1186 bool bottomoniums = settings.flag("Bottomonium:all");
1187 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[3S1(1)]g")) {
1188 sigmaPtr = new Sigma2gg2QQbar3S11g(5, 501);
1189 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1191 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[3P0(1)]g")) {
1192 sigmaPtr = new Sigma2gg2QQbar3PJ1g(5, 0, 502);
1193 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1195 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[3P1(1)]g")) {
1196 sigmaPtr = new Sigma2gg2QQbar3PJ1g(5, 1, 503);
1197 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1199 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[3P2(1)]g")) {
1200 sigmaPtr = new Sigma2gg2QQbar3PJ1g(5, 2, 504);
1201 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1203 if (bottomoniums || settings.flag("Bottomonium:qg2QQbar[3P0(1)]q")) {
1204 sigmaPtr = new Sigma2qg2QQbar3PJ1q(5, 0, 505);
1205 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1207 if (bottomoniums || settings.flag("Bottomonium:qg2QQbar[3P1(1)]q")) {
1208 sigmaPtr = new Sigma2qg2QQbar3PJ1q(5, 1, 506);
1209 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1211 if (bottomoniums || settings.flag("Bottomonium:qg2QQbar[3P2(1)]q")) {
1212 sigmaPtr = new Sigma2qg2QQbar3PJ1q(5, 2, 507);
1213 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1215 if (bottomoniums || settings.flag("Bottomonium:qqbar2QQbar[3P0(1)]g")) {
1216 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(5, 0, 508);
1217 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1219 if (bottomoniums || settings.flag("Bottomonium:qqbar2QQbar[3P1(1)]g")) {
1220 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(5, 1, 509);
1221 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1223 if (bottomoniums || settings.flag("Bottomonium:qqbar2QQbar[3P2(1)]g")) {
1224 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(5, 2, 510);
1225 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1227 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[3S1(8)]g")) {
1228 sigmaPtr = new Sigma2gg2QQbarX8g(5, 0, 511);
1229 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1231 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[1S0(8)]g")) {
1232 sigmaPtr = new Sigma2gg2QQbarX8g(5, 1, 512);
1233 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1235 if (bottomoniums || settings.flag("Bottomonium:gg2QQbar[3PJ(8)]g")) {
1236 sigmaPtr = new Sigma2gg2QQbarX8g(5, 2, 513);
1237 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1239 if (bottomoniums || settings.flag("Bottomonium:qg2QQbar[3S1(8)]q")) {
1240 sigmaPtr = new Sigma2qg2QQbarX8q(5, 0, 514);
1241 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1243 if (bottomoniums || settings.flag("Bottomonium:qg2QQbar[1S0(8)]q")) {
1244 sigmaPtr = new Sigma2qg2QQbarX8q(5, 1, 515);
1245 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1247 if (bottomoniums || settings.flag("Bottomonium:qg2QQbar[3PJ(8)]q")) {
1248 sigmaPtr = new Sigma2qg2QQbarX8q(5, 2, 516);
1249 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1251 if (bottomoniums || settings.flag("Bottomonium:qqbar2QQbar[3S1(8)]g")) {
1252 sigmaPtr = new Sigma2qqbar2QQbarX8g(5, 0, 517);
1253 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1255 if (bottomoniums || settings.flag("Bottomonium:qqbar2QQbar[1S0(8)]g")) {
1256 sigmaPtr = new Sigma2qqbar2QQbarX8g(5, 1, 518);
1257 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1259 if (bottomoniums || settings.flag("Bottomonium:qqbar2QQbar[3PJ(8)]g")) {
1260 sigmaPtr = new Sigma2qqbar2QQbarX8g(5, 2, 519);
1261 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1264 // Set up requested objects for top production
1265 bool tops = settings.flag("Top:all");
1266 if (tops || settings.flag("Top:gg2ttbar")) {
1267 sigmaPtr = new Sigma2gg2QQbar(6, 601);
1268 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1270 if (tops || settings.flag("Top:qqbar2ttbar")) {
1271 sigmaPtr = new Sigma2qqbar2QQbar(6, 602);
1272 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1274 if (tops || settings.flag("Top:qq2tq(t:W)")) {
1275 sigmaPtr = new Sigma2qq2QqtW(6, 603);
1276 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1278 if (tops || settings.flag("Top:ffbar2ttbar(s:gmZ)")) {
1279 sigmaPtr = new Sigma2ffbar2FFbarsgmZ(6, 604);
1280 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1282 if (tops || settings.flag("Top:ffbar2tqbar(s:W)")) {
1283 sigmaPtr = new Sigma2ffbar2FfbarsW(6, 0, 605);
1284 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1286 if (tops || settings.flag("Top:gmgm2ttbar")) {
1287 sigmaPtr = new Sigma2gmgm2ffbar(6, 606);
1288 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1291 // Set up requested objects for fourth-generation b' production
1292 bool bPrimes = settings.flag("FourthBottom:all");
1293 if (bPrimes || settings.flag("FourthBottom:gg2bPrimebPrimebar")) {
1294 sigmaPtr = new Sigma2gg2QQbar(7, 801);
1295 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1297 if (bPrimes || settings.flag("FourthBottom:qqbar2bPrimebPrimebar")) {
1298 sigmaPtr = new Sigma2qqbar2QQbar(7, 802);
1299 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1301 if (bPrimes || settings.flag("FourthBottom:qq2bPrimeq(t:W)")) {
1302 sigmaPtr = new Sigma2qq2QqtW(7, 803);
1303 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1305 if (bPrimes || settings.flag("FourthBottom:ffbar2bPrimebPrimebar(s:gmZ)")) {
1306 sigmaPtr = new Sigma2ffbar2FFbarsgmZ(7, 804);
1307 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1309 if (bPrimes || settings.flag("FourthBottom:ffbar2bPrimeqbar(s:W)")) {
1310 sigmaPtr = new Sigma2ffbar2FfbarsW(7, 0, 805);
1311 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1313 if (bPrimes || settings.flag("FourthBottom:ffbar2bPrimetbar(s:W)")) {
1314 sigmaPtr = new Sigma2ffbar2FfbarsW(7, 6, 806);
1315 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1318 // Set up requested objects for fourth-generation t' production
1319 bool tPrimes = settings.flag("FourthTop:all");
1320 if (tPrimes || settings.flag("FourthTop:gg2tPrimetPrimebar")) {
1321 sigmaPtr = new Sigma2gg2QQbar(8, 821);
1322 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1324 if (tPrimes || settings.flag("FourthTop:qqbar2tPrimetPrimebar")) {
1325 sigmaPtr = new Sigma2qqbar2QQbar(8, 822);
1326 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1328 if (tPrimes || settings.flag("FourthTop:qq2tPrimeq(t:W)")) {
1329 sigmaPtr = new Sigma2qq2QqtW(8, 823);
1330 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1332 if (tPrimes || settings.flag("FourthTop:ffbar2tPrimetPrimebar(s:gmZ)")) {
1333 sigmaPtr = new Sigma2ffbar2FFbarsgmZ(8, 824);
1334 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1336 if (tPrimes || settings.flag("FourthTop:ffbar2tPrimeqbar(s:W)")) {
1337 sigmaPtr = new Sigma2ffbar2FfbarsW(8, 0, 825);
1338 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1341 // Set up requested objects for two different fourth-generation fermions.
1342 if (bPrimes || tPrimes
1343 || settings.flag("FourthPair:ffbar2tPrimebPrimebar(s:W)")) {
1344 sigmaPtr = new Sigma2ffbar2FfbarsW(8, 7, 841);
1345 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1347 if (settings.flag("FourthPair:ffbar2tauPrimenuPrimebar(s:W)")) {
1348 sigmaPtr = new Sigma2ffbar2FfbarsW(17, 18, 842);
1349 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1352 // Flag for global choice between SM and BSM Higgses.
1353 bool useBSMHiggses = settings.flag("Higgs:useBSM");
1355 // Set up requested objects for Standard-Model Higgs production.
1356 if (!useBSMHiggses) {
1357 bool HiggsesSM = settings.flag("HiggsSM:all");
1358 if (HiggsesSM || settings.flag("HiggsSM:ffbar2H")) {
1359 sigmaPtr = new Sigma1ffbar2H(0);
1360 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1362 if (HiggsesSM || settings.flag("HiggsSM:gg2H")) {
1363 sigmaPtr = new Sigma1gg2H(0);
1364 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1366 if (HiggsesSM || settings.flag("HiggsSM:gmgm2H")) {
1367 sigmaPtr = new Sigma1gmgm2H(0);
1368 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1370 if (HiggsesSM || settings.flag("HiggsSM:ffbar2HZ")) {
1371 sigmaPtr = new Sigma2ffbar2HZ(0);
1372 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1374 if (HiggsesSM || settings.flag("HiggsSM:ffbar2HW")) {
1375 sigmaPtr = new Sigma2ffbar2HW(0);
1376 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1378 if (HiggsesSM || settings.flag("HiggsSM:ff2Hff(t:ZZ)")) {
1379 sigmaPtr = new Sigma3ff2HfftZZ(0);
1380 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1382 if (HiggsesSM || settings.flag("HiggsSM:ff2Hff(t:WW)")) {
1383 sigmaPtr = new Sigma3ff2HfftWW(0);
1384 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1386 if (HiggsesSM || settings.flag("HiggsSM:gg2Httbar")) {
1387 sigmaPtr = new Sigma3gg2HQQbar(6,0);
1388 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1390 if (HiggsesSM || settings.flag("HiggsSM:qqbar2Httbar")) {
1391 sigmaPtr = new Sigma3qqbar2HQQbar(6,0);
1392 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1395 // Further Standard-Model Higgs processes, not included in "all".
1396 if (settings.flag("HiggsSM:qg2Hq")) {
1397 sigmaPtr = new Sigma2qg2Hq(4,0);
1398 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1399 sigmaPtr = new Sigma2qg2Hq(5,0);
1400 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1402 if (settings.flag("HiggsSM:gg2Hbbbar")) {
1403 sigmaPtr = new Sigma3gg2HQQbar(5,0);
1404 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1406 if (settings.flag("HiggsSM:qqbar2Hbbbar")) {
1407 sigmaPtr = new Sigma3qqbar2HQQbar(5,0);
1408 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1410 if (settings.flag("HiggsSM:gg2Hg(l:t)")) {
1411 sigmaPtr = new Sigma2gg2Hglt(0);
1412 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1414 if (settings.flag("HiggsSM:qg2Hq(l:t)")) {
1415 sigmaPtr = new Sigma2qg2Hqlt(0);
1416 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1418 if (settings.flag("HiggsSM:qqbar2Hg(l:t)")) {
1419 sigmaPtr = new Sigma2qqbar2Hglt(0);
1420 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1424 // Common switch for the group of Higgs production BSM.
1425 if (useBSMHiggses) {
1426 bool HiggsesBSM = settings.flag("HiggsBSM:all");
1428 // Set up requested objects for BSM H1 production.
1429 bool HiggsesH1 = settings.flag("HiggsBSM:allH1");
1430 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:ffbar2H1")) {
1431 sigmaPtr = new Sigma1ffbar2H(1);
1432 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1434 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:gg2H1")) {
1435 sigmaPtr = new Sigma1gg2H(1);
1436 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1438 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:gmgm2H1")) {
1439 sigmaPtr = new Sigma1gmgm2H(1);
1440 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1442 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:ffbar2H1Z")) {
1443 sigmaPtr = new Sigma2ffbar2HZ(1);
1444 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1446 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:ffbar2H1W")) {
1447 sigmaPtr = new Sigma2ffbar2HW(1);
1448 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1450 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:ff2H1ff(t:ZZ)")) {
1451 sigmaPtr = new Sigma3ff2HfftZZ(1);
1452 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1454 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:ff2H1ff(t:WW)")) {
1455 sigmaPtr = new Sigma3ff2HfftWW(1);
1456 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1458 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:gg2H1ttbar")) {
1459 sigmaPtr = new Sigma3gg2HQQbar(6,1);
1460 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1462 if (HiggsesBSM || HiggsesH1 || settings.flag("HiggsBSM:qqbar2H1ttbar")) {
1463 sigmaPtr = new Sigma3qqbar2HQQbar(6,1);
1464 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1467 // Further BSM H1 processes, not included in "all".
1468 if (settings.flag("HiggsBSM:qg2H1q")) {
1469 sigmaPtr = new Sigma2qg2Hq(4,1);
1470 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1471 sigmaPtr = new Sigma2qg2Hq(5,1);
1472 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1474 if (settings.flag("HiggsBSM:gg2H1bbbar")) {
1475 sigmaPtr = new Sigma3gg2HQQbar(5,1);
1476 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1478 if (settings.flag("HiggsBSM:qqbar2H1bbbar")) {
1479 sigmaPtr = new Sigma3qqbar2HQQbar(5,1);
1480 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1482 if (settings.flag("HiggsBSM:gg2H1g(l:t)")) {
1483 sigmaPtr = new Sigma2gg2Hglt(1);
1484 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1486 if (settings.flag("HiggsBSM:qg2H1q(l:t)")) {
1487 sigmaPtr = new Sigma2qg2Hqlt(1);
1488 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1490 if (settings.flag("HiggsBSM:qqbar2H1g(l:t)")) {
1491 sigmaPtr = new Sigma2qqbar2Hglt(1);
1492 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1495 // Set up requested objects for BSM H2 production.
1496 bool HiggsesH2 = settings.flag("HiggsBSM:allH2");
1497 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:ffbar2H2")) {
1498 sigmaPtr = new Sigma1ffbar2H(2);
1499 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1501 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:gg2H2")) {
1502 sigmaPtr = new Sigma1gg2H(2);
1503 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1505 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:gmgm2H2")) {
1506 sigmaPtr = new Sigma1gmgm2H(2);
1507 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1509 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:ffbar2H2Z")) {
1510 sigmaPtr = new Sigma2ffbar2HZ(2);
1511 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1513 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:ffbar2H2W")) {
1514 sigmaPtr = new Sigma2ffbar2HW(2);
1515 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1517 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:ff2H2ff(t:ZZ)")) {
1518 sigmaPtr = new Sigma3ff2HfftZZ(2);
1519 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1521 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:ff2H2ff(t:WW)")) {
1522 sigmaPtr = new Sigma3ff2HfftWW(2);
1523 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1525 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:gg2H2ttbar")) {
1526 sigmaPtr = new Sigma3gg2HQQbar(6,2);
1527 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1529 if (HiggsesBSM || HiggsesH2 || settings.flag("HiggsBSM:qqbar2H2ttbar")) {
1530 sigmaPtr = new Sigma3qqbar2HQQbar(6,2);
1531 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1534 // Further BSM H2 processes, not included in "all".
1535 if (settings.flag("HiggsBSM:qg2H2q")) {
1536 sigmaPtr = new Sigma2qg2Hq(4,2);
1537 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1538 sigmaPtr = new Sigma2qg2Hq(5,2);
1539 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1541 if (settings.flag("HiggsBSM:gg2H2bbbar")) {
1542 sigmaPtr = new Sigma3gg2HQQbar(5,2);
1543 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1545 if (settings.flag("HiggsBSM:qqbar2H2bbbar")) {
1546 sigmaPtr = new Sigma3qqbar2HQQbar(5,2);
1547 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1549 if (settings.flag("HiggsBSM:gg2H2g(l:t)")) {
1550 sigmaPtr = new Sigma2gg2Hglt(2);
1551 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1553 if (settings.flag("HiggsBSM:qg2H2q(l:t)")) {
1554 sigmaPtr = new Sigma2qg2Hqlt(2);
1555 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1557 if (settings.flag("HiggsBSM:qqbar2H2g(l:t)")) {
1558 sigmaPtr = new Sigma2qqbar2Hglt(2);
1559 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1562 // Set up requested objects for BSM A3 production.
1563 bool HiggsesA3 = settings.flag("HiggsBSM:allA3");
1564 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:ffbar2A3")) {
1565 sigmaPtr = new Sigma1ffbar2H(3);
1566 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1568 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:gg2A3")) {
1569 sigmaPtr = new Sigma1gg2H(3);
1570 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1572 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:gmgm2A3")) {
1573 sigmaPtr = new Sigma1gmgm2H(3);
1574 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1576 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:ffbar2A3Z")) {
1577 sigmaPtr = new Sigma2ffbar2HZ(3);
1578 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1580 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:ffbar2A3W")) {
1581 sigmaPtr = new Sigma2ffbar2HW(3);
1582 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1584 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:ff2A3ff(t:ZZ)")) {
1585 sigmaPtr = new Sigma3ff2HfftZZ(3);
1586 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1588 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:ff2A3ff(t:WW)")) {
1589 sigmaPtr = new Sigma3ff2HfftWW(3);
1590 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1592 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:gg2A3ttbar")) {
1593 sigmaPtr = new Sigma3gg2HQQbar(6,3);
1594 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1596 if (HiggsesBSM || HiggsesA3 || settings.flag("HiggsBSM:qqbar2A3ttbar")) {
1597 sigmaPtr = new Sigma3qqbar2HQQbar(6,3);
1598 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1601 // Further BSM A3 processes, not included in "all".
1602 if (settings.flag("HiggsBSM:qg2A3q")) {
1603 sigmaPtr = new Sigma2qg2Hq(4,3);
1604 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1605 sigmaPtr = new Sigma2qg2Hq(5,3);
1606 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1608 if (settings.flag("HiggsBSM:gg2A3bbbar")) {
1609 sigmaPtr = new Sigma3gg2HQQbar(5,3);
1610 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1612 if (settings.flag("HiggsBSM:qqbar2A3bbbar")) {
1613 sigmaPtr = new Sigma3qqbar2HQQbar(5,3);
1614 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1616 if (settings.flag("HiggsBSM:gg2A3g(l:t)")) {
1617 sigmaPtr = new Sigma2gg2Hglt(3);
1618 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1620 if (settings.flag("HiggsBSM:qg2A3q(l:t)")) {
1621 sigmaPtr = new Sigma2qg2Hqlt(3);
1622 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1624 if (settings.flag("HiggsBSM:qqbar2A3g(l:t)")) {
1625 sigmaPtr = new Sigma2qqbar2Hglt(3);
1626 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1629 // Set up requested objects for Charged Higgs production
1630 bool HiggsesChg = settings.flag("HiggsBSM:allH+-");
1631 if (HiggsesBSM || HiggsesChg || settings.flag("HiggsBSM:ffbar2H+-")) {
1632 sigmaPtr = new Sigma1ffbar2Hchg;
1633 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1635 if (HiggsesBSM || HiggsesChg || settings.flag("HiggsBSM:bg2H+-t")) {
1636 sigmaPtr = new Sigma2qg2Hchgq(6, 1062, "b g -> H+- t");
1637 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1640 // Set up requested objects for Higgs pair-production
1641 bool HiggsesPairs = settings.flag("HiggsBSM:allHpair");
1642 if (HiggsesBSM || HiggsesPairs || settings.flag("HiggsBSM:ffbar2A3H1")) {
1643 sigmaPtr = new Sigma2ffbar2A3H12(1);
1644 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1646 if (HiggsesBSM || HiggsesPairs || settings.flag("HiggsBSM:ffbar2A3H2")) {
1647 sigmaPtr = new Sigma2ffbar2A3H12(2);
1648 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1650 if (HiggsesBSM || HiggsesPairs || settings.flag("HiggsBSM:ffbar2H+-H1")) {
1651 sigmaPtr = new Sigma2ffbar2HchgH12(1);
1652 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1654 if (HiggsesBSM || HiggsesPairs || settings.flag("HiggsBSM:ffbar2H+-H2")) {
1655 sigmaPtr = new Sigma2ffbar2HchgH12(2);
1656 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1658 if (HiggsesBSM || HiggsesPairs || settings.flag("HiggsBSM:ffbar2H+H-")) {
1659 sigmaPtr = new Sigma2ffbar2HposHneg();
1660 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1664 // Set up requested objects for SUSY pair processes.
1665 if(couplings->isSUSY){
1666 CoupSUSY* coupSUSY = (CoupSUSY *) couplings;
1668 bool SUSYs = settings.flag("SUSY:all");
1669 bool nmssm = settings.flag("SLHA:NMSSM");
1671 // Preselected SUSY codes
1672 int codeA = max( abs(settings.mode("SUSY:idA")),
1673 abs(settings.mode("SUSY:idB")));
1674 int codeB = min( abs(settings.mode("SUSY:idA")),
1675 abs(settings.mode("SUSY:idB")));
1677 // MSSM: 4 neutralinos
1679 if (nmssm) nNeut = 5;
1682 if (SUSYs || settings.flag("SUSY:gg2gluinogluino")) {
1683 // Skip if specific codes not asked for
1684 if (codeA == 0 || codeA == 1000021) {
1685 if (codeB == 0 || codeB == 1000021 ) {
1686 sigmaPtr = new Sigma2gg2gluinogluino();
1687 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1691 if (SUSYs || settings.flag("SUSY:qqbar2gluinogluino")) {
1692 // Skip if specific codes not asked for
1693 if (codeA == 0 || codeA == 1000021) {
1694 if (codeB == 0 || codeB == 1000021 ) {
1695 sigmaPtr = new Sigma2qqbar2gluinogluino();
1696 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1702 if (SUSYs || settings.flag("SUSY:qg2squarkgluino")) {
1704 for (int idx = 1; idx <= 6; ++idx) {
1705 for (int iso = 1; iso <= 2; ++iso) {
1707 int id3 = iso + ((idx <= 3)
1708 ? 1000000+2*(idx-1) : 2000000+2*(idx-4));
1710 // Skip if specific codes not asked for
1711 if (codeA != 0 && codeA != abs(id3) && codeA != abs(id4)) continue;
1712 if (codeB != 0 && ( codeA != max(abs(id3),abs(id4))
1713 || codeB != min(abs(id3),abs(id4)) ) ) continue;
1714 sigmaPtr = new Sigma2qg2squarkgluino(id3,iproc);
1715 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1720 // Squark-antisquark (gg initiated)
1721 if (SUSYs || settings.flag("SUSY:gg2squarkantisquark")) {
1723 for (int idx = 1; idx <= 6; ++idx) {
1724 for (int iso = 1; iso <= 2; ++iso) {
1726 int id = iso + ((idx <= 3)
1727 ? 1000000+2*(idx-1) : 2000000+2*(idx-4));
1728 // Skip if specific codes not asked for
1729 if (codeA != 0 && codeA != abs(id)) continue;
1730 if (codeA != 0 && codeB != 0 && codeB != abs(id)) continue;
1731 sigmaPtr = new Sigma2gg2squarkantisquark(id,iproc);
1732 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1737 // Squark-antisquark (qqbar initiated)
1738 if (SUSYs || settings.flag("SUSY:qqbar2squarkantisquark")) {
1740 for (int idx = 1; idx <= 6; ++idx) {
1741 for (int iso = 1; iso <= 2; ++iso) {
1742 for (int jso = iso; jso >= 1; --jso) {
1743 for (int jdx = 1; jdx <= 6; ++jdx) {
1744 if (iso == jso && jdx < idx) continue;
1745 int id1 = iso + ((idx <= 3) ? 1000000+2*(idx-1)
1746 : 2000000+2*(idx-4));
1747 int id2 = jso + ((jdx <= 3) ? 1000000+2*(jdx-1)
1748 : 2000000+2*(jdx-4));
1749 // Update process number counter (for ~q~q, +2 if not self-conj)
1750 //if (iproc == 1302) iproc=1310;
1752 if (iso == jso && id1 != id2) iproc++;
1753 // Skip if specific codes not asked for
1754 if (codeA != 0 && codeA != abs(id1)
1755 && codeA != abs(id2)) continue;
1756 if (codeB != 0 && ( codeA != max(abs(id1),abs(id2))
1757 || codeB != min(abs(id1),abs(id2)) ) ) continue;
1758 if (iso == jso && id1 != id2) {
1759 sigmaPtr = new Sigma2qqbar2squarkantisquark(id1,-id2,iproc-1);
1760 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1761 sigmaPtr = new Sigma2qqbar2squarkantisquark(id2,-id1,iproc);
1762 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1764 sigmaPtr = new Sigma2qqbar2squarkantisquark(id1,-id2,iproc);
1765 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1774 if (SUSYs || settings.flag("SUSY:qq2squarksquark")) {
1776 for (int idx = 1; idx <= 6; ++idx) {
1777 for (int iso = 1; iso <= 2; ++iso) {
1778 for (int jso = iso; jso >= 1; jso--) {
1779 for (int jdx = 1; jdx <= 6; ++jdx) {
1780 if (iso == jso && jdx < idx) continue;
1782 int id1 = iso + ((idx <= 3)
1783 ? 1000000+2*(idx-1) : 2000000+2*(idx-4));
1784 int id2 = jso + ((jdx <= 3)
1785 ? 1000000+2*(jdx-1) : 2000000+2*(jdx-4));
1786 // Skip if specific codes not asked for
1787 if (codeA != 0 && codeA != abs(id1) && codeA != abs(id2))
1789 if (codeB != 0 && ( codeA != max(abs(id1),abs(id2))
1790 || codeB != min(abs(id1),abs(id2)) ) ) continue;
1791 sigmaPtr = new Sigma2qq2squarksquark(id1,id2,iproc);
1792 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1799 // Neutralino + squark
1800 if (SUSYs || settings.flag("SUSY:qg2chi0squark")) {
1802 for (int iNeut= 1; iNeut <= nNeut; iNeut++) {
1803 for (int idx = 1; idx <= 6; idx++) {
1805 for (int iso = 1; iso <= 2; iso++) {
1806 if (iso == 2) isUp = true;
1808 int id3 = coupSUSY->idNeut(iNeut);
1809 int id4 = iso + ((idx <= 3)
1810 ? 1000000+2*(idx-1) : 2000000+2*(idx-4));
1811 // Skip if specific codes not asked for
1812 if (codeA != 0 && codeA != abs(id3) && codeA != abs(id4)) continue;
1813 if (codeB != 0 && codeB != min(abs(id3),abs(id4)) ) continue;
1814 if (codeA != 0 && codeA == codeB) continue;
1815 sigmaPtr = new Sigma2qg2chi0squark(iNeut,idx,isUp,iproc);
1816 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1822 // Chargino + squark
1823 if (SUSYs || settings.flag("SUSY:qg2chi+-squark")) {
1825 for (int iChar = 1; iChar <= 2; iChar++) {
1826 for (int idx = 1; idx <= 6; idx++) {
1828 for (int iso = 1; iso <= 2; iso++) {
1829 if (iso == 2) isUp = true;
1831 int id3 = coupSUSY->idChar(iChar);
1832 int id4 = iso + ((idx <= 3)
1833 ? 1000000+2*(idx-1) : 2000000+2*(idx-4));
1834 // Skip if specific codes not asked for
1835 if (codeA != 0 && codeA != abs(id3) && codeA != abs(id4)) continue;
1836 if (codeB != 0 && codeB != min(abs(id3),abs(id4)) ) continue;
1837 if (codeA != 0 && codeA == codeB) continue;
1838 sigmaPtr = new Sigma2qg2charsquark(iChar,idx,isUp,iproc);
1839 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1846 if (SUSYs || settings.flag("SUSY:qqbar2chi0chi0")) {
1848 for (int iNeut2 = 1; iNeut2 <= nNeut; iNeut2++) {
1849 for (int iNeut1 = 1; iNeut1 <= iNeut2; iNeut1++) {
1851 if (codeA != 0 && codeA != abs(coupSUSY->idNeut(iNeut1)) &&
1852 codeA != abs(coupSUSY->idNeut(iNeut2))) continue;
1853 if (codeB != 0 && (codeA != max(abs(coupSUSY->idNeut(iNeut1)),
1854 abs(coupSUSY->idNeut(iNeut2)))
1855 || codeB != min(abs(coupSUSY->idNeut(iNeut1)),
1856 abs(coupSUSY->idNeut(iNeut2)))) ) continue;
1857 sigmaPtr = new Sigma2qqbar2chi0chi0(iNeut1, iNeut2,iproc);
1858 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1863 // Neutralino-Chargino
1864 if (SUSYs || settings.flag("SUSY:qqbar2chi+-chi0")) {
1866 for (int iNeut = 1; iNeut <= nNeut; iNeut++) {
1867 for (int iChar = 1; iChar <= 2; ++iChar) {
1869 if (codeA != 0 && codeA != coupSUSY->idNeut(iNeut)
1870 && codeA != coupSUSY->idChar(iChar)) continue;
1872 && ( codeA != max(coupSUSY->idNeut(iNeut),coupSUSY->idChar(iChar))
1873 || codeB != min(coupSUSY->idNeut(iNeut),coupSUSY->idChar(iChar))
1875 sigmaPtr = new Sigma2qqbar2charchi0( iChar, iNeut, iproc-1);
1876 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1877 sigmaPtr = new Sigma2qqbar2charchi0(-iChar, iNeut, iproc);
1878 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1883 // Chargino-Chargino
1884 if (SUSYs || settings.flag("SUSY:qqbar2chi+chi-")) {
1886 for (int i = 1; i <= 2; ++i) {
1887 for (int j = 1; j <= 2; ++j) {
1889 if (codeA != 0 && codeA != abs(coupSUSY->idChar(i))
1890 && codeA != abs(coupSUSY->idChar(j))) continue;
1892 && ( codeA != max(coupSUSY->idChar(i),coupSUSY->idChar(j))
1893 || codeB != min(coupSUSY->idChar(i),coupSUSY->idChar(j)) ) )
1895 sigmaPtr = new Sigma2qqbar2charchar( i,-j, iproc);
1896 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1901 // RPV squark production
1902 if(SUSYs || settings.flag("SUSY:qq2antisquark")) {
1903 for (int idx = 1; idx <= 6; ++idx) {
1904 for (int iso = 1; iso <= 2; ++iso) {
1905 int id1 = iso + ((idx <= 3) ? 1000000+2*(idx-1) : 2000000+2*(idx-4));
1906 if(codeA !=0 && codeA != abs(id1)) continue;
1907 sigmaPtr = new Sigma1qq2antisquark(id1);
1908 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1915 // Set up requested objects for New-Gauge-Boson processes.
1916 if (settings.flag("NewGaugeBoson:ffbar2gmZZprime")) {
1917 sigmaPtr = new Sigma1ffbar2gmZZprime();
1918 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1920 if (settings.flag("NewGaugeBoson:ffbar2Wprime")) {
1921 sigmaPtr = new Sigma1ffbar2Wprime();
1922 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1924 if (settings.flag("NewGaugeBoson:ffbar2R0")) {
1925 sigmaPtr = new Sigma1ffbar2Rhorizontal();
1926 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1929 // Set up requested objects for Left-Right-Symmetry processes.
1930 bool leftrights = settings.flag("LeftRightSymmmetry:all");
1931 if (leftrights || settings.flag("LeftRightSymmmetry:ffbar2ZR")) {
1932 sigmaPtr = new Sigma1ffbar2ZRight();
1933 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1935 if (leftrights || settings.flag("LeftRightSymmmetry:ffbar2WR")) {
1936 sigmaPtr = new Sigma1ffbar2WRight();
1937 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1939 if (leftrights || settings.flag("LeftRightSymmmetry:ll2HL")) {
1940 sigmaPtr = new Sigma1ll2Hchgchg(1);
1941 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1943 if (leftrights || settings.flag("LeftRightSymmmetry:lgm2HLe")) {
1944 sigmaPtr = new Sigma2lgm2Hchgchgl(1, 11);
1945 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1947 if (leftrights || settings.flag("LeftRightSymmmetry:lgm2HLmu")) {
1948 sigmaPtr = new Sigma2lgm2Hchgchgl(1, 13);
1949 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1951 if (leftrights || settings.flag("LeftRightSymmmetry:lgm2HLtau")) {
1952 sigmaPtr = new Sigma2lgm2Hchgchgl(1, 15);
1953 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1955 if (leftrights || settings.flag("LeftRightSymmmetry:ff2HLff")) {
1956 sigmaPtr = new Sigma3ff2HchgchgfftWW(1);
1957 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1959 if (leftrights || settings.flag("LeftRightSymmmetry:ffbar2HLHL")) {
1960 sigmaPtr = new Sigma2ffbar2HchgchgHchgchg(1);
1961 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1963 if (leftrights || settings.flag("LeftRightSymmmetry:ll2HR")) {
1964 sigmaPtr = new Sigma1ll2Hchgchg(2);
1965 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1967 if (leftrights || settings.flag("LeftRightSymmmetry:lgm2HRe")) {
1968 sigmaPtr = new Sigma2lgm2Hchgchgl(2, 11);
1969 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1971 if (leftrights || settings.flag("LeftRightSymmmetry:lgm2HRmu")) {
1972 sigmaPtr = new Sigma2lgm2Hchgchgl(2, 13);
1973 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1975 if (leftrights || settings.flag("LeftRightSymmmetry:lgm2HRtau")) {
1976 sigmaPtr = new Sigma2lgm2Hchgchgl(2, 15);
1977 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1979 if (leftrights || settings.flag("LeftRightSymmmetry:ff2HRff")) {
1980 sigmaPtr = new Sigma3ff2HchgchgfftWW(2);
1981 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1983 if (leftrights || settings.flag("LeftRightSymmmetry:ffbar2HRHR")) {
1984 sigmaPtr = new Sigma2ffbar2HchgchgHchgchg(2);
1985 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1988 // Set up requested objects for leptoquark LQ processes.
1989 bool leptoquarks = settings.flag("LeptoQuark:all");
1990 if (leptoquarks || settings.flag("LeptoQuark:ql2LQ")) {
1991 sigmaPtr = new Sigma1ql2LeptoQuark;
1992 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1994 if (leptoquarks || settings.flag("LeptoQuark:qg2LQl")) {
1995 sigmaPtr = new Sigma2qg2LeptoQuarkl;
1996 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
1998 if (leptoquarks || settings.flag("LeptoQuark:gg2LQLQbar")) {
1999 sigmaPtr = new Sigma2gg2LQLQbar;
2000 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2002 if (leptoquarks || settings.flag("LeptoQuark:qqbar2LQLQbar")) {
2003 sigmaPtr = new Sigma2qqbar2LQLQbar;
2004 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2007 // Set up requested objects for excited-fermion processes.
2008 bool excitedfermions = settings.flag("ExcitedFermion:all");
2009 if (excitedfermions || settings.flag("ExcitedFermion:dg2dStar")) {
2010 sigmaPtr = new Sigma1qg2qStar(1);
2011 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2013 if (excitedfermions || settings.flag("ExcitedFermion:ug2uStar")) {
2014 sigmaPtr = new Sigma1qg2qStar(2);
2015 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2017 if (excitedfermions || settings.flag("ExcitedFermion:sg2sStar")) {
2018 sigmaPtr = new Sigma1qg2qStar(3);
2019 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2021 if (excitedfermions || settings.flag("ExcitedFermion:cg2cStar")) {
2022 sigmaPtr = new Sigma1qg2qStar(4);
2023 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2025 if (excitedfermions || settings.flag("ExcitedFermion:bg2bStar")) {
2026 sigmaPtr = new Sigma1qg2qStar(5);
2027 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2029 if (excitedfermions || settings.flag("ExcitedFermion:egm2eStar")) {
2030 sigmaPtr = new Sigma1lgm2lStar(11);
2031 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2033 if (excitedfermions || settings.flag("ExcitedFermion:mugm2muStar")) {
2034 sigmaPtr = new Sigma1lgm2lStar(13);
2035 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2037 if (excitedfermions || settings.flag("ExcitedFermion:taugm2tauStar")) {
2038 sigmaPtr = new Sigma1lgm2lStar(15);
2039 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2041 if (excitedfermions || settings.flag("ExcitedFermion:qq2dStarq")) {
2042 sigmaPtr = new Sigma2qq2qStarq(1);
2043 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2045 if (excitedfermions || settings.flag("ExcitedFermion:qq2uStarq")) {
2046 sigmaPtr = new Sigma2qq2qStarq(2);
2047 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2049 if (excitedfermions || settings.flag("ExcitedFermion:qq2sStarq")) {
2050 sigmaPtr = new Sigma2qq2qStarq(3);
2051 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2053 if (excitedfermions || settings.flag("ExcitedFermion:qq2cStarq")) {
2054 sigmaPtr = new Sigma2qq2qStarq(4);
2055 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2057 if (excitedfermions || settings.flag("ExcitedFermion:qq2bStarq")) {
2058 sigmaPtr = new Sigma2qq2qStarq(5);
2059 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2061 if (excitedfermions || settings.flag("ExcitedFermion:qqbar2eStare")) {
2062 sigmaPtr = new Sigma2qqbar2lStarlbar(11);
2063 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2065 if (excitedfermions || settings.flag("ExcitedFermion:qqbar2nueStarnue")) {
2066 sigmaPtr = new Sigma2qqbar2lStarlbar(12);
2067 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2069 if (excitedfermions || settings.flag("ExcitedFermion:qqbar2muStarmu")) {
2070 sigmaPtr = new Sigma2qqbar2lStarlbar(13);
2071 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2073 if (excitedfermions || settings.flag("ExcitedFermion:qqbar2numuStarnumu")) {
2074 sigmaPtr = new Sigma2qqbar2lStarlbar(14);
2075 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2077 if (excitedfermions || settings.flag("ExcitedFermion:qqbar2tauStartau")) {
2078 sigmaPtr = new Sigma2qqbar2lStarlbar(15);
2079 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2082 || settings.flag("ExcitedFermion:qqbar2nutauStarnutau")) {
2083 sigmaPtr = new Sigma2qqbar2lStarlbar(16);
2084 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2087 // Set up requested objects for contact interaction processes.
2088 if (settings.flag("ContactInteractions:QCqq2qq")) {
2089 sigmaPtr = new Sigma2QCqq2qq();
2090 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2092 if (settings.flag("ContactInteractions:QCqqbar2qqbar")) {
2093 sigmaPtr = new Sigma2QCqqbar2qqbar();
2094 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2096 if (settings.flag("ContactInteractions:QCffbar2eebar")) {
2097 sigmaPtr = new Sigma2QCffbar2llbar(11, 4003);
2098 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2100 if (settings.flag("ContactInteractions:QCffbar2mumubar")) {
2101 sigmaPtr = new Sigma2QCffbar2llbar(13, 4004);
2102 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2104 if (settings.flag("ContactInteractions:QCffbar2tautaubar")) {
2105 sigmaPtr = new Sigma2QCffbar2llbar(15, 4005);
2106 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2109 // Set spin of particles in the Hidden Valley scenario.
2110 int spinFv = settings.mode("HiddenValley:spinFv");
2111 for (int i = 1; i < 7; ++i) {
2112 if (particleDataPtr->spinType( 4900000 + i) != spinFv + 1)
2113 particleDataPtr->spinType( 4900000 + i, spinFv + 1);
2114 if (particleDataPtr->spinType( 4900010 + i) != spinFv + 1)
2115 particleDataPtr->spinType( 4900010 + i, spinFv + 1);
2118 if (particleDataPtr->spinType( 4900101) != 2)
2119 particleDataPtr->spinType( 4900101, 2);
2121 int spinqv = settings.mode("HiddenValley:spinqv");
2122 if (particleDataPtr->spinType( 4900101) != 2 * spinqv + 1)
2123 particleDataPtr->spinType( 4900101, 2 * spinqv + 1);
2126 // Set up requested objects for HiddenValley processes.
2127 bool hiddenvalleys = settings.flag("HiddenValley:all");
2128 if (hiddenvalleys || settings.flag("HiddenValley:gg2DvDvbar")) {
2129 sigmaPtr = new Sigma2gg2qGqGbar( 4900001, 4901, spinFv,
2131 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2133 if (hiddenvalleys || settings.flag("HiddenValley:gg2UvUvbar")) {
2134 sigmaPtr = new Sigma2gg2qGqGbar( 4900002, 4902, spinFv,
2136 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2138 if (hiddenvalleys || settings.flag("HiddenValley:gg2SvSvbar")) {
2139 sigmaPtr = new Sigma2gg2qGqGbar( 4900003, 4903, spinFv,
2141 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2143 if (hiddenvalleys || settings.flag("HiddenValley:gg2CvCvbar")) {
2144 sigmaPtr = new Sigma2gg2qGqGbar( 4900004, 4904, spinFv,
2146 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2148 if (hiddenvalleys || settings.flag("HiddenValley:gg2BvBvbar")) {
2149 sigmaPtr = new Sigma2gg2qGqGbar( 4900005, 4905, spinFv,
2151 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2153 if (hiddenvalleys || settings.flag("HiddenValley:gg2TvTvbar")) {
2154 sigmaPtr = new Sigma2gg2qGqGbar( 4900006, 4906, spinFv,
2156 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2158 if (hiddenvalleys || settings.flag("HiddenValley:qqbar2DvDvbar")) {
2159 sigmaPtr = new Sigma2qqbar2qGqGbar( 4900001, 4911, spinFv,
2160 "q qbar -> Dv Dvbar");
2161 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2163 if (hiddenvalleys || settings.flag("HiddenValley:qqbar2UvUvbar")) {
2164 sigmaPtr = new Sigma2qqbar2qGqGbar( 4900002, 4912, spinFv,
2165 "q qbar -> Uv Uvbar");
2166 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2168 if (hiddenvalleys || settings.flag("HiddenValley:qqbar2SvSvbar")) {
2169 sigmaPtr = new Sigma2qqbar2qGqGbar( 4900003, 4913, spinFv,
2170 "q qbar -> Sv Svbar");
2171 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2173 if (hiddenvalleys || settings.flag("HiddenValley:qqbar2CvCvbar")) {
2174 sigmaPtr = new Sigma2qqbar2qGqGbar( 4900004, 4914, spinFv,
2175 "q qbar -> Cv Cvbar");
2176 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2178 if (hiddenvalleys || settings.flag("HiddenValley:qqbar2BvBvbar")) {
2179 sigmaPtr = new Sigma2qqbar2qGqGbar( 4900005, 4915, spinFv,
2180 "q qbar -> Bv Bvbar");
2181 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2183 if (hiddenvalleys || settings.flag("HiddenValley:qqbar2TvTvbar")) {
2184 sigmaPtr = new Sigma2qqbar2qGqGbar( 4900006, 4916, spinFv,
2185 "q qbar -> Tv Tvbar");
2186 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2188 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2DvDvbar")) {
2189 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900001, 4921, spinFv,
2190 "f fbar -> Dv Dvbar");
2191 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2193 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2UvUvbar")) {
2194 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900002, 4922, spinFv,
2195 "f fbar -> Uv Uvbar");
2196 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2198 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2SvSvbar")) {
2199 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900003, 4923, spinFv,
2200 "f fbar -> Sv Svbar");
2201 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2203 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2CvCvbar")) {
2204 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900004, 4924, spinFv,
2205 "f fbar -> Cv Cvbar");
2206 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2208 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2BvBvbar")) {
2209 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900005, 4925, spinFv,
2210 "f fbar -> Bv Bvbar");
2211 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2213 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2TvTvbar")) {
2214 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900006, 4926, spinFv,
2215 "f fbar -> Tv Tvbar");
2216 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2218 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2EvEvbar")) {
2219 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900011, 4931, spinFv,
2220 "f fbar -> Ev Evbar");
2221 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2223 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2nuEvnuEvbar")) {
2224 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900012, 4932, spinFv,
2225 "f fbar -> nuEv nuEvbar");
2226 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2228 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2MUvMUvbar")) {
2229 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900013, 4933, spinFv,
2230 "f fbar -> MUv MUvbar");
2231 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2233 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2nuMUvnuMUvbar")) {
2234 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900014, 4934, spinFv,
2235 "f fbar -> nuMUv nuMUvbar");
2236 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2238 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2TAUvTAUvbar")) {
2239 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900015, 4935, spinFv,
2240 "f fbar -> TAUv TAUvbar");
2241 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2243 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2nuTAUvnuTAUvbar")) {
2244 sigmaPtr = new Sigma2ffbar2fGfGbar( 4900016, 4936, spinFv,
2245 "f fbar -> nuTAUv nuTAUvbar");
2246 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2248 if (hiddenvalleys || settings.flag("HiddenValley:ffbar2Zv")) {
2249 sigmaPtr = new Sigma1ffbar2Zv();
2250 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2253 // Set up requested objects for RS extra-dimensional G* processes.
2254 bool extraDimGstars = settings.flag("ExtraDimensionsG*:all");
2255 if (extraDimGstars || settings.flag("ExtraDimensionsG*:gg2G*")) {
2256 sigmaPtr = new Sigma1gg2GravitonStar;
2257 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2259 if (extraDimGstars || settings.flag("ExtraDimensionsG*:ffbar2G*")) {
2260 sigmaPtr = new Sigma1ffbar2GravitonStar;
2261 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2263 if (settings.flag("ExtraDimensionsG*:gg2G*g")) {
2264 sigmaPtr = new Sigma2gg2GravitonStarg;
2265 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2267 if (settings.flag("ExtraDimensionsG*:qg2G*q")) {
2268 sigmaPtr = new Sigma2qg2GravitonStarq;
2269 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2271 if (settings.flag("ExtraDimensionsG*:qqbar2G*g")) {
2272 sigmaPtr = new Sigma2qqbar2GravitonStarg;
2273 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2276 // Set up requested objects for RS extra-dimensional KKgluon processes.
2277 if (settings.flag("ExtraDimensionsG*:qqbar2KKgluon*")) {
2278 sigmaPtr = new Sigma1qqbar2KKgluonStar;
2279 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2282 // NOAM: Set up requested objects for TEV extra-dimensional processes.
2283 if (settings.flag("ExtraDimensionsTEV:ffbar2ddbar")) {
2284 sigmaPtr = new Sigma2ffbar2TEVffbar(1, 5061);
2285 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2287 if (settings.flag("ExtraDimensionsTEV:ffbar2uubar")) {
2288 sigmaPtr = new Sigma2ffbar2TEVffbar(2, 5062);
2289 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2291 if (settings.flag("ExtraDimensionsTEV:ffbar2ssbar")) {
2292 sigmaPtr = new Sigma2ffbar2TEVffbar(3, 5063);
2293 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2295 if (settings.flag("ExtraDimensionsTEV:ffbar2ccbar")) {
2296 sigmaPtr = new Sigma2ffbar2TEVffbar(4, 5064);
2297 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2299 if (settings.flag("ExtraDimensionsTEV:ffbar2bbbar")) {
2300 sigmaPtr = new Sigma2ffbar2TEVffbar(5, 5065);
2301 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2303 if (settings.flag("ExtraDimensionsTEV:ffbar2ttbar")) {
2304 sigmaPtr = new Sigma2ffbar2TEVffbar(6, 5066);
2305 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2307 if (settings.flag("ExtraDimensionsTEV:ffbar2e+e-")) {
2308 sigmaPtr = new Sigma2ffbar2TEVffbar(11, 5071);
2309 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2311 if (settings.flag("ExtraDimensionsTEV:ffbar2nuenuebar")) {
2312 sigmaPtr = new Sigma2ffbar2TEVffbar(12, 5072);
2313 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2315 if (settings.flag("ExtraDimensionsTEV:ffbar2mu+mu-")) {
2316 sigmaPtr = new Sigma2ffbar2TEVffbar(13, 5073);
2317 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2319 if (settings.flag("ExtraDimensionsTEV:ffbar2numunumubar")) {
2320 sigmaPtr = new Sigma2ffbar2TEVffbar(14, 5074);
2321 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2323 if (settings.flag("ExtraDimensionsTEV:ffbar2tau+tau-")) {
2324 sigmaPtr = new Sigma2ffbar2TEVffbar(15, 5075);
2325 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2327 if (settings.flag("ExtraDimensionsTEV:ffbar2nutaunutaubar")) {
2328 sigmaPtr = new Sigma2ffbar2TEVffbar(16, 5076);
2329 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2332 // Set up requested objects for large extra-dimensional G processes.
2333 bool extraDimLEDmono = settings.flag("ExtraDimensionsLED:monojet");
2334 if (extraDimLEDmono || settings.flag("ExtraDimensionsLED:gg2Gg")) {
2335 sigmaPtr = new Sigma2gg2LEDUnparticleg( true );
2336 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2338 if (extraDimLEDmono || settings.flag("ExtraDimensionsLED:qg2Gq")) {
2339 sigmaPtr = new Sigma2qg2LEDUnparticleq( true );
2340 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2342 if (extraDimLEDmono || settings.flag("ExtraDimensionsLED:qqbar2Gg")) {
2343 sigmaPtr = new Sigma2qqbar2LEDUnparticleg( true );
2344 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2346 if (settings.flag("ExtraDimensionsLED:ffbar2GZ")) {
2347 sigmaPtr = new Sigma2ffbar2LEDUnparticleZ( true );
2348 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2350 if (settings.flag("ExtraDimensionsLED:ffbar2Ggamma")) {
2351 sigmaPtr = new Sigma2ffbar2LEDUnparticlegamma( true );
2352 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2354 if (settings.flag("ExtraDimensionsLED:ffbar2gammagamma")) {
2355 sigmaPtr = new Sigma2ffbar2LEDgammagamma( true );
2356 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2358 if (settings.flag("ExtraDimensionsLED:gg2gammagamma")) {
2359 sigmaPtr = new Sigma2gg2LEDgammagamma( true );
2360 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2362 if (settings.flag("ExtraDimensionsLED:ffbar2llbar")) {
2363 sigmaPtr = new Sigma2ffbar2LEDllbar( true );
2364 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2366 if (settings.flag("ExtraDimensionsLED:gg2llbar")) {
2367 sigmaPtr = new Sigma2gg2LEDllbar( true );
2368 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2370 bool extraDimLEDdij = settings.flag("ExtraDimensionsLED:dijets");
2371 if (extraDimLEDdij || settings.flag("ExtraDimensionsLED:gg2DJgg")) {
2372 sigmaPtr = new Sigma2gg2LEDgg;
2373 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2375 if (extraDimLEDdij || settings.flag("ExtraDimensionsLED:gg2DJqqbar")) {
2376 sigmaPtr = new Sigma2gg2LEDqqbar;
2377 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2379 if (extraDimLEDdij || settings.flag("ExtraDimensionsLED:qg2DJqg")) {
2380 sigmaPtr = new Sigma2qg2LEDqg;
2381 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2383 if (extraDimLEDdij || settings.flag("ExtraDimensionsLED:qq2DJqq")) {
2384 sigmaPtr = new Sigma2qq2LEDqq;
2385 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2387 if (extraDimLEDdij || settings.flag("ExtraDimensionsLED:qqbar2DJgg")) {
2388 sigmaPtr = new Sigma2qqbar2LEDgg;
2389 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2391 if (extraDimLEDdij || settings.flag("ExtraDimensionsLED:qqbar2DJqqbarNew")) {
2392 sigmaPtr = new Sigma2qqbar2LEDqqbarNew;
2393 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2396 // Set up requested objects for unparticle processes.
2397 bool extraDimUnpartmono = settings.flag("ExtraDimensionsUnpart:monojet");
2398 if (extraDimUnpartmono || settings.flag("ExtraDimensionsUnpart:gg2Ug")) {
2399 sigmaPtr = new Sigma2gg2LEDUnparticleg( false );
2400 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2402 if (extraDimUnpartmono || settings.flag("ExtraDimensionsUnpart:qg2Uq")) {
2403 sigmaPtr = new Sigma2qg2LEDUnparticleq( false );
2404 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2406 if (extraDimUnpartmono || settings.flag("ExtraDimensionsUnpart:qqbar2Ug")) {
2407 sigmaPtr = new Sigma2qqbar2LEDUnparticleg( false );
2408 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2410 if (settings.flag("ExtraDimensionsUnpart:ffbar2UZ")) {
2411 sigmaPtr = new Sigma2ffbar2LEDUnparticleZ( false );
2412 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2414 if (settings.flag("ExtraDimensionsUnpart:ffbar2Ugamma")) {
2415 sigmaPtr = new Sigma2ffbar2LEDUnparticlegamma( false );
2416 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2418 if (settings.flag("ExtraDimensionsUnpart:ffbar2gammagamma")) {
2419 sigmaPtr = new Sigma2ffbar2LEDgammagamma( false );
2420 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2422 if (settings.flag("ExtraDimensionsUnpart:gg2gammagamma")) {
2423 sigmaPtr = new Sigma2gg2LEDgammagamma( false );
2424 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2426 if (settings.flag("ExtraDimensionsUnpart:ffbar2llbar")) {
2427 sigmaPtr = new Sigma2ffbar2LEDllbar( false );
2428 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2430 if (settings.flag("ExtraDimensionsUnpart:gg2llbar")) {
2431 sigmaPtr = new Sigma2gg2LEDllbar( false );
2432 containerPtrs.push_back( new ProcessContainer(sigmaPtr) );
2440 //--------------------------------------------------------------------------
2442 // Routine to initialize list of second hard processes.
2444 bool SetupContainers::init2(vector<ProcessContainer*>& container2Ptrs,
2445 Settings& settings) {
2447 // Reset process list, if filled in previous subrun.
2448 if (container2Ptrs.size() > 0) {
2449 for (int i = 0; i < int(container2Ptrs.size()); ++i)
2450 delete container2Ptrs[i];
2451 container2Ptrs.clear();
2453 SigmaProcess* sigmaPtr;
2455 // Two hard QCD jets.
2456 if (settings.flag("SecondHard:TwoJets")) {
2457 sigmaPtr = new Sigma2gg2gg;
2458 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2459 sigmaPtr = new Sigma2gg2qqbar;
2460 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2461 sigmaPtr = new Sigma2qg2qg;
2462 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2463 sigmaPtr = new Sigma2qq2qq;
2464 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2465 sigmaPtr = new Sigma2qqbar2gg;
2466 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2467 sigmaPtr = new Sigma2qqbar2qqbarNew;
2468 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2469 sigmaPtr = new Sigma2gg2QQbar(4, 121);
2470 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2471 sigmaPtr = new Sigma2qqbar2QQbar(4, 122);
2472 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2473 sigmaPtr = new Sigma2gg2QQbar(5, 123);
2474 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2475 sigmaPtr = new Sigma2qqbar2QQbar(5, 124);
2476 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2479 // A prompt photon and a hard jet.
2480 if (settings.flag("SecondHard:PhotonAndJet")) {
2481 sigmaPtr = new Sigma2qg2qgamma;
2482 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2483 sigmaPtr = new Sigma2qqbar2ggamma;
2484 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2485 sigmaPtr = new Sigma2gg2ggamma;
2486 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2489 // Two prompt photons.
2490 if (settings.flag("SecondHard:TwoPhotons")) {
2491 sigmaPtr = new Sigma2ffbar2gammagamma;
2492 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2493 sigmaPtr = new Sigma2gg2gammagamma;
2494 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2497 // Charmonium production.
2498 if (settings.flag("SecondHard:Charmonium")) {
2499 sigmaPtr = new Sigma2gg2QQbar3S11g(4, 401);
2500 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2501 sigmaPtr = new Sigma2gg2QQbar3PJ1g(4, 0, 402);
2502 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2503 sigmaPtr = new Sigma2gg2QQbar3PJ1g(4, 1, 403);
2504 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2505 sigmaPtr = new Sigma2gg2QQbar3PJ1g(4, 2, 404);
2506 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2507 sigmaPtr = new Sigma2qg2QQbar3PJ1q(4, 0, 405);
2508 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2509 sigmaPtr = new Sigma2qg2QQbar3PJ1q(4, 1, 406);
2510 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2511 sigmaPtr = new Sigma2qg2QQbar3PJ1q(4, 2, 407);
2512 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2513 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(4, 0, 408);
2514 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2515 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(4, 1, 409);
2516 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2517 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(4, 2, 410);
2518 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2519 sigmaPtr = new Sigma2gg2QQbarX8g(4, 0, 411);
2520 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2521 sigmaPtr = new Sigma2gg2QQbarX8g(4, 1, 412);
2522 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2523 sigmaPtr = new Sigma2gg2QQbarX8g(4, 2, 413);
2524 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2525 sigmaPtr = new Sigma2qg2QQbarX8q(4, 0, 414);
2526 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2527 sigmaPtr = new Sigma2qg2QQbarX8q(4, 1, 415);
2528 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2529 sigmaPtr = new Sigma2qg2QQbarX8q(4, 2, 416);
2530 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2531 sigmaPtr = new Sigma2qqbar2QQbarX8g(4, 0, 417);
2532 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2533 sigmaPtr = new Sigma2qqbar2QQbarX8g(4, 1, 418);
2534 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2535 sigmaPtr = new Sigma2qqbar2QQbarX8g(4, 2, 419);
2536 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2539 // Bottomonium production.
2540 if (settings.flag("SecondHard:Bottomonium")) {
2541 sigmaPtr = new Sigma2gg2QQbar3S11g(5, 501);
2542 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2543 sigmaPtr = new Sigma2gg2QQbar3PJ1g(5, 0, 502);
2544 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2545 sigmaPtr = new Sigma2gg2QQbar3PJ1g(5, 1, 503);
2546 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2547 sigmaPtr = new Sigma2gg2QQbar3PJ1g(5, 2, 504);
2548 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2549 sigmaPtr = new Sigma2qg2QQbar3PJ1q(5, 0, 505);
2550 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2551 sigmaPtr = new Sigma2qg2QQbar3PJ1q(5, 1, 506);
2552 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2553 sigmaPtr = new Sigma2qg2QQbar3PJ1q(5, 2, 507);
2554 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2555 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(5, 0, 508);
2556 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2557 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(5, 1, 509);
2558 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2559 sigmaPtr = new Sigma2qqbar2QQbar3PJ1g(5, 2, 510);
2560 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2561 sigmaPtr = new Sigma2gg2QQbarX8g(5, 0, 511);
2562 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2563 sigmaPtr = new Sigma2gg2QQbarX8g(5, 1, 512);
2564 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2565 sigmaPtr = new Sigma2gg2QQbarX8g(5, 2, 513);
2566 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2567 sigmaPtr = new Sigma2qg2QQbarX8q(5, 0, 514);
2568 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2569 sigmaPtr = new Sigma2qg2QQbarX8q(5, 1, 515);
2570 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2571 sigmaPtr = new Sigma2qg2QQbarX8q(5, 2, 516);
2572 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2573 sigmaPtr = new Sigma2qqbar2QQbarX8g(5, 0, 517);
2574 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2575 sigmaPtr = new Sigma2qqbar2QQbarX8g(5, 1, 518);
2576 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2577 sigmaPtr = new Sigma2qqbar2QQbarX8g(5, 2, 519);
2578 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2581 // A single gamma*/Z0.
2582 if (settings.flag("SecondHard:SingleGmZ")) {
2583 sigmaPtr = new Sigma1ffbar2gmZ;
2584 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2588 if (settings.flag("SecondHard:SingleW")) {
2589 sigmaPtr = new Sigma1ffbar2W;
2590 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2593 // A gamma*/Z0 and a hard jet.
2594 if (settings.flag("SecondHard:GmZAndJet")) {
2595 sigmaPtr = new Sigma2qqbar2gmZg;
2596 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2597 sigmaPtr = new Sigma2qg2gmZq;
2598 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2601 // A W+- and a hard jet.
2602 if (settings.flag("SecondHard:WAndJet")) {
2603 sigmaPtr = new Sigma2qqbar2Wg;
2604 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2605 sigmaPtr = new Sigma2qg2Wq;
2606 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2609 // Top pair production.
2610 if (settings.flag("SecondHard:TopPair")) {
2611 sigmaPtr = new Sigma2gg2QQbar(6, 601);
2612 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2613 sigmaPtr = new Sigma2qqbar2QQbar(6, 602);
2614 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2615 sigmaPtr = new Sigma2ffbar2FFbarsgmZ(6, 604);
2616 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2619 // Single top production.
2620 if (settings.flag("SecondHard:SingleTop")) {
2621 sigmaPtr = new Sigma2qq2QqtW(6, 603);
2622 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2623 sigmaPtr = new Sigma2ffbar2FfbarsW(6, 0, 605);
2624 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2627 // Two b jets - already part of TwoJets sample above.
2628 if (settings.flag("SecondHard:TwoBJets")) {
2629 sigmaPtr = new Sigma2gg2QQbar(5, 123);
2630 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2631 sigmaPtr = new Sigma2qqbar2QQbar(5, 124);
2632 container2Ptrs.push_back( new ProcessContainer(sigmaPtr) );
2640 //==========================================================================
2642 } // end namespace Pythia8