bachelor_thesis/prog/python/qqgg/analytical_xs.org

#+PROPERTY: header-args :exports both

* Init
** Required Modules
#+NAME: e988e3f2-ad1f-49a3-ad60-bedba3863283
#+BEGIN_SRC ein-python :session :exports both
import numpy as np
import matplotlib.pyplot as plt
#+END_SRC

#+RESULTS: e988e3f2-ad1f-49a3-ad60-bedba3863283

** Utilities
#+NAME: 53548778-a4c1-461a-9b1f-0f401df12b08
#+BEGIN_SRC ein-python :session :exports both :results raw drawer
%run ../utility.py
#+END_SRC

#+RESULTS: 53548778-a4c1-461a-9b1f-0f401df12b08
:RESULTS:
:END:

* Implementation
#+NAME: 777a013b-6c20-44bd-b58b-6a7690c21c0e
#+BEGIN_SRC ein-python :session :exports both :results raw drawer :exports code :tangle xs.py
  """
  Implementation of the analytical cross section for q q_bar ->
  gamma gamma

  Author: Valentin Boettcher <hiro@protagon.space>
  """

  import numpy as np
  from scipy.constants import alpha

  # NOTE: a more elegant solution would be a decorator
  def energy_factor(charge, esp):
      """
      Calculates the factor common to all other values in this module

      Arguments:
      esp -- center of momentum energy in GeV
      charge -- charge of the particle in units of the elementary charge
      """

      return charge**4*(alpha/esp)**2/6


  def diff_xs(theta, charge, esp):
      """
      Calculates the differential cross section as a function of the
      azimuth angle theta in units of 1/GeV^2.

      Arguments:
      theta -- azimuth angle
      esp -- center of momentum energy in GeV
      charge -- charge of the particle in units of the elementary charge
      """

      f = energy_factor(charge, esp)
      return f*((np.cos(theta)**+1)/np.sin(theta)**2)

  def diff_xs_eta(eta, charge, esp):
      """
      Calculates the differential cross section as a function of the
      pseudo rapidity of the photons in units of 1/GeV^2.

      Arguments:
      eta -- pseudo rapidity
      esp -- center of momentum energy in GeV
      charge -- charge of the particle in units of the elementary charge
      """

      f = energy_factor(charge, esp)
      return f*(2*np.cosh(eta)**2 - 1)

  def total_xs_eta(eta, charge, esp):
      """
      Calculates the total cross section as a function of the pseudo
      rapidity of the photons in units of 1/GeV^2.  If the rapditiy is
      specified as a tuple, it is interpreted as an interval.  Otherwise
      the interval [-eta, eta] will be used.

      Arguments:
      eta -- pseudo rapidity (tuple or number)
      esp -- center of momentum energy in GeV
      charge -- charge of the particle in units of the elementar charge
      """

      f = energy_factor(charge, esp)
      if not isinstance(eta, tuple):
          eta = (-eta, eta)

      if len(eta) != 2:
          raise ValueError('Invalid eta cut.')

      def F(x):
          return np.tanh(x) - 2*x

      return 2*np.pi*f*(F(eta[0]) - F(eta[1]))
#+END_SRC

#+RESULTS: 777a013b-6c20-44bd-b58b-6a7690c21c0e
:RESULTS:
:END:

* Calculations
** XS qq -> gamma gamma
First, set up the input parameters.
#+NAME: 7e62918a-2935-41ac-94e0-f0e7c3af8e0d
#+BEGIN_SRC ein-python :session :exports both :results raw drawer
eta = 2.5
charge = 1/3
esp = 200  # GeV
#+END_SRC

#+RESULTS: 7e62918a-2935-41ac-94e0-f0e7c3af8e0d
:RESULTS:
:END:

And now calculate the cross section in picobarn.
#+NAME: cf853fb6-d338-482e-bc55-bd9f8e796495
#+BEGIN_SRC ein-python :session :exports both :results raw drawer
xs_gev = total_xs_eta(eta, charge, esp)
xs_pb = gev_to_pb(xs_gev)
xs_pb
#+END_SRC

#+RESULTS: cf853fb6-d338-482e-bc55-bd9f8e796495
:RESULTS:
0.053793289459925515
:END:

Compared to sherpa, it's pretty close.
#+NAME: 81b5ed93-0312-45dc-beec-e2ba92e22626
#+BEGIN_SRC ein-python :session :exports both :results raw drawer
  sherpa = 0.0538009
  xs_pb/sherpa
#+END_SRC

#+RESULTS: 81b5ed93-0312-45dc-beec-e2ba92e22626
:RESULTS:
0.9998585425137037
:END:

I had to set the runcard option ~EW_SCHEME: alpha0~ to use the pure
QED coupling constant.