Statistical approaches in high energy physics using ROOT data analysis framework

Samgeeth P

IISER Trivandrum, Maruthamala P.O, Vithura, Thiruvananthapuram, Kerala, India 695551

Dr. Prolay K Mal

School of Physics, NISER Bhubaneswar, P.O. Jatni, Khurda 752050, Odisha, India

Abstract

This project was an attempt to provide an idea of how the software ROOT is used to statistically analyse the events that occur in high energy physics with an example where the Bhabha scattering cross-section was computed. The approach used for this study included both acquisition of theoretical and computational knowledge and application of the techniques learned on various tasks. The theory background gathered comprised of statistics, probability, relativity and particle physics. The fundamentals of particles and their relativistic aspects were learned along with information about particle collisions, important particle colliders, particle detectors and their components. A method in statistics and probability, which was of significant importance to the project, was Monte Carlo integration, which makes use of random numbers to find solutions for problems with enormous amount of data. The random numbers were generated using ROOT data analysis framework. Two techniques, namely the inversion technique and accept-reject method helped in confining the random numbers to a desired sample space. Out of the two, feasibility of one technique over another was studied. Methods of sampling like important sampling and stratified sampling were learned and their use in increasing the efficiency of Monte Carlo integration was also noted. All the knowledge acquired was put together for the computation of Bhabha scattering cross-section and plotting the relation between differential cross-section and the scattering angle. The task ultimately was successful, as the value of total cross-section was obtained and the plot of the differential cross-section versus scattering angle matched the expected nature of the plot. More background information, methodologies used and deeper insight to the results are provided as we go further into the report.

Keywords: Bhabha scattering cross-section, particle collisions, Monte Carlo technique, random number generation, differential cross-section, scattering angle

Abbreviations

CMS	Compact Muon Solenoid
MC	Monte Carlo
CERN	European Organization for Nuclear Research
NISER	National Institute of Science Education and Research
PRNG	Pseudo Random Number Generator
pdf	Probability Distribution Function
cdf	Cumulative Distribution Function
CM	Center of Mass
HEP	High Energy Physics

INTRODUCTION

Random Numbers

What are random numbers?

We all live in a universe that always increases its randomness. From particles of the size of an electron, to huge clusters with billions of stars, everything tends to increase their randomness over time. Which is why simulations with very systematic arrangements of events might not always yield data that can perfectly explain the actual events that occured. In order to simulate an event with minimal deviations from the actual experimental data, we need the simulation itself to be random. This is where random numbers come into the picture. Random numbers are a set of numbers with no pattern at all. The outcome of a simulation using these numbers provides a much better depiction of an an actual event.

Pseudo-random numbers

Suppose a blind-folded person picks numbers from a box full of numbers. He can start from anywhere. He just needs to keep on picking with no pattern till he took the last one. In other words, the numbers he has are completely random. But, the case changes if we generate numbers using a computer. The systems we use cannot just blindly pick numbers like that, because the machines themselves are much more systematic than us. Which is why it always has to follow a pattern in picking the numbers. The scientists have come up with algorithms that maximize the randomness of this pattern. But after that limit, it will repeat the pattern. However, numbers generated in such a way can be considered random for almost all practical puposes. Such numbers are called pseudo-random numbers. Since all computers generate random numbers based on some algorithm, all computer generated random numbers are pseudo-random numbers. We will learn methods to generate these random numbers in section 2.1

Monte Carlo

Monte Carlo technique is a statistical technique for finding the integral of functions using random numbers. Instead of taking stepwise increasing inervals, the MC technique randomly chooses an interval and adds it to the sum. Which is why random numbers are used for MC integration.

Why Monte Carlo?

In a normal one dimensional MC integral, the error of the estimate is proportional to (1/√n). But for numerical methods (even the simplest ones), the accuracy is proportional to 1/n² which makes one wonder if MC is good enough for finding accurate values of integrals. But the benefit of choosing MC over other techniques is when the number of variables (dimensions) is more. The general formula for error in simple numerical methods is actually proportional to n^-2/d (where d is the number of dimensions), while for MC, the error remains inversely proportional to the square root of n, which means that for higher dimensions (more than 4, which can be frequent in real life scenarios), MC becomes a better method with less error.

