66 {
67
68
69 if (!inResonance) return false;
70
71
72 int iTop = partonSystemsPtr->getInRes(iSys);
73 if (iTop == 0 || event[iTop].idAbs() != 6) return false;
74
75
76 int sizeOut = partonSystemsPtr->sizeOut(iSys);
77 if (sizeOut == 2) return false;
78
79
80 int iRad = sizeOld;
81 int iEmt = sizeOld + 1;
82 int iRec = sizeOld + 2;
83
84
85
87 iRad = event[iRad].mother1();
88 iRec = event[iRec].mother1();
89 }
90
91
92 if (event[iEmt].id() != 21) return false;
93 if (event[iTop].id() == 6) {
94 if (event[iEmt].col() != event[iTop].col()) return false;
95 } else {
96 if (event[iEmt].acol() != event[iTop].acol()) return false;
97 }
98
99
100 if (event[iRec].idAbs() != 24) {
101 cout << " ERROR: recoiler is " << event[iRec].id() << endl;
102 return false;
103 }
104
105
106 double pRadRec = event[iRad].p() * event[iRec].p();
107 double pRadEmt = event[iRad].p() * event[iEmt].p();
108 double pRecEmt = event[iRec].p() * event[iEmt].p();
109 double wtW = 2. * pRadRec / (pRadEmt * pRecEmt)
110 -
pow2(event[iRad].
m() / pRadEmt);
111
113
114
115 double pRadTop = event[iRad].p() * event[iTop].p();
116 double pTopEmt = event[iTop].p() * event[iEmt].p();
117 double wtT = 2. * pRadTop / (pRadEmt * pTopEmt)
118 -
pow2(event[iRad].
m() / pRadEmt) -
pow2(event[iTop].
m() / pTopEmt);
119
120
122
123
125 std::ios oldState(nullptr);
126 std::cout << "\n now event with sizeOld = " << sizeOld << ", iSys = "
127 << iSys << ", sizeOut = " << sizeOut << scientific
128 << setprecision(3)
129 << ", weight with W = " << wtW << " and with t = " << wtT << std::endl;
130 partonSystemsPtr->list();
131 event.list();
132 std::cout.copyfmt(oldState);
133 }
134
135
136 return (wtT < wtW * rndmPtr->flat());
137 }
static const std::map< unsigned int, unsigned int > pow2
1.14.0