[1] | 1 | <html> |
---|
| 2 | <head> |
---|
| 3 | <title>Resonance Decays</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 | <script language=javascript type=text/javascript> |
---|
| 10 | function stopRKey(evt) { |
---|
| 11 | var evt = (evt) ? evt : ((event) ? event : null); |
---|
| 12 | var node = (evt.target) ? evt.target :((evt.srcElement) ? evt.srcElement : null); |
---|
| 13 | if ((evt.keyCode == 13) && (node.type=="text")) |
---|
| 14 | {return false;} |
---|
| 15 | } |
---|
| 16 | |
---|
| 17 | document.onkeypress = stopRKey; |
---|
| 18 | </script> |
---|
| 19 | <?php |
---|
| 20 | if($_POST['saved'] == 1) { |
---|
| 21 | if($_POST['filepath'] != "files/") { |
---|
| 22 | echo "<font color='red'>SETTINGS SAVED TO FILE</font><br/><br/>"; } |
---|
| 23 | else { |
---|
| 24 | echo "<font color='red'>NO FILE SELECTED YET.. PLEASE DO SO </font><a href='SaveSettings.php'>HERE</a><br/><br/>"; } |
---|
| 25 | } |
---|
| 26 | ?> |
---|
| 27 | |
---|
| 28 | <form method='post' action='ResonanceDecays.php'> |
---|
| 29 | |
---|
| 30 | <h2>Resonance Decays</h2> |
---|
| 31 | |
---|
| 32 | The <code>ResonanceDecays</code> class performs the sequential decays of |
---|
| 33 | all resonances formed in the hard process. Note the important distinction |
---|
| 34 | between "resonances" and other "particles" made in PYTHIA. |
---|
| 35 | <ul> |
---|
| 36 | <li> |
---|
| 37 | The list of resonances contains <i>gamma^*/Z^0</i>, <i>W^+-</i>, top, |
---|
| 38 | the Higgs, and essentially all new particles of Beyond-the-Standard-Model |
---|
| 39 | physics: further Higgses, sfermions, gauginos, techniparticles, and so on. |
---|
| 40 | The partial widths to different decay channels are perturbatively |
---|
| 41 | calculable, given the parameters of the respective model, and branching |
---|
| 42 | ratios may be allowed to vary across a (reasonably broad) resonance peak. |
---|
| 43 | Usually resonances are short-lived, and therefore it makes sense to consider |
---|
| 44 | their decays immediately after the primary hard process has been set up. |
---|
| 45 | Furthermore, in several cases the decay angular distributions are encoded |
---|
| 46 | as part of the specific process, e.g. the <i>W</i> decays differently in |
---|
| 47 | <i>f fbar -> W^+-</i>, <i>f fbar -> W^+ W^-</i> and |
---|
| 48 | <i>h^0 -> W^+ W^- </i>. All of these particles are (in PYTHIA) only |
---|
| 49 | produced as part of the hard process itself, i.e. they are not produced |
---|
| 50 | in showers or hadronization processes. Therefore the restriction to |
---|
| 51 | specific decay channels can be consistently taken into account as a |
---|
| 52 | corresponding reduction in the cross section of a process. Finally, note |
---|
| 53 | that all of these resonances have an on-shell mass above 20 GeV, with the |
---|
| 54 | exception of some hypothetical weakly interacting and stable particles |
---|
| 55 | such as the gravitino. |
---|
| 56 | </li> |
---|
| 57 | <li> |
---|
| 58 | The other particles include normal hadrons and the Standard-Model leptons, |
---|
| 59 | including the <i>tau^+-</i>. These can be produced in the normal |
---|
| 60 | hadronization and decay description, which involve unknown nonperturbative |
---|
| 61 | parameters and multistep chains that cannot be predicted beforehand: |
---|
| 62 | a hard process like <i>g g -> g g</i> can develop a shower with a |
---|
| 63 | <i>g -> b bbar</i> branching, where the <i>b</i> hadronizes to a |
---|
| 64 | <i>B^0bar</i> that oscillates to a <i>B^0</i> that decays to a |
---|
| 65 | <i>tau^+</i>. Therefore any change of branching ratios - most of which |
---|
| 66 | are determined from data rather than from first principles anyway - |
---|
| 67 | will not be taken into account in the cross section of a process. |
---|
| 68 | Exceptions exist, but most particles in this class are made to decay |
---|
| 69 | isotropically. Finally, note that all of these particles have a mass |
---|
| 70 | below 20 GeV. |
---|
| 71 | </li> |
---|
| 72 | </ul> |
---|
| 73 | |
---|
| 74 | There is one ambiguous case in this classification, namely the photon. |
---|
| 75 | The <i>gamma^*/Z^0</i> combination contains a low-mass peak when |
---|
| 76 | produced in a hard process. On the other hand, photons can participate |
---|
| 77 | in shower evolution, and therefore a photon originally assumed |
---|
| 78 | massless can be assigned an arbitrarily high mass when it is allowed |
---|
| 79 | to branch into a fermion pair. In some cases this could lead to |
---|
| 80 | doublecounting, e.g. between processes such as |
---|
| 81 | <i>f fbar -> (gamma^*/Z^0) (gamma^*/Z^0)</i>, |
---|
| 82 | <i>f fbar -> (gamma^*/Z^0) gamma</i> and |
---|
| 83 | <i>f fbar -> gamma gamma</i>. Here it make sense to limit the |
---|
| 84 | lower mass allowed for the <i>gamma^*/Z^0</i> combination, |
---|
| 85 | in <code>23:mMin</code>, to be the same as the upper limit allowed |
---|
| 86 | for an off-shell photon in the shower evolution, in |
---|
| 87 | <code>TimeShower:mMaxGamma</code>. By default this matching is done |
---|
| 88 | at 10 GeV. |
---|
| 89 | |
---|
| 90 | <p/> |
---|
| 91 | In spite of the above-mentioned differences, the resonances and the |
---|
| 92 | other particles are all stored in one common |
---|
| 93 | <?php $filepath = $_GET["filepath"]; |
---|
| 94 | echo "<a href='ParticleData.php?filepath=".$filepath."' target='page'>";?>particle data table</a>, so as to offer a |
---|
| 95 | uniform interface to <?php $filepath = $_GET["filepath"]; |
---|
| 96 | echo "<a href='ParticleDataScheme.php?filepath=".$filepath."' target='page'>";?>setting and |
---|
| 97 | getting</a> properties such as name, mass, charge and decay modes, |
---|
| 98 | also for the <?php $filepath = $_GET["filepath"]; |
---|
| 99 | echo "<a href='ParticleProperties.php?filepath=".$filepath."' target='page'>";?>particle properties</a> |
---|
| 100 | in the event record. Some methods are specific to resonances, however, |
---|
| 101 | in particular for the calculation of partial widths and thereby of |
---|
| 102 | branching ratio. For resonances these can be calculated dynamically, |
---|
| 103 | set up at initialization for the nominal mass and then updated to the |
---|
| 104 | current mass when these are picked according to a Breit-Wigner resonance |
---|
| 105 | shape. |
---|
| 106 | |
---|
| 107 | <h3>Resonance Decays and Cross Sections</h3> |
---|
| 108 | |
---|
| 109 | As already hinted above, you have the possibility to set the allowed |
---|
| 110 | decay channels of resonances, see |
---|
| 111 | <?php $filepath = $_GET["filepath"]; |
---|
| 112 | echo "<a href='ParticleDataScheme.php?filepath=".$filepath."' target='page'>";?>Particle Data Scheme</a> description. |
---|
| 113 | For instance, if you study the process <i>q qbar -> H^0 Z^0</i> |
---|
| 114 | you could specify that the <i>Z^0</i> should decay only to |
---|
| 115 | lepton pairs, the <i>H^0</i> only to <i>W^+ W^-</i>, the |
---|
| 116 | <i>W^+</i> only to a muon and a neutrino, while the <i>W^-</i> |
---|
| 117 | can decay to anything. Unfortunately there are limits to the |
---|
| 118 | flexibility: you cannot set a resonance to have different properties |
---|
| 119 | in different places of a process, e.g. if instead |
---|
| 120 | <i>H^0 -> Z^0 Z^0</i> in the above process then the three |
---|
| 121 | <i>Z^0</i>'s would all obey the same rules. |
---|
| 122 | |
---|
| 123 | <p/> |
---|
| 124 | The restrictions on the allowed final states of a process is directly |
---|
| 125 | reflected in the cross section of it. That is, if some final states |
---|
| 126 | are excluded then the cross section is reduced accordingly. Such |
---|
| 127 | restrictions are built up recursively in cases of sequential decay |
---|
| 128 | chains. The restrictions are also reflected in the compositions of |
---|
| 129 | those events that actually do get to be generated. For instance, |
---|
| 130 | the relative rates of <i>H^0 -> W^+ W^-</i> and |
---|
| 131 | <i>H^0 -> Z^0 Z^0</i> are shifted when the allowed sets of |
---|
| 132 | <i>W^+-</i> and <i>Z^0</i> decay channels are changed. |
---|
| 133 | |
---|
| 134 | <p/> |
---|
| 135 | We remind that only those particles that Pythia treat as resonances |
---|
| 136 | enjoy this property, and only those that are considered as part of the |
---|
| 137 | hard process and its assocaited resonance decays. |
---|
| 138 | |
---|
| 139 | <p/> |
---|
| 140 | There is one key restriction on resonances: |
---|
| 141 | <br/><br/><table><tr><td><strong>ResonanceWidths:minWidth </td><td></td><td> <input type="text" name="1" value="1e-20" size="20"/> (<code>default = <strong>1e-20</strong></code>; <code>minimum = 1e-30</code>)</td></tr></table> |
---|
| 142 | Minimal allowed width of a resonance, in GeV. If the width falls below |
---|
| 143 | this number the resonance is considered stable and will not be allowed |
---|
| 144 | to decay. This is mainly intended as a technical parameter, to avoid |
---|
| 145 | disasters in cases where no open decay channels exists at all. It could |
---|
| 146 | be used for real-life decisions as well, however, but then typically |
---|
| 147 | would have to be much bigger than the default value. Special caution |
---|
| 148 | would be needed if coloured resonance particles were made stable, since |
---|
| 149 | the program would not necessarily know how to hadronize them, and |
---|
| 150 | therefore fail at that stage. |
---|
| 151 | |
---|
| 152 | |
---|
| 153 | <p/> |
---|
| 154 | In spite of this technical parameter choice, it is possible to set |
---|
| 155 | a lifetime for a resonance, and thereby to obtain displaced vertices. |
---|
| 156 | If a resonance is allowed to decay it will do so, irrespective of |
---|
| 157 | the location of the decay vertex. This is unlike |
---|
| 158 | <?php $filepath = $_GET["filepath"]; |
---|
| 159 | echo "<a href='ParticleDecays.php?filepath=".$filepath."' target='page'>";?>normal particle decays</a>, |
---|
| 160 | where it is possible to define some region around the primary |
---|
| 161 | vertex within which all decays should happen, with particles |
---|
| 162 | leaving that region considered stable. The logic is that resonances |
---|
| 163 | as a rule are too short-lived for secondary vertices, |
---|
| 164 | so if you pick a scenario with a long-lived but unstable resonance |
---|
| 165 | it is because you <i>want</i> to study secondary vertices. |
---|
| 166 | How to interface those decays to a detector simulation program then |
---|
| 167 | is another story, to be solved separately. Do note that a special |
---|
| 168 | treatment is needed for coloured long-lived resonances, that form |
---|
| 169 | <?php $filepath = $_GET["filepath"]; |
---|
| 170 | echo "<a href='Rhadrons.php?filepath=".$filepath."' target='page'>";?>R-hadrons</a>, and where charge and flavour |
---|
| 171 | may change between the production and decay vertices. |
---|
| 172 | |
---|
| 173 | <h3>Special properties and methods for resonances</h3> |
---|
| 174 | |
---|
| 175 | The method <code>ParticleData::isResonance(id)</code> allows you to |
---|
| 176 | query whether a given particle species is considered a resonance or not. |
---|
| 177 | You can also change the default value of this flag in the normal way, |
---|
| 178 | e.g. <code>pythia.readString("id:isResonance = true")</code>. |
---|
| 179 | |
---|
| 180 | <p/> |
---|
| 181 | An option with a forced width can be set with the |
---|
| 182 | <code>id:doForceWidth</code> flag as above, and queried with |
---|
| 183 | <code>ParticleData::doForceWidth(id)</code>. It is by default |
---|
| 184 | <code>off</code>, and should normally so remain. If switched |
---|
| 185 | <code>on</code> then the width stored in <code>id:mWidth</code> is |
---|
| 186 | strictly used to describe the Breit-Wigner of the resonance. This is |
---|
| 187 | unlike the normal behaviour of standard resonances such as the |
---|
| 188 | <i>Z^0</i>, <i>W^+-</i>, <i>t</i> or <i>h^0</i>, which have |
---|
| 189 | explicit decay-widths formulae encoded, in classes derived from the |
---|
| 190 | <code><?php $filepath = $_GET["filepath"]; |
---|
| 191 | echo "<a href='SemiInternalResonances.php?filepath=".$filepath."' target='page'>";?>ResonanceWidths</a></code> |
---|
| 192 | base class. These formulae are used, e.g., to derive all the Higgs partial |
---|
| 193 | widths as a function of the Higgs mass you choose, and at initialization |
---|
| 194 | overwrites the existing total width value. The reason for forcing the |
---|
| 195 | width to another value specified by you would normally more have to do |
---|
| 196 | with experimental issues than with physics ones, e.g. how sensitive your |
---|
| 197 | detector would be to changes in the Higgs width by a factor of two. |
---|
| 198 | A warning is that such a rescaling could modify the cross section of |
---|
| 199 | a process correspondingly for some processes, while leaving it |
---|
| 200 | (essentially) unchanged for others (as would seem most logical), |
---|
| 201 | depending on how these were encoded. A further warning is that, |
---|
| 202 | if you use this facility for <i>Z^0</i> or <i>Z'^0</i> with |
---|
| 203 | <i>gamma^*/Z^0</i> or <i>gamma^*/Z^0/Z'^0</i> interference on, |
---|
| 204 | then also the handling of this interference is questionable. |
---|
| 205 | So, if you need to use the width-rescaling option, be extremely cautios. |
---|
| 206 | |
---|
| 207 | <p/> |
---|
| 208 | If a resonance does not have a class of its own, with hardcoded equations |
---|
| 209 | for all relevant partial widths, then a simpler object will be created |
---|
| 210 | at initialization. This object will take the total width and branching |
---|
| 211 | ratios as is (with the optional variations explained in the next section), |
---|
| 212 | and thus the rescaling approach brings no further freedom. |
---|
| 213 | |
---|
| 214 | <p/> |
---|
| 215 | Mainly for internal usage, the |
---|
| 216 | <code><?php $filepath = $_GET["filepath"]; |
---|
| 217 | echo "<a href='ParticleDataScheme.php?filepath=".$filepath."' target='page'>";?>ParticleData</a></code> contain |
---|
| 218 | some special methods that are only meaningful for resonances: |
---|
| 219 | <ul> |
---|
| 220 | <li><code>resInit(...)</code> to initialize a resonance, possibly |
---|
| 221 | including a recalculation of the nominal width to match the nominal |
---|
| 222 | mass;</li> |
---|
| 223 | <li><code>resWidth(...)</code> to calculate the partial and total widths |
---|
| 224 | at the currently selected mass;</li> |
---|
| 225 | <li><code>resWidthOpen(...)</code> to calculate the partial and total |
---|
| 226 | widths of those channels left open by user switches, at the currently |
---|
| 227 | selected mass;</li> |
---|
| 228 | <li><code>resWidthStore(...)</code> to calculate the partial and total |
---|
| 229 | widths of those channels left open by user switches, at the currently |
---|
| 230 | selected mass, and store those as input for a subsequent selection of |
---|
| 231 | decay channel;</li> |
---|
| 232 | <li><code>resOpenFrac(...)</code> to return the fraction of the total |
---|
| 233 | width that is open by the decay channel selection made by users (based on |
---|
| 234 | the choice of <code><?php $filepath = $_GET["filepath"]; |
---|
| 235 | echo "<a href='ParticleDataScheme.php?filepath=".$filepath."' target='page'>";?>onMode</a></code> |
---|
| 236 | for the various decay channels, recursively calculated for sequential |
---|
| 237 | decays);</li> |
---|
| 238 | <li><code>resWidthRescaleFactor(...)</code> returns the factor by which |
---|
| 239 | the internally calculated PYTHIA width has to be rescaled to give the |
---|
| 240 | user-enforced width;</li> |
---|
| 241 | <li><code>resWidthChan(...)</code> to return the width for one particular |
---|
| 242 | channel (currently only used for Higgs decays, to obtain instate coupling |
---|
| 243 | from outstate width).</li> |
---|
| 244 | </ul> |
---|
| 245 | These methods actually provide an interface to the classes derived from |
---|
| 246 | the <code>ResonanceWidths</code> base class, to describe various |
---|
| 247 | resonances. |
---|
| 248 | |
---|
| 249 | <h3>Modes for Matrix Element Processing</h3> |
---|
| 250 | |
---|
| 251 | The <code>meMode()</code> value for a decay mode is used to specify |
---|
| 252 | <?php $filepath = $_GET["filepath"]; |
---|
| 253 | echo "<a href='ParticleDecays.php?filepath=".$filepath."' target='page'>";?>nonisotropic decays or the conversion of |
---|
| 254 | a parton list into a set of hadrons</a> in some channels of normal |
---|
| 255 | particles. For resonances it can also take a third function, namely |
---|
| 256 | to describe how the branching ratios and widths of a resonance should |
---|
| 257 | be rescaled as a function of the current mass of the decaying resonance. |
---|
| 258 | The rules are especially useful when new channels are added to an |
---|
| 259 | existing particle, or a completely new resonance added. |
---|
| 260 | |
---|
| 261 | <ul> |
---|
| 262 | <li>0 : channels for which hardcoded partial-width expressions are |
---|
| 263 | expected to exist in the derived class of the respective resonance. |
---|
| 264 | Should no such code exist then the partial width defaults to zero. |
---|
| 265 | </li> |
---|
| 266 | <li>1 - 99 : same as 0, but normally not used for resonances.</li> |
---|
| 267 | <li>100 : calculate the partial width of the channel from its stored |
---|
| 268 | branching ratio times the stored total width. This value remains unchanged |
---|
| 269 | when the resonance fluctuates in mass. Specifically there are no |
---|
| 270 | threshold corrections. That is, if the resonance fluctuates down in |
---|
| 271 | mass, to below the nominal threshold, it is assumed that one of the |
---|
| 272 | daughters could also fluctuate down to keep the channel open. (If not, |
---|
| 273 | there may be problems later on.) |
---|
| 274 | </li> |
---|
| 275 | <li>101 : calculate the partial width of the channel from its stored |
---|
| 276 | branching ratio times the stored total width. Multiply by a step threshold, |
---|
| 277 | i.e. the channel is switched off when the sum of the daughter on-shell |
---|
| 278 | masses is above the current mother mass.</li> |
---|
| 279 | <li>102 : calculate the partial width of the channel from its stored |
---|
| 280 | branching ratio times the stored total width. Multiply by a smooth |
---|
| 281 | threshold factor |
---|
| 282 | <i>beta = sqrt( (1 - m_1^2/m_2 - m_2^2/m^2)^2 - 4 m_1^2 m_2^2/m^4)</i> |
---|
| 283 | for two-body decays and <i>sqrt(1 - Sum_i m_i / m)</i> for multibody |
---|
| 284 | ones. The former correctly encodes the size of the phase space but |
---|
| 285 | misses out on any nontrivial matrix-element behaviour, while the latter |
---|
| 286 | obviously is a very crude simplification of the correct phase-space |
---|
| 287 | expression. Specifically, it is thereby assumed that the stored branching |
---|
| 288 | ratio and total width did not take into account such a factor.</li> |
---|
| 289 | <li>103 : use the same kind of behaviour and threshold factor as for |
---|
| 290 | 102 above, but assume that such a threshold factor has been used when |
---|
| 291 | the default branching ratio and total width were calculated, so that one |
---|
| 292 | should additionally divide by the on-shell threshold factor. Specifically, |
---|
| 293 | this will give back the stored branching ratios for on-shell mass, |
---|
| 294 | unlike the 102 option. To avoid division by zero, or in general |
---|
| 295 | unreasonably big rescaling factors, a lower limit |
---|
| 296 | <code>minThreshold</code> (see below) on the value of the on-shell |
---|
| 297 | threshold factor is imposed. (In cases where a big rescaling is |
---|
| 298 | intentional, code 102 would be more appropriate.) </li> |
---|
| 299 | </ul> |
---|
| 300 | |
---|
| 301 | <br/><br/><table><tr><td><strong>ResonanceWidths:minThreshold </td><td></td><td> <input type="text" name="2" value="0.1" size="20"/> (<code>default = <strong>0.1</strong></code>; <code>minimum = 0.01</code>)</td></tr></table> |
---|
| 302 | Used uniquely for <code>meMode = 103</code> to set the minimal value |
---|
| 303 | assumed for the threshold factor, |
---|
| 304 | <i>sqrt( (1 - m_1^2/m_2 - m_2^2/m^2)^2 - 4 m_1^2 m_2^2/m^4)</i> |
---|
| 305 | for two-body decays and <i>sqrt(1 - Sum_i m_i / m)</i> for multibody |
---|
| 306 | ones. Thus the inverse of this number sets an upper limit for how |
---|
| 307 | much the partial width of a channel can increase from the on-shell |
---|
| 308 | value to the value for asymptotically large resonance masses. Is mainly |
---|
| 309 | intended as a safety measure, to avoid unintentionally large rescalings. |
---|
| 310 | |
---|
| 311 | |
---|
| 312 | <p/> |
---|
| 313 | All of these <code>meMode</code>'s may coexist for the same resonance. |
---|
| 314 | This would be the case e.g. if you want to add a few new channels to an |
---|
| 315 | already existing resonance, where the old partial widths come hardcoded |
---|
| 316 | while the new ones are read in from an external file. The typical example |
---|
| 317 | would be an MSSM Higgs sector, where partial widths to SM particles are |
---|
| 318 | already encoded, <code>meMode = 0</code>, while decay rates to sparticles |
---|
| 319 | are read in from some external calculation and maybe would be best |
---|
| 320 | approximated by using <code>meMode = 103</code>. Indeed the default |
---|
| 321 | particle table in PYTHIA uses 103 for all channels that are expected |
---|
| 322 | to be provided by external input. |
---|
| 323 | |
---|
| 324 | <p/> |
---|
| 325 | Some further clarification may be useful. At initialization the existing |
---|
| 326 | total width and on-shell branching ratios will be updated. For channels |
---|
| 327 | with <code>meMode < 100</code> the originally stored branching ratios |
---|
| 328 | are irrelevant, since the existing code will anyway be used to calculate |
---|
| 329 | the partial widths from scratch. For channels with <code>meMode = 100</code> |
---|
| 330 | or bigger, instead the stored branching ratio is used together with the |
---|
| 331 | originally stored total width to define the correct on-shell partial width. |
---|
| 332 | The sum of partial widths then gives the new total width, and from there |
---|
| 333 | new branching ratios are defined. |
---|
| 334 | |
---|
| 335 | <p/> |
---|
| 336 | In these operations the original sum of branching ratios need not be |
---|
| 337 | normalized to unity. For instance, you may at input have a stored total |
---|
| 338 | width of 1 GeV and a sum of branching ratios of 2. After initialization |
---|
| 339 | the width will then have been changed to 2 GeV and the sum of branching |
---|
| 340 | ratios rescaled to unity. This might happen e.g. if you add a few channels |
---|
| 341 | to an existing resonance, without changing the branching ratios of the |
---|
| 342 | existing channels or the total width of the resonance. |
---|
| 343 | |
---|
| 344 | <p/> |
---|
| 345 | In order to simulate the Breit-Wigner shape correctly, it is important |
---|
| 346 | that all channels that contribute to the total width are included in the |
---|
| 347 | above operations. This must be kept separate from the issue of which |
---|
| 348 | channels you want to have switched on for a particular study, to be |
---|
| 349 | considered next. |
---|
| 350 | |
---|
| 351 | <p/> |
---|
| 352 | |
---|
| 353 | In the event-generation process, when an off-shell resonance mass has been |
---|
| 354 | selected, the width and branching ratios are re-evaluated for this new mass. |
---|
| 355 | At this stage also the effects of restrictions on allowed decay modes are |
---|
| 356 | taken into account, as set by the <code>onMode</code> switch for each |
---|
| 357 | separate decay channel. Thus a channel may be on or off, with different |
---|
| 358 | choices of open channels between the particle and its antiparticle. |
---|
| 359 | In addition, even when a channel is on, the decay may be into another |
---|
| 360 | resonance with its selection of allowed channels. It is these kinds of |
---|
| 361 | restrictions that lead to the <i>Gamma_out</i> possibly being |
---|
| 362 | smaller than <i>Gamma_tot</i>. As a reminder, the Breit-Wigner for |
---|
| 363 | decays behaves like <i>Gamma_out / ((s - m^2)^2 + s * Gamma_tot^2)</i>, |
---|
| 364 | where the width in the numerator is only to those channels being studied, |
---|
| 365 | but the one in the denominator to all channels of the particle. These |
---|
| 366 | ever-changing numbers are not directly visible to the user, but are only |
---|
| 367 | stored in a work area. |
---|
| 368 | |
---|
| 369 | <input type="hidden" name="saved" value="1"/> |
---|
| 370 | |
---|
| 371 | <?php |
---|
| 372 | echo "<input type='hidden' name='filepath' value='".$_GET["filepath"]."'/>"?> |
---|
| 373 | |
---|
| 374 | <table width="100%"><tr><td align="right"><input type="submit" value="Save Settings" /></td></tr></table> |
---|
| 375 | </form> |
---|
| 376 | |
---|
| 377 | <?php |
---|
| 378 | |
---|
| 379 | if($_POST["saved"] == 1) |
---|
| 380 | { |
---|
| 381 | $filepath = $_POST["filepath"]; |
---|
| 382 | $handle = fopen($filepath, 'a'); |
---|
| 383 | |
---|
| 384 | if($_POST["1"] != "1e-20") |
---|
| 385 | { |
---|
| 386 | $data = "ResonanceWidths:minWidth = ".$_POST["1"]."\n"; |
---|
| 387 | fwrite($handle,$data); |
---|
| 388 | } |
---|
| 389 | if($_POST["2"] != "0.1") |
---|
| 390 | { |
---|
| 391 | $data = "ResonanceWidths:minThreshold = ".$_POST["2"]."\n"; |
---|
| 392 | fwrite($handle,$data); |
---|
| 393 | } |
---|
| 394 | fclose($handle); |
---|
| 395 | } |
---|
| 396 | |
---|
| 397 | ?> |
---|
| 398 | </body> |
---|
| 399 | </html> |
---|
| 400 | |
---|
| 401 | <!-- Copyright (C) 2012 Torbjorn Sjostrand --> |
---|
| 402 | |
---|