Another fact is that for real life cases, the boundary of the multi-dimensional region may be so poorly defined that it might become harder to perform a simple integration over it, while MC only needs to check whether the point is inside the curve or not to either accept or reject the sample (remember accept-reject method from section 2.1.3)

Also, each time we evaluate an integral, MC provides an answer slightly different than before, which implies it is very suitable to simulate real life scenarios, because real life examples almost always contains random events.

ROOT

About ROOT

ROOT is an object oriented program and library developed by CERN. It is a modular scientific toolkit originally designed for particle physics data analysis purposes with features specifically designed for this field. But now it is also being used in other fields like astronomy and data mining. It is mainly written in C++, but can be integrated with other languages like Python and R.

Being a framework, ROOT provides many utilities like input outpit functions and graphics. Also, being a framework primarily designed for HEP, there are utilities such as histograms and fitting. Its like the utility belt of Batman that lets him do a lot of things that he couldn't do before. This makes the codes less bulky, more robust, more consistant and much easier.

Another fact that makes ROOT easier is the autoload of libraries. While reading a code, ROOT automatically loads the libraries so that even if you forget to include a library, ROOT gets the job done for you. This can be very useful when codes that need many libraries to be loaded are compiled.

An interesting fact about ROOT is that it uses a command line interpreter method. Which means the executable need not be made always to run the code that was written.

We shall now discuss a few utilities of ROOT that we will be using in this project. Note that the codes are shown only in appendix part.

TCanvas

TCanvas is the class for creating and editing canvases. It can be used to add a canvas, divide it into rows of columns to draw multiple objects, draw objects with different scales, overlay objects etc.

Histograms

TH1,TH2 and TH3 are a set of classes representing one dimensional, two dimensional and three dimensional histograms respectively. Each class can be further divided into five depending on the bin content. For example, TH1C, TH1S, TH1I, TH1F and TH1D are the five one dimensional histogram classes.

Histograms can have constant bin width, or we can set the width of each bin seperately. For both, the number of bins must be defined and for the constant bin width, we need to give only the starting and ending point of the histogram, whereas all the bin edges should be made in an array and given to make a histogram with variable bin width.

Filling histograms can be done in two ways. First way is to fill each value the desired number of times using Fill. Another way is to use SetBinContent to give the value of each bin content.

TGraph, TGraph2D and TGraphErrors

TGraph classes are used to plot graphs. Graphs can be created by defining a set of coordinates as arrays or just creating a graph first and adding points later. TGraphErrors comes with the ability to assign error bars to each point. There are a few Draw options by which the graph can be styled as we want.

Functions and fitting

Functions can be defined using the classes TF1, TF2 and TF3. There are a few pre defined functions such as the gaussian or the poisson distributions. These functions can later be used to fit histograms with Fit to find out the distribution the histograms are following. Both pre-defined functions and functions we defined can be used to fit histograms.

gRandom

This is a call for a class TRandom3, which is basically a random number generator. It uses an algorithm to generate pseudo random numbers. More information regarding this is given in section 2.1.1

LITERATURE REVIEW

Random Number Generation

PRNG algorithms

As we saw before, the concept of pseudo random numbers are based on algorithms. The values generated by these algorithms are usually based on a fixed number called a seed. One of the most common algorithms used for this is the linear congruential generator, given by the relation,

$\;\;\;\;\;\;\;\;\;\;X_{n+1}=(aX_n+b)\mod m$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><msub><mi>X</mi><mrow><mi>n</mi><mo>+</mo><mn>1</mn></mrow></msub><mo>=</mo><mo>(</mo><mi>a</mi><msub><mi>X</mi><mi>n</mi></msub><mo>+</mo><mi>b</mi><mo>)</mo><mi>mod</mi><mi>m</mi></math>\;\;\;\;\;\;\;\;\;\;X_{n+1}=(aX_n+b)\mod m

where X is a series of random numbers with n^th element X_n and a,b and m are large integers. The sequence can continue generating random numbers upto n=m-1.

The PRNG we are using, however, is the gRandom generator in ROOT, which uses the Mersenne Twister generator. The algorithm can generate upto 10⁶⁰⁰⁰ random numbers and is very fast.

The PRNG sequences are usually programmed to generate a uniform distribution of numbers between 0 and 1. We will now look at the methods by which we can span this over a desired interval with a desired probability distribution.

