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