{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Calculate the Fourier coefficient $v_2$ from two-particle correlations "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 76,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import math\n",
    "from itertools import combinations"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Read data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "metadata": {},
   "outputs": [],
   "source": [
    "event, phi = np.loadtxt(\"dndphi_events.csv\", delimiter=',', skiprows=1, unpack=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define function that calculates $v_2$ for a given event\n",
    "\n",
    "One can use [``itertools.combinations``](https://docs.python.org/3/library/itertools.html#itertools.combinations) to get all pairs for a given 1d array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "def v2(phi_vals):\n",
    "    '''\n",
    "    phi_vals: 1d numpy array with phi values\n",
    "    return v2\n",
    "    '''\n",
    "    c=list(combinations(phi_vals, 2))\n",
    "    v2=0\n",
    "    for x in c:\n",
    "        v2 = v2 + math.cos(2*(x[0]-x[1]))/len(c)\n",
    "    return v2\n",
    "    # cumulants method, averaging the cosine of twice the difference of the pairs from a given event"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Remember: $<v_2^2> = <<e^{i 2 (\\varphi_1-\\varphi_2)}>>$."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Loop over all events and determine $v_2$ averaged over all events"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {},
   "outputs": [],
   "source": [
    "v2val = np.array([]) # array with v2 values for each event, can use numpy.append() to append a value\n",
    "\n",
    "nevt = 100\n",
    "v22=0\n",
    "for i in range(nevt):\n",
    "    phi_vals = phi[event == i]\n",
    "    v22 = v22 + v2(phi_vals)/nevt\n",
    "# average of the v2s for all the events"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We just calculated $<v_2^2>$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.1462940675740694\n"
     ]
    }
   ],
   "source": [
    "v2f=math.sqrt(v22)\n",
    "print(v2f)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our final estimate (ignoring non-flow effects) for $v_2$ is then $v_2= 0.146$."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}