Inversion technique

Inversion technique makes use of the probability density function (pdf) of a given distribution to generate random numbers according to that distribution. Pdf is nothing but the function itself with a normalisation factor such that the integral of the function over the whole domain is equal to 1. It is denoted by p(x)

$\;\;\;\;\;\;\;\;\;\;\int\limits_Xp(x)dx=1$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><msub><mo>∫</mo><mi>x</mi></msub><mi>p</mi><mo>(</mo><mi>x</mi><mo>)</mo><mi>d</mi><mi>x</mi><mo>=</mo><mn>1</mn></math>\;\;\;\;\;\;\;\;\;\;\int\limits_Xp(x)dx=1

Here, the technique is to integrate the pdf from the lower limit to an arbitrary value x. This integral is called the cumulative distribution function (cdf), denoted by F(x).

$\;\;\;\;\;\;\;\;\;\;F(x)=\int\limits_{-\infty}^xp(x)dx$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mi>F</mi><mo>(</mo><mi>x</mi><mo>)</mo><mo>=</mo><msubsup><mi>\intlimits</mi><mrow><mo>-</mo><mo>∞</mo></mrow><mi>x</mi></msubsup><mi>p</mi><mo>(</mo><mi>x</mi><mo>)</mo><mi>d</mi><mi>x</mi></math>\;\;\;\;\;\;\;\;\;\;F(x)=\int\limits_{-\infty}^xp(x)dx

This cdf is the uniform random number we are generating (say, z). So, if we denote x as a function of z, that is, x=F^-1(z), we can generate x randomly with the probability distribution.

Box Muller Method:

This is a special case of the inversion technique, used specifically to generate Gaussian distributed numbers with mean 0 and standard deviation 1. Instead of defining the desired random number as a function of one random number, Box and Muller generated two random numbers and transformed them into two Gaussian distributed functions. The equations for the Gaussian distributed random numbers are given by,

$\begin{array}{l}\;\;\;\;\;\;\;\;\;\;x_1=\sqrt{-2\ln(z_1)}\cos(2\mathrm\pi z_2)\\\;\;\;\;\;\;\;\;\;\;x_2=\sqrt{-2\ln({\mathrm z}_1)}\sin(2\mathrm\pi z_2)\end{array}$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mtable columnalign="left"><mtr><mtd><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><msub><mi>x</mi><mn>1</mn></msub><mo>=</mo><msqrt><mo>-</mo><mn>2</mn><mi>ln</mi><mo>(</mo><msub><mi>z</mi><mn>1</mn></msub><mo>)</mo></msqrt><mi>cos</mi><mo>(</mo><mn>2</mn><mi mathvariant="normal">π</mi><msub><mi>z</mi><mn>2</mn></msub><mo>)</mo></mtd></mtr><mtr><mtd><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><msub><mi>x</mi><mn>2</mn></msub><mo>=</mo><msqrt><mo>-</mo><mn>2</mn><mi>ln</mi><mo>(</mo><msub><mi mathvariant="normal">z</mi><mn>1</mn></msub><mo>)</mo></msqrt><mi>sin</mi><mo>(</mo><mn>2</mn><mi mathvariant="normal">π</mi><msub><mi>z</mi><mn>2</mn></msub><mo>)</mo></mtd></mtr></mtable></math>\begin{array}{l}\;\;\;\;\;\;\;\;\;\;x_1=\sqrt{-2\ln(z_1)}\cos(2\mathrm\pi z_2)\\\;\;\;\;\;\;\;\;\;\;x_2=\sqrt{-2\ln({\mathrm z}_1)}\sin(2\mathrm\pi z_2)\end{array}

The benefit of using inversion technique is that we can use all the random numbers generated, as all are coming according to the desired distribution. But the catch is that if the cdf is not invertible, then the whole idea of inversion technique becomes useless. Also, if the exact equation of the function is unknown, then again the method will become useless.

Accept-reject method

As we saw before, the inversion technique can only generate functions whose cdf are invertible. So, we need a way to generate other functions as well. Accept-reject method is based on a trial and error approach where two random numbers x and y are generated simultaneously and if y is less than the desired value of the function of x (ie, if y<f(x)), then the value x is accepted.

