PYTHIA8/pythia8145/src/BeamShape.cxx

   1 // BeamShape.cc is a part of the PYTHIA event generator.
   2 // Copyright (C) 2010 Torbjorn Sjostrand.
   3 // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
   4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
   5
   6 // Function definitions (not found in the header) for the BeamShape class.
   7
   8 #include "BeamShape.h"
   9
  10 namespace Pythia8 {
  11
  12 //==========================================================================
  13
  14 // The BeamShape class.
  15
  16 //--------------------------------------------------------------------------
  17
  18 // Initialize beam parameters.
  19
  20   void BeamShape::init( Settings& settings, Rndm* rndmPtrIn) {
  21
  22   // Save pointer.
  23   rndmPtr             = rndmPtrIn;
  24
  25   // Main flags.
  26   allowMomentumSpread = settings.flag("Beams:allowMomentumSpread");
  27   allowVertexSpread   = settings.flag("Beams:allowVertexSpread");
  28
  29   // Parameters for beam A momentum spread.
  30   sigmaPxA            = settings.parm("Beams:sigmaPxA");
  31   sigmaPyA            = settings.parm("Beams:sigmaPyA");
  32   sigmaPzA            = settings.parm("Beams:sigmaPzA");
  33   maxDevA             = settings.parm("Beams:maxDevA");
  34
  35   // Parameters for beam B momentum spread.
  36   sigmaPxB            = settings.parm("Beams:sigmaPxB");
  37   sigmaPyB            = settings.parm("Beams:sigmaPyB");
  38   sigmaPzB            = settings.parm("Beams:sigmaPzB");
  39   maxDevB             = settings.parm("Beams:maxDevB");
  40
  41   // Parameters for beam vertex spread.
  42   sigmaVertexX        = settings.parm("Beams:sigmaVertexX");
  43   sigmaVertexY        = settings.parm("Beams:sigmaVertexY");
  44   sigmaVertexZ        = settings.parm("Beams:sigmaVertexZ");
  45   maxDevVertex        = settings.parm("Beams:maxDevVertex");
  46   sigmaTime           = settings.parm("Beams:sigmaTime");
  47   maxDevTime          = settings.parm("Beams:maxDevTime");
  48
  49   // Parameters for beam vertex offset.
  50   offsetX             = settings.parm("Beams:offsetVertexX");
  51   offsetY             = settings.parm("Beams:offsetVertexY");
  52   offsetZ             = settings.parm("Beams:offsetVertexZ");
  53   offsetT             = settings.parm("Beams:offsetTime");
  54
  55 }
  56
  57 //--------------------------------------------------------------------------
  58
  59 // Set the two beam momentum deviations and the beam vertex.
  60
  61 void BeamShape::pick() {
  62
  63   // Reset all values.
  64   deltaPxA = deltaPyA = deltaPzA = deltaPxB = deltaPyB = deltaPzB
  65     = vertexX = vertexY = vertexZ = vertexT = 0.;
  66
  67   // Set beam A momentum deviation by a three-dimensional Gaussian.
  68   if (allowMomentumSpread) {
  69     double totalDev, gauss;
  70     do {
  71       totalDev = 0.;
  72       if (sigmaPxA > 0.) {
  73         gauss     = rndmPtr->gauss();
  74         deltaPxA  = sigmaPxA * gauss;
  75         totalDev += gauss * gauss;
  76       }
  77       if (sigmaPyA > 0.) {
  78         gauss     = rndmPtr->gauss();
  79         deltaPyA  = sigmaPyA * gauss;
  80         totalDev += gauss * gauss;
  81       }
  82       if (sigmaPzA > 0.) {
  83         gauss     = rndmPtr->gauss();
  84         deltaPzA  = sigmaPzA * gauss;
  85         totalDev += gauss * gauss;
  86       }
  87     } while (totalDev > maxDevA * maxDevA);
  88
  89     // Set beam B momentum deviation by a three-dimensional Gaussian.
  90     do {
  91       totalDev = 0.;
  92       if (sigmaPxB > 0.) {
  93         gauss     = rndmPtr->gauss();
  94         deltaPxB  = sigmaPxB * gauss;
  95         totalDev += gauss * gauss;
  96       }
  97       if (sigmaPyB > 0.) {
  98         gauss     = rndmPtr->gauss();
  99         deltaPyB  = sigmaPyB * gauss;
 100         totalDev += gauss * gauss;
 101       }
 102       if (sigmaPzB > 0.) {
 103         gauss     = rndmPtr->gauss();
 104         deltaPzB  = sigmaPzB * gauss;
 105         totalDev += gauss * gauss;
 106       }
 107     } while (totalDev > maxDevB * maxDevB);
 108   }
 109
 110   // Set beam vertex location by a three-dimensional Gaussian.
 111   if (allowVertexSpread) {
 112     double totalDev, gauss;
 113     do {
 114       totalDev = 0.;
 115       if (sigmaVertexX > 0.) {
 116         gauss     = rndmPtr->gauss();
 117         vertexX   = sigmaVertexX * gauss;
 118         totalDev += gauss * gauss;
 119       }
 120       if (sigmaVertexY > 0.) {
 121         gauss     = rndmPtr->gauss();
 122         vertexY   = sigmaVertexY * gauss;
 123         totalDev += gauss * gauss;
 124       }
 125       if (sigmaVertexZ > 0.) {
 126         gauss     = rndmPtr->gauss();
 127         vertexZ   = sigmaVertexZ * gauss;
 128         totalDev += gauss * gauss;
 129       }
 130     } while (totalDev > maxDevVertex * maxDevVertex);
 131
 132     // Set beam collision time by a Gaussian.
 133     if (sigmaTime > 0.) {
 134       do gauss    = rndmPtr->gauss();
 135       while (abs(gauss) > maxDevTime);
 136       vertexT     = sigmaTime * gauss;
 137     }
 138
 139     // Add offset to beam vertex.
 140     vertexX      += offsetX;
 141     vertexY      += offsetY;
 142     vertexZ      += offsetZ;
 143     vertexT      += offsetT;
 144   }
 145
 146 }
 147
 148 //==========================================================================
 149
 150 } // end namespace Pythia8
 151