]>
Commit | Line | Data |
---|---|---|
9419eeef | 1 | <chapter name="Semi-Internal Processes"> |
2 | ||
3 | <h2>Semi-Internal Processes</h2> | |
4 | ||
5 | Normally users are expected to implement new processes via the | |
6 | <aloc href="LesHouchesAccord">Les Houches Accord</aloc>. Then | |
7 | you do all flavour, colour and phase-space selection externally, | |
8 | before your process-level events are input for further processing | |
9 | by PYTHIA. However, it is also possible to implement a | |
10 | new process in exactly the same way as the internal PYTHIA | |
11 | ones, thus making use of the internal phase space selection machinery | |
12 | to sample an externally provided cross-section expression. | |
13 | This page gives a brief summary how to do that. If you additionally | |
14 | want to introduce a new resonance species, with its own internal | |
15 | width calculations, you will find further instructions | |
16 | <aloc href="SemiInternalResonances">here</aloc>. | |
17 | ||
18 | <p/> | |
19 | Should you actually go ahead, it is strongly recommended to shop around | |
20 | for a similar process that has already been implemented, and to use that | |
21 | existing code as a template. Look for processes with the same combinations | |
22 | of incoming flavours and colour flows, rather than the shape of the | |
23 | cross section itself. With a reasonable such match the task should be | |
24 | of medium difficulty, without it more demanding. | |
25 | ||
26 | <p/> | |
27 | PYTHIA is rather good at handling the phase space of | |
28 | <ei>2 -> 1</ei> and <ei>2 -> 2</ei> processes, is more primitive for | |
29 | <ei>2 -> 3</ei> ones and does not at all address higher multiplicities. | |
30 | This limits the set of processes that you can implement in this | |
31 | framework. The produced particles may be resonances, however, so it is | |
32 | possible to end up with bigger "final" multiplicities through sequential | |
33 | decays, and to include further matrix-element weighting in those decays. | |
34 | ||
35 | <p/> | |
36 | There are three steps involved in implementing a process: | |
37 | <br/>1) making use of the PYTHIA-provided kinematics information to | |
38 | calculate the relevant cross section, | |
39 | <br/>2) writing a new class, where the matrix elements are implemented, | |
40 | including information on incoming and outgoing flavours and colours, and | |
41 | <br/>3) making the process available. | |
42 | <br/>We consider these aspects in turn. An example where it all comes | |
43 | together is found in <code>main25.cc</code>. | |
44 | ||
45 | <h3>The Cross Section Calculation</h3> | |
46 | ||
47 | The key method for the cross section calculation is | |
48 | <code>SigmaProcess::sigmaHat()</code>, described below. At the point when | |
49 | it is called, the kinematics has already been set up, and from these | |
50 | phase space variables the differential cross section is to be calculated. | |
51 | ||
52 | <p/> | |
53 | For a <ei>2 -> 1</ei> process, the returned value should be | |
54 | <ei>sigmaHat(sHat)</ei>, where <code>mH</code> (= <ei>mHat</ei>), | |
55 | <code>sH</code> (= <ei>sHat</ei>) and <code>sH2</code> (= <ei>sHat^2</ei>) | |
56 | are available to be used. Incoming partons are massless. Overload the | |
57 | <code>convertME()</code> method below if you instead plan to return | |
58 | <ei>|M|^2</ei>. | |
59 | ||
60 | <p/> | |
61 | For a <ei>2 -> 2</ei> process, instead <ei>d(sigmaHat)/d(tHat)</ei> | |
62 | should be returned, based on provided | |
63 | <code>mH, sH, sH2, tH, tH2, uH, uH2, m3, s3, m4, s4</code> and | |
64 | <code>pT2</code> values (<code>s3 = m3*m3</code> etc.). Incoming | |
65 | partons are massless. Overload the <code>convertME()</code> method | |
66 | below if you instead plan to return <ei>|M|^2</ei>. | |
67 | ||
68 | <p/> | |
69 | For a <ei>2 -> 3</ei> process, instead <ei>|M|^2</ei> should be | |
70 | returned, with normalization such that <ei>|M|^2 / (2 sHat)</ei> integrated | |
71 | over the three-body phase space gives the cross section. Here no standard | |
72 | set of Mandelstam-style variables exists. Instead the obvious ones, | |
73 | <code>mH, sH, m3, s3, m4, s4, m5, s5</code>, are complemented by the | |
74 | four-vectors <code>p3cm, p4cm, p5cm</code>, from which further invariants | |
75 | may be calculated. The four-vectors are defined in the CM frame of the | |
76 | subcollision, with massless incoming partons along the <ei>+-z</ei> axis. | |
77 | ||
78 | <p/> | |
79 | In either case, <ei>alpha_s</ei> and <ei>alpha_em</ei> have already | |
80 | been calculated, and are stored in <code>alpS</code> and <code>alpEM</code>. | |
81 | Also other standard variables may be used, like | |
82 | <code>CoupEW::sin2thetaW()</code>, and related flavour-dependent | |
83 | vector and axial couplings in <code>CoupEW</code> and CKM combinations | |
84 | in <code>VCKM</code>. | |
85 | ||
86 | <p/> | |
87 | In case some of the final-state particles are resonances, their | |
88 | squared masses have already been selected according to a Breit-Wigner | |
89 | with a linearly running width <ei>Gamma(m) = Gamma(m_0) * m / m_0</ei>. | |
90 | More precisely, the mass spectrum is weighted according to | |
91 | <ei>w_BW(m^2) d(m^2)</ei>, where | |
92 | <eq> | |
93 | w_BW(m^2) = (1/pi) * (m * Gamma(m)) / ( (m^2 - m_0^2)^2 + (m * Gamma(m))^2 ) . | |
94 | </eq> | |
95 | If you would like to have another expression, the above weights are stored | |
96 | in <code>runBW3</code>, <code>runBW4</code> and <code>runBW5</code>, | |
97 | respectively. If you divide out one of these factors, you just remain with | |
98 | a phase space selection <ei>d(m^2)</ei> for this particle, | |
99 | and can multiply on your desired shape factor instead. Unfortunately, the | |
100 | Monte Carlo efficiency will drop if your new mass distribution differs | |
101 | dramatically from the input one. Therefore it does make sense to adjust the | |
102 | database value of the width to be slightly (but not too much) broader | |
103 | than the distribution you have in mind. Also note that, already by default, | |
104 | the wings of the Breit-Wigner are oversampled (with a compensating lower | |
105 | internal weight) by partly sampling like <ei>(a + b/m^2 + c/m^4) d(m^2)</ei>, | |
106 | where the last term is only used for <ei>gamma^*/Z^0</ei>. | |
107 | ||
108 | <p/> | |
109 | As alternative to the kinematics variables defined above, also the two | |
110 | arrays <code>mME[5]</code> and <code>pME[5]</code>, for masses and | |
111 | four-momenta, respectively, can be used for cross-section calculations. | |
112 | Here indices 0 and 1 are the two incoming beams, and 2 and onwards the | |
113 | outgoing particles. Note that this differs by one step from the normal | |
114 | internal labelling, where slot 0 is left empty. The four-momenta are | |
115 | defined in the rest frame of the subcollision, with the incoming partons | |
116 | along the <ei>+-z</ei> direction. The kinematics need not agree with the | |
117 | "correct" one stored in the event record, for three reasons. | |
118 | <br/>1) Gauge invariance forces matrix-element calculations to use | |
119 | the same masses for incoming as outgoing legs of a particle species, | |
120 | say <ei>b</ei> quarks. Therefore the kinematics of the two incoming | |
121 | partons is recalculated, relative to the normal event record, to put | |
122 | the partons on the mass shell. (Note that initial masses is a technical | |
123 | issue, not the correct physics picture: the incoming partons are likely | |
124 | to be spacelike virtual rather than on the mass shell.) | |
125 | <br/>2) In principle each fermion flavour has to be treated separately, | |
126 | owing to a different mass. However, in many cases fermions can be | |
127 | assumed massless, which speeds up the calculations, and further gains | |
128 | occur if then different flavours can use the same cross-section | |
129 | expression. In MadGraph the default is that fermions up to and including | |
130 | the <ei>c</ei> quark and the <ei>mu</ei> lepton are considered massless, | |
131 | while the <ei>b</ei> quark and the <ei>tau</ei> lepton are considered | |
132 | massive. This can be modified however, and below we provide four flags | |
133 | that can be used to consider the "borderline" fermions either as | |
134 | massless or as massive when matrix elements are evaluated, to match the | |
135 | assumptions made for the matrix elements themselves. | |
136 | <br/>3) For <ei>2 -> 2</ei> and <ei>2 -> 3</ei> processes of massive | |
137 | identical particles (or antiparticles) in the final state, such as | |
138 | <ei>t tbar</ei> or <ei>W^+ W^-</ei>, the kinematics is here adjusted | |
139 | so that the two or three particles have the same mass, formed as a | |
140 | suitable average of the actual Breit-Wigner-distributed masses. This | |
141 | allows the evaluation of matrix-element expressions that only have | |
142 | meaning if the two/three have the same mass. | |
143 | <br/>Thus the mass array <code>mME[5]</code> and the four-momentum array | |
144 | <code>pME[5]</code> present values both for initial- and final-state | |
145 | particles based on these mass principles suited for matrix-element input. | |
146 | Note that these variables therefore differ from the kinematics stored in | |
147 | the event record proper, where incoming fermions are always massless and | |
148 | outgoing resonances have independent Breit-Wigner mass distributions. | |
149 | <br/>The conversion from the normal to the special kinematics is done | |
150 | by calling the <code>setupForME()</code> method. This you have to do | |
151 | yourself in the <code>SigmaHat()</code> member of your derived class. | |
152 | Alternatively it could be done in <code>SigmaKin()</code>, i.e. before | |
153 | the loop over incoming flavours, but then these would be considered | |
154 | massless. The identity of final-state particles is obtained from the | |
155 | <code>id3Mass()</code>, <code>id4Mass()</code> and <code>id5Mass()</code> | |
156 | methods. Should the conversion to <code>mME[5]</code> and | |
157 | <code>pME[5]</code> not work, <code>setupForME()</code> will return | |
158 | <code>false</code>, and then the cross section should be put zero. | |
159 | ||
160 | <flag name="SigmaProcess:cMassiveME" default="off"> | |
161 | Let the <ei>c</ei> quark be massive or not in the kinematics set up for | |
162 | external matrix-element evaluation. | |
163 | </flag> | |
164 | ||
165 | <flag name="SigmaProcess:bMassiveME" default="on"> | |
166 | Let the <ei>b</ei> quark be massive or not in the kinematics set up for | |
167 | external matrix-element evaluation. | |
168 | </flag> | |
169 | ||
170 | <flag name="SigmaProcess:muMassiveME" default="off"> | |
171 | Let the <ei>mu</ei> lepton be massive or not in the kinematics set up for | |
172 | external matrix-element evaluation. | |
173 | </flag> | |
174 | ||
175 | <flag name="SigmaProcess:tauMassiveME" default="on"> | |
176 | Let the <ei>tau</ei> lepton be massive or not in the kinematics set up for | |
177 | external matrix-element evaluation. | |
178 | </flag> | |
179 | ||
180 | ||
181 | <h3>The Cross Section Class</h3> | |
182 | ||
183 | The matrix-element information has to be encoded in a new class. | |
184 | The relevant code could either be put before the main program in the | |
185 | same file, or be stored separately, e.g. in a matched pair | |
186 | of <code>.h</code> and <code>.cc</code> files. The latter may be more | |
187 | convenient, in particular if the cross sections are lengthy, or if you | |
188 | intend to build up your own little process library, but of course | |
189 | requires that these additional files are correctly compiled and linked. | |
190 | ||
191 | <p/> | |
192 | The class has to be derived either from | |
193 | <code>Sigma1Process</code>, for <ei>2 -> 1</ei> processes, from | |
194 | <code>Sigma2Process</code>, for <ei>2 -> 2</ei> ones, or from | |
195 | <code>Sigma3Process</code>, for <ei>2 -> 3</ei> ones. (The | |
196 | <code>Sigma0Process</code> class is used for elastic, diffractive | |
197 | and minimum-bias events, and is not recommended for use beyond that.) | |
198 | These are in their turn derived from the <code>SigmaProcess</code> | |
199 | base class. | |
200 | ||
201 | <p/> | |
202 | The class can implement a number of methods. Some of these are | |
203 | compulsory, others strongly recommended, and the rest are to be | |
204 | used only when the need arises to override the default behaviour. | |
205 | The methods are: | |
206 | ||
207 | <p/> | |
208 | A <b>constructor</b> for the derived class obviously must be available. | |
209 | Here you are quite free to allow a list of arguments, to set | |
210 | the parameters of your model, or even to create a set of closely | |
211 | related but distinct processes. For instance, <ei>g g -> Q Qbar</ei>, | |
212 | <ei>Q = c</ei> or <ei>b</ei>, is only coded once, and then the | |
213 | constructor takes the quark code (4 or 5) as argument, | |
214 | to allow the proper amount of differentiation. | |
215 | ||
216 | <p/> | |
217 | A <b>destructor</b> is only needed if you plan to delete the process | |
218 | before the natural end of the run, and require some special behaviour | |
219 | at that point. If you call such a destructor you will leave a pointer | |
220 | dangling inside the <code>Pythia</code> object you gave it in to, | |
221 | if that still exists. | |
222 | ||
223 | <method name="void SigmaProcess::initProc()"> | |
224 | is called once during initalization, and can then be used to set up | |
225 | parameters, such as masses and couplings, and perform calculations | |
226 | that need not be repeated for each new event, thereby saving time. | |
227 | This method needs not be implemented, since in principle all | |
228 | calculations can be done in <code>sigmaHat</code> below. | |
229 | </method> | |
230 | ||
231 | <method name="void SigmaProcess::sigmaKin()"> | |
232 | is called once a kinematical configuration has been determined, but | |
233 | before the two incoming flavours are known. This routine can therefore | |
234 | be used to perform calculations that otherwise might have to be repeated | |
235 | over and over again in <code>sigmaHat</code> below. For instance | |
236 | a flavour-independent cross section calculation for a <ei>q g</ei> | |
237 | initial state would be repeated 20 times in <code>sigmaHat</code>, | |
238 | five times for the five quark flavours allowed in the incoming beams, | |
239 | times twice to include antiquarks, times twice since the (anti)quark | |
240 | could be in either of the two beams. You could therefore calculate the | |
241 | result once only and store it as a private data member of the class. | |
242 | It is optional whether you want to use this method, however, or put | |
243 | everything in <code>sigmaHat</code>. | |
244 | </method> | |
245 | ||
246 | <method name="double SigmaProcess::sigmaHat()"> | |
247 | is the key method for cross section calculations and returns a cross section | |
248 | value, as described in the previous section. It is called when also a | |
249 | preliminary set of incoming flavours has been picked, in addition to the | |
250 | kinematical ones already available for <code>sigmaKin</code>. | |
251 | Typically <code>sigmaHat</code> is called inside a loop over all allowed | |
252 | incoming flavour combinations, stored in <code>id1</code> and | |
253 | <code>id2</code>, with fixed kinematics, as already illustrated above. | |
254 | The sum over the different flavour combinations provides the total | |
255 | cross section, while their relative size is used to make a selection of | |
256 | a specific incomimg state. | |
257 | </method> | |
258 | ||
259 | <method name="bool SigmaProcess::setupForME()"> | |
260 | to be called by the user from inside <code>sigmaHat()</code> | |
261 | (or possibly <code>sigmaKin()</code>) to setup alternative kinematics | |
262 | in the <code>mME[5]</code> and <code>pME[5]</code> arrays, better | |
263 | suited for matrix-element calculations. See the end of the previous | |
264 | section for a more detailed description. Should the method return | |
265 | <code>false</code> then the conversion did not work, and | |
266 | <code>sigmaHat()</code> (or <code>sigmaKin()</code>) should be set to | |
267 | vanish. | |
268 | </method> | |
269 | ||
270 | <method name="void SigmaProcess::setIdColAcol()"> | |
271 | is called only once an initial state and a kinematical configuration has | |
272 | been picked. This routine must set the complete flavour information and | |
273 | the colour flow of the process. This may involve further random choices, | |
274 | between different possible final-state flavours or between possible | |
275 | competing colour flows. Private data members of the class may be used to | |
276 | retain some information from the previous steps above. | |
277 | <br/>When this routine is called the two incoming flavours have already | |
278 | been selected and are available in <code>id1</code> and <code>id2</code>, | |
279 | whereas the one, two or three outgoing ones either are fixed for a given | |
280 | process or can be determined from the instate (e.g. whether a <ei>W^+</ei> | |
281 | or <ei>W^-</ei> was produced). There is also a standard method in | |
282 | <code>VCKM</code> to pick a final flavour from an initial one with CKM | |
283 | mixing. Once you have figured out the value of | |
284 | <code>id3</code> and, the case being, <code>id4</code> and | |
285 | <code>id5</code>, you store these values permanently by a call | |
286 | <code>setId( id1, id2, id3, id4, id5)</code>, where the last two may be | |
287 | omitted if irrelevant. | |
288 | <br/>Correspondingly, the colours are stored with | |
289 | <code>setColAcol( col1, acol1, col2, acol2, col3, acol3, col4, acol4, | |
290 | col5, acol5)</code>, where the final ones may be omitted if irrelevant. | |
291 | Les Houches style colour tags are used, but starting with number 1 | |
292 | (and later shifted by the currently requested offset). The | |
293 | input is grouped particle by particle, with the colour index before the | |
294 | anticolour one. You may need to select colour flow dynamically, depending | |
295 | on the kinematics, when several distinct possibilities exist. Trivial | |
296 | operations, like swapping colours and anticolours, can be done with | |
297 | existing methods. | |
298 | <br/>When the <code>id3Mass()</code> and <code>id4Mass()</code> | |
299 | methods have been used, the order of the outgoing particles may be | |
300 | inconsistent with the way the <ei>tHat</ei> and <ei>uHat</ei> | |
301 | variables have been defined. A typical example would be a process like | |
302 | <ei>q g -> q' W</ei> with <ei>tHat</ei> defined between incoming and | |
303 | outgoing quark, but where <code>id3Mass() = 24</code> and so the | |
304 | process is to be stored as <ei>q g -> W q'</ei>. One should then put | |
305 | the variable <code>swapTU = true</code> in <code>setIdColAcol()</code> | |
306 | for each event where the <ei>tHat</ei> and <ei>uHat</ei> variables | |
307 | should be swapped before the event kinematics is reconstructed. This | |
308 | variable is automatically restored to <code>false</code> for each new | |
309 | event. | |
310 | </method> | |
311 | ||
312 | <method name="double SigmaProcess::weightDecayFlav( Event& process)"> | |
313 | is called to allow a reweighting of the simultaneous flavour choices of | |
314 | resonance decay products. Is currently only used for the | |
315 | <ei>q qbar -> gamma*/Z^0 gamma*/Z^0</ei> process, and will likely not | |
316 | be of interest for you. | |
317 | </method> | |
318 | ||
319 | <method name="double SigmaProcess::weightDecay( Event& process, | |
320 | int iResBeg, int iResEnd)"> | |
321 | is called when the basic process has one or several resonances, after each | |
322 | set of related resonances in <code>process[i]</code>, | |
323 | <code>iResBeg</code> <= <code>i </code> <= <code>iResEnd</code>, | |
324 | has been allowed to decay. The calculated weight, to be normalized | |
325 | to the range between 0 and 1, is used to decide whether to accept the | |
326 | decay(s) or try for a new decay configuration. The base-class version of | |
327 | this method returns unity, i.e. gives isotropic decays by default. | |
328 | This method may be called repeatedly for a single event. For instance, in | |
329 | <ei>q qbar -> H^0 Z^0</ei> with <ei>H^0 -> W^+ W^-</ei>, a first call | |
330 | would be made after the <ei>H^0</ei> and <ei>Z^0</ei> decays, and then | |
331 | depend only on the <ei>Z^0</ei> decay angles since the <ei>H^0</ei> | |
332 | decays isotropically. The second call would be after the <ei>W^+ W^-</ei> | |
333 | decays and then involve correlations between the four daughter fermions. | |
334 | </method> | |
335 | ||
336 | <method name="string SigmaProcess::name()"> | |
337 | returns the name of the process, as you want it to be shown in listings. | |
338 | </method> | |
339 | ||
340 | <method name="int SigmaProcess::code()"> | |
341 | returns an integer identifier of the process. This has no internal function, | |
342 | but is only intended as a service for the user to rapidly (and hopefully | |
343 | uniquely) identify which process occured in a given event. Numbers below | |
344 | 10000 are reserved for internal PYTHIA use. | |
345 | </method> | |
346 | ||
347 | <method name="string SigmaProcess::inFlux()"> | |
348 | this string specifies the combinations of incoming partons that are | |
349 | allowed for the process under consideration, and thereby which incoming | |
350 | flavours <code>id1</code> and <code>id2</code> the <code>sigmaHat()</code> | |
351 | calls will be looped over. It is always possible to pick a wider flavour | |
352 | selection than strictly required and then put to zero cross sections in | |
353 | the superfluous channels, but of course this may cost some extra execution | |
354 | time. Currently allowed options are: | |
355 | <br/>* <code>gg</code>: two gluons. | |
356 | <br/>* <code>qg</code>: one (anti)quark and one gluon. | |
357 | <br/>* <code>qq</code>: any combination of two quarks, two antiquarks or | |
358 | a quark and an antiquark. | |
359 | <br/>* <code>qqbarSame</code>: a quark and its antiquark; | |
360 | this is a subset of the above <code>qq</code> option. | |
361 | <br/>* <code>ff</code>: any combination of two fermions, two antifermions | |
362 | or a fermion and an antifermion; is the same as <code>qq</code> for | |
363 | hadron beams but also allows processes to work with lepton beams. | |
364 | <br/>* <code>ffbarSame</code>: a fermion and its antifermion; is the | |
365 | same as <code>qqbarSame</code> for hadron beams but also allows processes | |
366 | to work with lepton beams. | |
367 | <br/>* <code>ffbarChg</code>: a fermion and an antifermion that combine | |
368 | to give charge +-1. | |
369 | <br/>* <code>fgm</code>: a fermion and a photon (gamma). | |
370 | <br/>* <code>ggm</code>: a gluon and a photon. | |
371 | <br/>* <code>gmgm</code>: two photons. | |
372 | </method> | |
373 | ||
374 | <method name="bool SigmaProcess::convert2mb()"> | |
375 | it is assumed that cross sections normally come in dimensions such that | |
376 | they, when integrated over the relevant phase space, obtain the dimension | |
377 | GeV^-2, and therefore need to be converted to mb. If the cross section | |
378 | is already encoded as mb then <code>convert2mb()</code> should be | |
379 | overloaded to instead return <code>false</code>. | |
380 | </method> | |
381 | ||
382 | <method name="bool SigmaProcess::convertME()"> | |
383 | it is assumed that <ei>2 -> 1</ei> cross sections are encoded as | |
384 | <ei>sigmaHat(sHat)</ei>, and <ei>2 -> 2</ei> ones as | |
385 | <ei>d(sigmaHat)/d(tHat)</ei> in the <code>SigmaProcess::sigmaHat()</code> | |
386 | methods. If <code>convertME()</code> is overloaded to instead return | |
387 | <code>true</code> then the return value is instead assumed to be the | |
388 | squared matrix element <ei>|M|^2</ei>, and | |
389 | <code>SigmaProcess::sigmaHatWrap(...)</code> converts to | |
390 | <ei>sigmaHat(sHat)</ei> or <ei>d(sigmaHat)/d(tHat)</ei>, respectively. | |
391 | This switch has no effect on <ei>2 -> 3</ei> processes, where | |
392 | <ei>|M|^2</ei> is the only allowed input anyway. | |
393 | </method> | |
394 | ||
395 | <method name="int SigmaProcess::id3Mass()"> | |
396 | </method> | |
397 | <methodmore name="int SigmaProcess::id4Mass()"> | |
398 | </methodmore> | |
399 | <methodmore name="int SigmaProcess::id5Mass()"> | |
400 | are the one, two or three final-state flavours, where masses are to be | |
401 | selected before the matrix elements are evaluated. Only the absolute value | |
402 | should be given. For massless particles, like gluons and photons, one need | |
403 | not give anything, i.e. one defaults to 0. The same goes for normal light | |
404 | quarks, where masses presumably are not implemented in the matrix elements. | |
405 | Later on, these quarks can still (automatically) obtain constituent masses, | |
406 | once a <ei>u</ei>, <ei>d</ei> or <ei>s</ei> flavour has been selected. | |
407 | </methodmore> | |
408 | ||
409 | <method name="int SigmaProcess::resonanceA()"> | |
410 | </method> | |
411 | <methodmore name="int SigmaProcess::resonanceB()"> | |
412 | are the codes of up to two <ei>s</ei>-channel resonances contributing to | |
413 | the matrix elements. These are used by the program to improve the phase-space | |
414 | selection efficiency, by partly sampling according to the relevant | |
415 | Breit-Wigners. Massless resonances (the gluon and photon) need not be | |
416 | specified. | |
417 | </methodmore> | |
418 | ||
419 | <method name="bool SigmaProcess::isSChannel()"> | |
420 | normally the choice of renormalization and factorization scales in | |
421 | <ei>2 -> 2</ei> and <ei>2 -> 3</ei> processes is based on the assumption | |
422 | that <ei>t</ei>- and <ei>u</ei>-channel exchanges dominates the | |
423 | cross section. In cases such as <ei>f fbar -> gamma* -> f' fbar'</ei> a | |
424 | <ei>2 -> 2</ei> process actually ought to be given scales as a | |
425 | <ei>2 -> 1</ei> one, in the sense that it proceeds entirely through | |
426 | an <ei>s</ei>-channel resonance. This can be achieved if you override the | |
427 | default <code>false</code> to return <code>true</code>. See further the | |
428 | page on <aloc href="CouplingsAndScales">couplings and scales</aloc>. | |
429 | </method> | |
430 | ||
431 | <method name="int SigmaProcess::idSChannel()"> | |
432 | normally no intermediate state is shown in the event record for | |
433 | <ei>2 -> 2</ei> and <ei>2 -> 3</ei> processes. However, in case | |
434 | that <code>idSChannel</code> is overloaded to return a nonzero value, | |
435 | an intermediate particle with that identity code is inserted into the | |
436 | event record, to make it a <ei>2 -> 1 -> 2</ei> or <ei>2 -> 1 -> 3</ei> | |
437 | process. Thus if both <code>isSChannel</code> and <code>idSChannel</code> | |
438 | are overloaded, a process will behave and look like it proceeded through | |
439 | a resonance. The one difference is that the implementation of the | |
440 | matrix element is not based on the division into a production and a | |
441 | decay of an intermediate resonance, but is directly describing the | |
442 | transition from the initial to the final state. | |
443 | </method> | |
444 | ||
445 | <method name="int SigmaProcess::isQCD3body()"> | |
446 | there are two different 3-body phase-space selection machineries, | |
447 | of which the non-QCD one is default. If you overload this method | |
448 | instead the QCD-inspired machinery will be used. The differences | |
449 | between these two is related to which | |
450 | <aloc href="PhaseSpaceCuts">phase space cuts</aloc> | |
451 | can be set, and also that the QCD machinery assumes (almost) massless | |
452 | outgoing partons. | |
453 | </method> | |
454 | ||
455 | <method name="int SigmaProcess::idTchan1()"> | |
456 | </method> | |
457 | <methodmore name="int SigmaProcess::idTchan2()"> | |
458 | the non-QCD <ei>2 -> 3</ei> phase space selection machinery is rather | |
459 | primitive, as already mentioned. The efficiency can be improved in | |
460 | processes that proceed though <ei>t</ei>-channel exchanges, such as | |
461 | <ei>q qbar' -> H^0 q qbar'</ei> via <ei>Z^0 Z^0</ei> fusion, if the identity | |
462 | of the <ei>t</ei>-channel-exchanged particles on the two side of the | |
463 | event are provided. Only the absolute value is of interest. | |
464 | </methodmore> | |
465 | ||
466 | <method name="double SigmaProcess::tChanFracPow1()"> | |
467 | </method> | |
468 | <methodmore name="double SigmaProcess::tChanFracPow2()"> | |
469 | in the above kind of <ei>2 -> 3</ei> phase-space selection, the | |
470 | sampling of <ei>pT^2</ei> is done with one part flat, one part weighted | |
471 | like <ei>1 / (pT^2 + m_R^2)</ei> and one part like | |
472 | <ei>1 / (pT^2 + m_R^2)^2</ei>. The above values provide the relative | |
473 | amount put in the latter two channels, respectively, with the first | |
474 | obtaining the rest. Thus the sum of <code>tChanFracPow1()</code> and | |
475 | <code>tChanFracPow2()</code> must be below unity. The final results | |
476 | should be independent of these numbers, but the Monte Carlo efficiency | |
477 | may be quite low for a bad choice. Here <ei>m_R</ei> is the mass of the | |
478 | exchanged resonance specified by <code>idTchan1()</code> or | |
479 | <code>idTchan2()</code>. Note that the order of the final-state | |
480 | listing is important in the above <ei>q qbar' -> H^0 q qbar'</ei> example, | |
481 | i.e. the <ei>H^0</ei> must be returned by <code>id3Mass()</code>, | |
482 | since it is actually the <ei>pT^2</ei> of the latter two that are | |
483 | selected independently, with the first <ei>pT</ei> then fixed | |
484 | by transverse-momentum conservation. | |
485 | </methodmore> | |
486 | ||
487 | <method name="bool SigmaProcess::useMirrorWeight()"> | |
488 | in <ei>2 -> 3</ei> processes the phase space selection used here | |
489 | involves a twofold ambiguity basically corresponding to a flipping of | |
490 | the positions of last two outgoing particles. These are assumed equally | |
491 | likely by default, <code>false</code>, but for processes proceeding entirely | |
492 | through <ei>t</ei>-channel exchange the Monte Carlo efficiency can be | |
493 | improved by making a preselection based on the relative propagator | |
494 | weights, <code>true</code>. | |
495 | </method> | |
496 | ||
497 | <method name="int SigmaProcess::gmZmode()"> | |
498 | allows a possibility to override the global mode | |
499 | <code><aloc href="ElectroweakProcesses">WeakZ0:gmZmode</aloc></code> | |
500 | for a specific process. The global mode normally is used to switch off | |
501 | parts of the <ei>gamma^*/Z^0</ei> propagator for test purposes. The | |
502 | above local mode is useful for processes where a <ei>Z^0</ei> really is | |
503 | that and nothing more, such as <ei>q qbar -> H^0 Z^0</ei>. The default | |
504 | value -1 returned by <code>gmZmode()</code> ensures that the global | |
505 | mode is used, while 0 gives full <ei>gamma^*/Z^0</ei> interference, | |
506 | 1 <ei>gamma^*</ei> only and 2 <ei>Z^0</ei> only. | |
507 | </method> | |
508 | ||
509 | <h3>Access to a process</h3> | |
510 | ||
511 | Once you have implemented a class, it is straightforward to make use of | |
512 | it in a run. Assume you have written a new class <code>MySigma</code>, | |
513 | which inherits from <code>Sigma1Process</code>, <code>Sigma2Process</code> | |
514 | or <code>Sigma3Process</code>, which in their turn inherit from | |
515 | <code>SigmaProcess</code>. You then create an instance of this class | |
516 | and hand it in to a <code>pythia</code> object with | |
517 | <pre> | |
518 | SigmaProcess* mySigma = new MySigma(); | |
519 | pythia.setSigmaPtr( mySigma); | |
520 | </pre> | |
521 | If you have several processes you can repeat the procedure any number | |
522 | of times. When <code>pythia.init(...)</code> is called these processes | |
523 | are initialized along with any internal processes you may have switched on, | |
524 | and treated in exactly the same manner. The <code>pythia.next()</code> | |
525 | will therefore generate a mix of the different kinds of processes without | |
526 | distinction. See also the <aloc href="ProgramFlow">Program Flow</aloc> | |
527 | description. | |
528 | ||
529 | <p/> | |
530 | If the code should be of good quality and general usefulness, it would | |
531 | be simple to include it as a permanently available process in the | |
532 | standard program distribution. The final step of that integration ought to | |
533 | be left for the PYTHIA authors, but here is a description of what is | |
534 | required. | |
535 | ||
536 | <p/> | |
537 | A flag has to be defined, that allows the process to be switched on; | |
538 | by default it should always be off. The name of the flag should be | |
539 | chosen of the type <code>model:process</code>. Here the | |
540 | <code>model</code> would be related to the general scenario considered, | |
541 | e.g. <code>Compositeness</code>, while <code>process</code> would | |
542 | specify instate and outstate, separated by a 2 (= to), e.g. | |
543 | <code>ug2u*g</code>. | |
544 | When several processes are implemented and "belong together" it is | |
545 | also useful to define a <code>model:all</code> switch that affects | |
546 | all the separate processes. | |
547 | ||
548 | <p/> | |
549 | The flags should normally be stored in the <code>ProcessSelection.xml</code> | |
550 | file or one of its daughters for a specific kind of processes. This is to | |
551 | make them easily found by users. You could create and use your own | |
552 | <code>.xml</code> file, so long as you then add that name to the | |
553 | list of files in the <code>Index.xml</code> file. (If not, | |
554 | the flags would never be created and the program would not work.) | |
555 | ||
556 | <p/> | |
557 | In the <code>ProcessContainer.c</code> file, the | |
558 | <code>SetupContainers::init()</code> method needs to be expanded to | |
559 | create instances of the processes switched on. This code is fairly | |
560 | repetitive, and should be easy to copy and modify from the code | |
561 | already there. The basic structure is | |
562 | <br/>(i) check whether a process is requested by the user and, if so, | |
563 | <br/>(ii) create an instance of the matrix-element class, | |
564 | <br/>(iii)create a container for the matrix element and its associated | |
565 | phase-space handling, and | |
566 | <br>(iv) add the container to the existing process list. | |
567 | ||
568 | <p/> | |
569 | Two minor variations are possible. One is that a set of related | |
570 | processes are lumped inside the the same initial check, i.e. are | |
571 | switched on all together. The second is that the matrix-element | |
572 | constructor may take arguments, as specified by you (see above). | |
573 | If so, the same basic matrix element may be recycled for a set of | |
574 | related processes, e.g. one for a composite <ei>u</ei> and one for | |
575 | a composite <ei>d</ei>. Obviously these variations may be combined. | |
576 | ||
577 | </chapter> | |
578 | ||
579 | <!-- Copyright (C) 2010 Torbjorn Sjostrand --> |