The benefit of accept-reject method is that it can produce random numbers in almost any distribution. Which means it can be used as an alternative for functions who aren't invertible. But the problem is, the rejected numbers are technically wasted. Which means if the region where y is generated is much bigger than the curve f(x), then a lot of numbers are lost in the process. This makes the accept-reject method a very wasteful process if proper sampling is not done.

The random number generation mentioned above can have many applications. One of which we are going to discuss here is Monte Carlo.

Techniques of Monte Carlo

Monte Carlo integration

The equation for simple MC integral is given by,

$\int\limits_{V} f dV = V \langle f \rangle \pm V \sqrt{{ \langle f^{2} \rangle - \langle f \rangle ^{2} } \over {n}}\;\;\;where\;\;\;\langle f \rangle = {{1}\over{n}}\sum\limits_{i=1}^{n} f_{i}\;\;;\;\;\langle f^{2} \rangle = {{1} \over {n}} \sum\limits_{i=1}^{n} f_{i}^{2}$ <math xmlns="http://www.w3.org/1998/Math/MathML"><munder><mo>∫</mo><mi>V</mi></munder><mi>f</mi><mi>d</mi><mi>V</mi><mo>=</mo><mi>V</mi><mo>⟨</mo><mi>f</mi><mo>⟩</mo><mo>±</mo><mi>V</mi><msqrt><mfrac><mrow><mo>⟨</mo><msup><mi>f</mi><mn>2</mn></msup><mo>⟩</mo><mo>-</mo><mo>⟨</mo><mi>f</mi><msup><mo>⟩</mo><mn>2</mn></msup></mrow><mi>n</mi></mfrac></msqrt><mo> </mo><mo> </mo><mo> </mo><mi>w</mi><mi>h</mi><mi>e</mi><mi>r</mi><mi>e</mi><mo> </mo><mo> </mo><mo> </mo><mo>⟨</mo><mi>f</mi><mo>⟩</mo><mo>=</mo><mfrac><mn>1</mn><mi>n</mi></mfrac><munderover><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>1</mn></mrow><mi>n</mi></munderover><msub><mi>f</mi><mi>i</mi></msub><mo> </mo><mo>,</mo><mo> </mo><mo> </mo><mo> </mo><mo>⟨</mo><msup><mi>f</mi><mn>2</mn></msup><mo>⟩</mo><mo>=</mo><mfrac><mn>1</mn><mi>n</mi></mfrac><munderover><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>1</mn></mrow><mi>n</mi></munderover><msubsup><mi>f</mi><mi>i</mi><mn>2</mn></msubsup></math>\int\limits_{V} f dV = V \langle f \rangle \pm V \sqrt{{ \langle f^{2} \rangle - \langle f \rangle ^{2} } \over {n}}\;\;\;where\;\;\;\langle f \rangle = {{1}\over{n}}\sum\limits_{i=1}^{n} f_{i}\;\;;\;\;\langle f^{2} \rangle = {{1} \over {n}} \sum\limits_{i=1}^{n} f_{i}^{2}

Here n is the number of events generated and V is the total volume of the region in which the numbers are generated.

Important sampling in Monte Carlo

The process of Monte Carlo is more or less like based on the accept-reject method. which means in order to avoid losing too many random numbers in the process of rejection, we need to define the region to generate random numbers more specifically. Here we use a technique called important sampling, where we define a region whose boundaries are closer to the boundaries of the function than a uniform distribution. The region is chosen in such a way that the probability distribution of the region still remains equal to 1. The idea can be represented mathematically by,

$\;\;\;\;\;\;\;\;\;\;\int\limits_{V} f dV = \int \big( {{f} \over {p}} \big) p dV = \big\langle {{f} \over {p}} \big\rangle \pm \sqrt{ { \langle f^{2}/p^{2} \rangle - \langle f/p \rangle^{2} } \over {n} }\;\;\;where\;\;\;\int\limits_{V} p dV = 1$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><munder><mo>∫</mo><mi>V</mi></munder><mi>f</mi><mi>d</mi><mi>V</mi><mo>=</mo><mo>∫</mo><mi>\big</mi><mo>(</mo><mfrac><mi>f</mi><mi>p</mi></mfrac><mi>\big</mi><mo>)</mo><mi>p</mi><mi>d</mi><mi>V</mi><mo>=</mo><mi>\big</mi><mo>⟨</mo><mfrac><mi>f</mi><mi>p</mi></mfrac><mi>\big</mi><mo>⟩</mo><mo>±</mo><msqrt><mfrac><mrow><mo>⟨</mo><msup><mi>f</mi><mn>2</mn></msup><mo>/</mo><msup><mi>p</mi><mn>2</mn></msup><mo>⟩</mo><mo>-</mo><mo>⟨</mo><mi>f</mi><mo>/</mo><mi>p</mi><msup><mo>⟩</mo><mn>2</mn></msup></mrow><mi>n</mi></mfrac></msqrt><mo> </mo><mo> </mo><mo> </mo><mi>w</mi><mi>h</mi><mi>e</mi><mi>r</mi><mi>e</mi><mo> </mo><mo> </mo><mo> </mo><munder><mo>∫</mo><mi>V</mi></munder><mi>p</mi><mi>d</mi><mi>V</mi><mo>=</mo><mn>1</mn></math>\;\;\;\;\;\;\;\;\;\;\int\limits_{V} f dV = \int \big( {{f} \over {p}} \big) p dV = \big\langle {{f} \over {p}} \big\rangle \pm \sqrt{ { \langle f^{2}/p^{2} \rangle - \langle f/p \rangle^{2} } \over {n} }\;\;\;where\;\;\;\int\limits_{V} p dV = 1

Stratified sampling in Monte Carlo

Just like important sampling, stratified sampling also enables us to minimize the wastage of random numbers due to rejection. As the name suggests, the region is divided into seperate strata and numbers are generated in each strata. This way, each strata can be assigned with an enclosing region with a suitable shape to minimize wastage. The equation follows directly from the integration formula,

$\begin{array}{l}\;\;\;\;\;\;\;\;\;\;\int\limits_{V}fdV=\int\limits_{V_{1}}fdV+\int\limits_{V_{2}}fdV+\int\limits_{V_{3}}fdV\\\;\;\;\;\;\;\;\;\;\;where\;\;\;V=V_{1} \cup V_{2} \cup V_{3}\;\;\;and\;\;\;V_{1} \cap V_{2} \cap V_{3}=\phi\end{array}$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mtable columnalign="left"><mtr><mtd><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><munder><mo>∫</mo><mi>V</mi></munder><mi>f</mi><mi>d</mi><mi>V</mi><mo>=</mo><munder><mo>∫</mo><msub><mi>V</mi><mn>1</mn></msub></munder><mi>f</mi><mi>d</mi><mi>V</mi><mo>+</mo><munder><mo>∫</mo><msub><mi>V</mi><mn>2</mn></msub></munder><mi>f</mi><mi>d</mi><mi>V</mi><mo>+</mo><munder><mo>∫</mo><msub><mi>V</mi><mn>3</mn></msub></munder><mi>f</mi><mi>d</mi><mi>V</mi></mtd></mtr><mtr><mtd><mi>w</mi><mi>h</mi><mi>e</mi><mi>r</mi><mi>e</mi><mo> </mo><mo> </mo><mo> </mo><mi>V</mi><mo>=</mo><msub><mi>V</mi><mn>1</mn></msub><mo>∪</mo><msub><mi>V</mi><mn>2</mn></msub><mo>∪</mo><msub><mi>V</mi><mn>3</mn></msub><mo> </mo><mo> </mo><mo> </mo><mi>a</mi><mi>n</mi><mi>d</mi><mo> </mo><mo> </mo><mo> </mo><msub><mi>V</mi><mn>1</mn></msub><mo>∩</mo><msub><mi>V</mi><mn>2</mn></msub><mo>∩</mo><msub><mi>V</mi><mn>3</mn></msub><mo>=</mo><mi>ϕ</mi></mtd></mtr></mtable></math>\begin{array}{l}\;\;\;\;\;\;\;\;\;\;\int\limits_{V}fdV=\int\limits_{V_{1}}fdV+\int\limits_{V_{2}}fdV+\int\limits_{V_{3}}fdV\\\;\;\;\;\;\;\;\;\;\;where\;\;\;V=V_{1} \cup V_{2} \cup V_{3}\;\;\;and\;\;\;V_{1} \cap V_{2} \cap V_{3}=\phi\end{array}

(*Examples of Monte carlo techniques generated will be given in methodologies part)

Bhabha Scattering

Scattering is a common process that happens when two particles interact with each other. Bhabha scattering is the process of electron-positron scattering In simple ways, it can be written as e⁺e^-→e⁺e^-

There are two ways Bhabha scattering can occur: s channel and t channel, where s and t are the Mandelstam variables.

s channel Bhabha scattering is when the electron and positron annihilate to produce a photon, which in turn produce another pair of electron and positron.

t channel Bhabha scattering is when the electron and positron simply interact with each other by the exchange of a photon (electromagnetic interaction).

s channel Bhabha scattering

t channel Bhabha scattering

Feynman diagrams for Bhabha scattering (horizontal axis shows the flow of time)

Figure shows the Feynman diagrams for both s and t channel Bhabha Scattering. Feynman diagrams are diagrams that represent the interactions between particles with the flow of time. One of the axes will be represent the flow of time. A line perpendicular to that axis will represent the particles and interactions that exist at that instant of time.

In our project, we consider both interactions to find out the total Bhabha scattering cross section.

Scattering Cross Section

Cross section

In order for two particles to scatter, the particles should come within a certain area. This area transverse to the motion of the particles within which they should meet in order to interact with each other is called cross section of the particles. Cross sections are usually denoted by σ and has the units of area. In particle physics, the most commonly used unit for cross section is barn (1b = 10^-24cm² or 0.389379mb = 1GeV^-2).

Differential cross section

Differential cross section is the cross section represented as a function of a final-state parameter (such as scattering angle or particle energy). In this project, we shall represent the differential cross section as the function of scattering angle θ. For Bhabha scattering, the differential cross section is given by,

$\;\;\;\;\;\;\;\;\;\;\big( {{d\sigma} \over {d\Omega}} \big) = {{\alpha^{2}} \over {2s}} \Bigg( {{1+cos^{4}\big({{θ} \over {2}} \big) } \over {sin^{4}\big({{θ}\over{2}}\big)}} - {{2cos^{4}\big({{θ} \over {2}} \big) } \over {sin^{2}\big({{θ}\over{2}}\big)}} +{{1+cos^{2}({θ}) } \over {2}} \Bigg)$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mi>\big</mi><mo>(</mo><mfrac><mrow><mi>d</mi><mi>σ</mi></mrow><mrow><mi>d</mi><mi>Ω</mi></mrow></mfrac><mi>\big</mi><mo>)</mo><mo>=</mo><mfrac><msup><mi>α</mi><mn>2</mn></msup><mrow><mn>2</mn><mi>s</mi></mrow></mfrac><mi>\Bigg</mi><mo>(</mo><mfrac><mrow><mn>1</mn><mo>+</mo><mi>c</mi><mi>o</mi><msup><mi>s</mi><mn>4</mn></msup><mi>\big</mi><mo>(</mo><mfrac><mi>θ</mi><mn>2</mn></mfrac><mi>\big</mi><mo>)</mo></mrow><mrow><mi>s</mi><mi>i</mi><msup><mi>n</mi><mn>4</mn></msup><mi>\big</mi><mo>(</mo><mfrac><mi>θ</mi><mn>2</mn></mfrac><mi>\big</mi><mo>)</mo></mrow></mfrac><mo>-</mo><mfrac><mrow><mn>2</mn><mi>c</mi><mi>o</mi><msup><mi>s</mi><mn>4</mn></msup><mi>\big</mi><mo>(</mo><mfrac><mi>θ</mi><mn>2</mn></mfrac><mi>\big</mi><mo>)</mo></mrow><mrow><mi>s</mi><mi>i</mi><msup><mi>n</mi><mn>2</mn></msup><mi>\big</mi><mo>(</mo><mfrac><mi>θ</mi><mn>2</mn></mfrac><mi>\big</mi><mo>)</mo></mrow></mfrac><mo>+</mo><mfrac><mrow><mn>1</mn><mo>+</mo><mi>c</mi><mi>o</mi><msup><mi>s</mi><mn>2</mn></msup><mo>(</mo><mi>θ</mi><mo>)</mo></mrow><mn>2</mn></mfrac><mi>\Bigg</mi><mo>)</mo></math>\;\;\;\;\;\;\;\;\;\;\big( {{d\sigma} \over {d\Omega}} \big) = {{\alpha^{2}} \over {2s}} \Bigg( {{1+cos^{4}\big({{θ} \over {2}} \big) } \over {sin^{4}\big({{θ}\over{2}}\big)}} - {{2cos^{4}\big({{θ} \over {2}} \big) } \over {sin^{2}\big({{θ}\over{2}}\big)}} +{{1+cos^{2}({θ}) } \over {2}} \Bigg)

Where α is the fine structure constant (which has a value of approx. 1/137), s is the square of center of mass energy and θ is the scattering angle. The unit of differential cross section is barn/steradians (b/sr).

Total cross section

When the differential cross section is integrated over the entire phase space, it is called the total cross section, denoted by σ_T. It is given by the integral,

$\;\;\;\;\;\;\;\;\;\;\sigma_T=\oint\limits_{\Omega} \big({d\sigma \over d\Omega}\big)d\Omega=\iint\big({d\sigma \over d\Omega}\big)\sin(\theta)d\theta d\varphi$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><mo> </mo><msub><mi>σ</mi><mi>T</mi></msub><mo>=</mo><munder><mo>∮</mo><mi>Ω</mi></munder><mi>\big</mi><mo>(</mo><mfrac><mrow><mi>d</mi><mi>σ</mi></mrow><mrow><mi>d</mi><mi>Ω</mi></mrow></mfrac><mi>\big</mi><mo>)</mo><mi>d</mi><mi>Ω</mi><mo>=</mo><mo>∬</mo><mi>\big</mi><mo>(</mo><mfrac><mrow><mi>d</mi><mi>σ</mi></mrow><mrow><mi>d</mi><mi>Ω</mi></mrow></mfrac><mi>\big</mi><mo>)</mo><mi>sin</mi><mo>(</mo><mi>θ</mi><mo>)</mo><mi>d</mi><mi>θ</mi><mi>d</mi><mi>φ</mi></math>\;\;\;\;\;\;\;\;\;\;\sigma_T=\oint\limits_{\Omega} \big({d\sigma \over d\Omega}\big)d\Omega=\iint\big({d\sigma \over d\Omega}\big)\sin(\theta)d\theta d\varphi

METHODOLOGY

Learning ROOT

The Basic methodologies used in ROOT were studied during the first few days of the project. The available classes were read upon and some codes were compiled.

Random number generation

Random numbers were generated with different distributions. Since the first ones were simple ones, inversion technique was used to generate them. The code is given in Appendix I

Random numbers generated from 0 to 1 uniformly

Random numbers generated from 1 to 2 with a distribution f(x)=kx²

Random numbers generated from 0 to 4×10^-6 with a distribution f(x)=ke^-x

Different trials where random numbers were generated using inversion technique

Generating A Gaussian

The Gaussian distribution was generated with Box-Muller method. Also Gaussean was generated using accept-reject method, importaint sampling and stratified sampling. Box-Muller generated histogram was plotted on the same canvas as the other three to see how much of a deviation is shown by the other three. The code is given in Appendix II.

Accept-reject method vs Box-Muller method

Importance sampling vs Box-Muller method

Stratified sampling vs Box-Muller method

Three types of accept-reject methods with Box-Muller method. Notice that the number of accepted entries in each case increases.

Sample Monte Carlo: Finding center of mass of A Torus

Before moving into the actual work on total cross section, a sample MC was done. A torus of outer radius 4 and inner radius 2 was considered. The portion of torus enclosed by the box 1≤x≤4, -3≤y≤4 and -1≤z≤1 was taken and its CM was found using MC integration. Also, the CM of a piece of torus of the same dimensions but a z-dependant density (ρ=e^5z) was also found out using MC integration. The figures for the pieces of tori are given below. The code is provided in Appendix III.

Torus with uniform density

torus with density ρ=e^5z

Pieces of tori of outer radius 4 and inner radius 2 enclosed by the box 1≤x≤4, -3≤y≤4 and -1≤z≤1

Finding Bhabha Scattering Cross Section

The equation from section 2.4.2 was taken as the differential cross section.

For values of: √s=208GeV, 0.4≤θ≤π-0.4 and 0≤φ≤2π ,

the differential cross section was integrated using MC techniques. The value of θ was taken from 0.4 to π-0.4 because we are trying to make the simulation as realistic as possible and for an actual Bhabha scattering occuring in a particle collider, the value of theta can only be within this range, because the angles beyond that is where the beam pipes are located. The particle detectors cannot detect electrons hitting there. Hence to make the simulation more agreeable to a real life scenario, we took the same range of θ for our simulation.

RESULTS AND DISCUSSION

The plot between differential cross section and scattering angle θ

The code to compute the total cross section and plot the differential cross section vs scattering angle is given in Appendix IV.

For the values of s, θ and φ mentioned in section 3.2, the total cross section was found to be,

$\;\;\;\;\;\;\;\;\;\;\sigma_{T}=120.739±0.022 \;pb$ <math xmlns="http://www.w3.org/1998/Math/MathML"><mi>\sgma</mi></math>\;\;\;\;\;\;\;\;\;\;\sigma_{T}=120.739±0.022 \;pb

A plot of the differental cross section w.r.t θ was plotted and it was seen that the differential cross section is the highest for lower values of θ. Which is reasonable because scattering in higher angles requires more change in momentum, ie, more force. So it should be relatively easier for smaller angles to scatter.
The differential cross section shows no dependance on φ, which is also justifiable because given equation has azimuthal symmetry.

CONCLUSION

In this project, our main objective was to learn about the statistical methods used in HEP, like the MC integration. Also, we planned to learn the functionalities of ROOT and how it is commonly used in HEP experiments. While we succeeded in finding out the Bhabha scattering cross section and explaining the dependance of differential cross section on the scattering angle, much is to be learned in both the physics part and the computational part of HEP experiments. This was just an introduction to the methods that are commonly used in these experiments. Different aspects of random number generation and MC integration was also studied. It was understood that although the experiment we did was basically a single variable integration, for which the accuracy of MC is comparitively less, this was purely meant to be an exercise for the actual experiments that are to be learned. This is just the dawn of the long day that lies ahead of us in the field of future HEP experiments. We hope to walk that path towards the mysteries yet to be solved in the world of HEP.

REFERENCES

Lyons, L. (1986). Statistics for Nuclear and Particle Physicists. Cambridge: Cambridge University Press. doi:10.1017/CBO9781139167710
Griffiths, D. (2008). Introduction to Elementary Particles. Weinheim: Wiley-VCH Verlag GmbH & Co. KgaA. doi:10.1002/9783527618460
Press, W. H. (1992). Numerical recipes in C: The art of scientific computing. Cambridge: Cambridge University Press.
Wikipedia contributors. (2019, June 25). Cross section (physics). In Wikipedia, The Free Encyclopedia.
Wikipedia contributors. (2018, December 4). Bhabha scattering. In Wikipedia, The Free Encyclopedia.
Wikipedia contributors. (2019, May 2). Barn (unit). In Wikipedia, The Free Encyclopedia.
Wikipedia contributors. (2019, June 25). ROOT. In Wikipedia, The Free Encyclopedia.

ACKNOWLEDGEMENTS

First of all, I would like to thank Indian Academy of Sciences for giving me an oppertunity to do this project. Next I thank Dr.Prolay K Mal for guiding me throughout the course of the project, frequently giving instructions. Then I thank all Post Doctoral and PhD fellows from CMS lab in NISER for being there for me whenever I had doubts. Also, I thank my lab-mates Satyabrata Sahoo and Vinay D R for making this internship more lively and enjoyable. The help of Dr.Sudakshina Prusty, the summer internship coordinator at NISER made this internship much more effective, so I thank her too. I use this occasion to thank NISER for providing me accomodation and other services during this period. Last, but not least, I thank everyone who is not mentioned here who made this internship a successful venture, without whose support I wouldn't have been able to complete the same.

APPENDICES

Appendix I

Code for uniform random number

Code for f(x)=kx²

Code for f(x)=ke^-x

Appendix II

Code for accept-reject method overlaid with Box-Muller method

Code for importance sampling overlaid with Box-Muller method

Code for stratified sampling overlaid with Box-Muller method

Appendix III

Code for piece of torus with constant density

Code for piece of torus with density ρ=e^5z

Appendix IV

Source

Fig 1a: https://medium.com/search?q=design
Fig 1 a: https://en.wikipedia.org/wiki/Bhabha_scattering
Fig 1 b: https://en.wikipedia.org/wiki/Bhabha_scattering