Lab2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# MOwNiT Laboratorium 2\n",
    "## Zadanie 1\n",
    "Napisz w dowolnym języku zwyczajną (wolną) dyskretną transformatę Fouriera."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Podstawowy wzór na dyskretną transformatę Fouriera prezentuje się w następujący sposób:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$X_k = \\sum_{n=0}^{N-1}x_n \\cdot e^{-i2\\pi kn/N}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Można jednak przedstawić go w wygodniejszej postaci macierzowej:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\vec{X} = M\\vec{x},\\;\\;\\textrm{gdzie}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\vec{x} =\n",
    "\\begin{bmatrix}\n",
    "x_0 \\\\\n",
    "x_1 \\\\\n",
    "\\vdots \\\\\n",
    "x_{N-1}\n",
    "\\end{bmatrix},     \n",
    "\\;\\;\\vec{X} =\n",
    "\\begin{bmatrix}\n",
    "X_0 \\\\\n",
    "X_1 \\\\\n",
    "\\vdots \\\\\n",
    "X_{N-1}\n",
    "\\end{bmatrix},\n",
    "\\;\\;M = \n",
    "\\begin{bmatrix}\n",
    "e^{-i2\\pi kn/N}\n",
    "\\end{bmatrix}_{k\\times n}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "def dft(x):\n",
    "    N = x.shape[0]\n",
    "    n = np.arange(N)\n",
    "    k = n.reshape((N, 1))\n",
    "    M = np.exp(-2j * np.pi * k * n/N)   \n",
    "    return  np.dot(M, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Aby upewnić się, że funkcja działa poprawnie, można porównać jej wynik z transformatą liczoną przez numpy:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = np.random.random(2048)\n",
    "np.allclose(dft(x), np.fft.fft(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Algorytm działa poprawnie, ale jest niewydajny - złożoność O(N²), co dla dużych N daje niesatysfakcjonujący czas wykonania."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Zadanie 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wykorzystaj implementację z zadania 1. do napisania szybkiej wersji transformaty (używając pomysłu z algorytmu Cooleya-Tukeya)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Algorytm Cooleya-Tukeya polega na wykorzystaniu metody divide and conquer do podzielenia problemu na mniejsze podproblemy i zredukowaniu obliczeń w każdym z nich o połowę (dzięki symetrii DFT). Warunkiem stopu w poniższej implementacji jest zmniejszenie wejściowego wektora na tyle, że użycie zwykłej wersji DFT jest korzystniejsze. W tym przypadku empirycznie wyznaczyłem, że dzieje się to dla N <= 16."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def fft(x):\n",
    "    N = x.shape[0]\n",
    "    if not np.log2(N).is_integer():\n",
    "        raise ValueError('Length of x must be a power of two')\n",
    "    if N <= 16:\n",
    "        return dft(x)\n",
    "    \n",
    "    n = np.arange(N)\n",
    "    even = fft(x[::2])\n",
    "    odd = fft(x[1::2])\n",
    "    exponent = np.exp(-2j * np.pi * n / N)\n",
    "    return np.concatenate([even + exponent[:N//2] * odd, even + exponent[N//2:] * odd])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Aby upewnić się, że funkcja działa poprawnie, można porównać jej wynik z transformatą liczoną przez numpy:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.allclose(fft(x), np.fft.fft(x))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Algorytm Cooleya-Tukeya pozwala zmniejszyć złożoność z O(N²) do O(NlogN). Porównanie czasów obu wersji:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "249 ms ± 8.16 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n",
      "5.24 ms ± 198 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n"
     ]
    }
   ],
   "source": [
    "%timeit dft(x)\n",
    "%timeit fft(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Zadanie 3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Przetestuj implementację z zadania 2. do wykonania analizy szeregu czasowego:\n",
    "    1. Znajdź dane przedstawiające jakiś szereg czasowy.\n",
    "    2. Załaduj je do programu.\n",
    "    3. Zobacz, czy wykonanie analizy Fouriera na tych danych ma sens -- być może trzeba pogrupować je w równe odstępy.\n",
    "    4. Narysuj wykres częstotliwości i postaraj się opisać jaka jest główna składowa."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Szereg czasowy:\n",
    "https://www.kaggle.com/robikscube/hourly-energy-consumption#AEP_hourly.csv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Datetime     object\n",
       "AEP_MW      float64\n",
       "dtype: object"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "\n",
    "df = pd.read_csv('resources/AEP_hourly.csv')\n",
    "df.dtypes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dataframe składa się z dwóch kolumn:\n",
    " - Daty danego pomiaru (object, czyli string w Pandas)\n",
    " - Konsumpcji energii w megawatach (float64)\n",
    " \n",
    "Po załadowaniu i zbadaniu danych można zająć się ich analizą:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3sAAAFICAYAAADpvOfJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3XmYpGV97//3t/fu2Vd2GBAUARFkUaMmLqi4QY4hJ/i79GhiYnJO0CTGhURDlByTEI0eEzFGjTFxQyERERAQBQVFZUb2fYCBGWD2tfeuqu/vj6ppe2Z6qZnp6uqueb+ua66p56n7qfpWzT3V9en7fu4nMhNJkiRJUmNpqncBkiRJkqTJZ9iTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JMkSZKkBjQjw15EfCki1kfEvVW0PTIiboqIOyLi7oh4/VTUKEmSJEn1NCPDHvBl4Owq234Y+FZmngqcD3y2VkVJkiRJ0nQxI8NeZv4Y2DxyX0Q8KyKui4gVEXFLRBy/szkwt3J7HvD0FJYqSZIkSXXRUu8CJtHngT/KzEci4oWUR/BeCXwEuCEi3g3MAs6qX4mSJEmSNDUaIuxFxGzg14DLI2Ln7vbK328BvpyZ/xgRLwa+EhEnZWapDqVKkiRJ0pRoiLBHeTrq1sw8ZZT73knl/L7MvC0iOoDFwPoprE+SJEmSptSMPGdvd5m5HXg8In4bIMqeX7n7SeBVlf3PBTqADXUpVJIkSZKmSGRmvWvYaxHxDeDllEfo1gF/DfwQ+BfgEKAVuCwzL46IE4AvALMpL9bygcy8oR51S5IkSdJUmZFhT5IkSZI0voaYxilJkiRJ2pVhT5IkSZIa0IxbjXPx4sW5bNmyepchSZIkSXWxYsWKjZm5ZKJ2My7sLVu2jOXLl9e7DEmSJEmqi4h4opp2TuOUJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNjTfslMbnlkA6VS1rsUSZIkSSMY9rRfvnfvWt72b7/gP29bVe9SJEmSJI1g2NN+eXprHwBPbu6rcyWSJEmSRjLsSZIkSVIDMuxJkiRJUgMy7EmSJElSAzLsSZIkSVIDMuxJkiRJUgMy7EmSJElSAzLsSZIkSVIDMuxJkiRJUgMy7EmSJElSAzLsSZIkSVIDMuxJkiRJUgMy7EmSJElSAzLsSZIkSVIDMuxpUiRZ7xIkSZIkjWDY036JiHqXIEmSJGkUhj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPUmSJElqQIY9SZIkSWpANQ17EXF2RDwUESsj4sJR7j8yIm6KiDsi4u6IeH0t65EkSZKkA0XNwl5ENAOXAq8DTgDeEhEn7Nbsw8C3MvNU4Hzgs7WqR7WVWe8KJEmSJI1Uy5G9M4GVmflYZg4ClwHn7tYmgbmV2/OAp2tYj2og6l2AJEmSpFG11PCxDwNWj9heA7xwtzYfAW6IiHcDs4CzaliPJEmSJB0wajmyN9qgz+6T/d4CfDkzDwdeD3wlIvaoKSLeFRHLI2L5hg0balCqJEmSJDWWWoa9NcARI7YPZ89pmu8EvgWQmbcBHcDi3R8oMz+fmadn5ulLliypUbmSJEmS1DhqGfZuB46LiKMjoo3yAixX7dbmSeBVABHxXMphz6E7SZIkSdpPNQt7mVkALgCuBx6gvOrmfRFxcUScU2n258AfRMRdwDeAd2S6rqMkSZIk7a9aLtBCZl4LXLvbvotG3L4feEkta5AkSZKkA1FNL6ouSZIkSaoPw94kyEzuXrO13mVIkiRJ0jDD3iT46s+f5JzP/ISbHlpf71IkSZIkCTDsTYqH1+4AYPXm3jpXIkmSJEllhj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPe2XiHpXIEmSJGk0hj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGvUmQZL1LkCRJkqRdGPYmkdcXlyRJkjRdGPYkSZIkqQEZ9iRJkiSpARn2JEmSJKkBGfYkSZIkqQEZ9jQpMl2RVJIkSZpODHuSJEmS1IAMe5oUEV54QpIkSZpODHuSJEmS1IAMe5IkSZLUgAx7kiRJktSADHuSJEmS1IAMe5PAqw5IkiRJmm4Me5PJFSklSZIkTROGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPUmSJElqQIa9yeSynJIkSZKmCcPeJHARTkmSJEnTjWFPkyId1ZQkSZKmFcOe9ouDmpIkSdL0ZNiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JsELkQpSZIkabox7E0mL7gnSZIkaZqoadiLiLMj4qGIWBkRF47R5n9GxP0RcV9EfL2W9UiSJEnSgaKlVg8cEc3ApcCrgTXA7RFxVWbeP6LNccBfAC/JzC0RsbRW9UiSJEnSgaSWI3tnAisz87HMHAQuA87drc0fAJdm5haAzFxfw3okSZIk6YBRy7B3GLB6xPaayr6Rng08OyJ+EhE/i4izR3ugiHhXRCyPiOUbNmyoUbmSJEmS1DhqGfZGW61k93UrW4DjgJcDbwG+GBHz9zgo8/OZeXpmnr5kyZJJL1SSJEmSGk0tw94a4IgR24cDT4/S5juZOZSZjwMPUQ5/kiRJkqT9UMuwdztwXEQcHRFtwPnAVbu1uRJ4BUBELKY8rfOxGtYkSZIkSQeEmoW9zCwAFwDXAw8A38rM+yLi4og4p9LsemBTRNwP3AS8PzM31aomSZIkSTpQ1OzSCwCZeS1w7W77LhpxO4H3Vv5oBtv9ZExJkiRJ9VXTi6qr8UWMtg6PJEmSpHoz7EmSJElSAzLsSZIkSVIDMuxNAs9XkyRJkjTdGPYmkWevSZIkSZouDHuSJEmS1IAMe5IkSZLUgAx7kiRJktSAqgp7UfbWiLiosn1kRJxZ29IkSZIkSfuq2pG9zwIvBt5S2d4BXFqTiiRJkiRJ+62lynYvzMwXRMQdAJm5JSLaaliXJEmSJGk/VDuyNxQRzVQuKRcRS4BSzaqSJEmSJO2XasPePwHfBpZGxMeAW4G/rVlVkiRJkqT9UtU0zsz8WkSsAF5F+drhv5mZD9S0MkmSJEnSPhs37EXEwhGb64FvjLwvMzfXqrCZKOtdQB3lgfziJUmSpGloopG9FZQzTABHAlsqt+cDTwJH17S6GSLqXUAdxYH84iVJkqRpbNxz9jLz6Mw8BrgeeFNmLs7MRcAbgf+eigIlSZIkSXuv2gVazsjMa3duZOb3gN+oTUkzjzMYJUmSJE031V5nb2NEfBj4KuVs81ZgU82qmqGc0ShJkiRpuqh2ZO8twBLKl1+4Elha2SdJkiRJmoaqvfTCZuBPalyLJEmSJGmSVBX2IuImRjk1LTNfOekVSZIkSZL2W7Xn7L1vxO0O4LeAwuSXI0mSJEmaDNVO41yx266fRMSPalCPJEmSJGkSVDuNc+GIzSbgNODgmlQkSZIkSdpv1U7jXEH5nL2gPH3zceCdtSpKkiRJkrR/qg17z83M/pE7IqK9BvVIkiRJkiZBtdfZ++ko+26bzEIkSZIkSZNn3JG9iDgYOAzojIhTKU/jBJgLdNW4NkmSJEnSPppoGudrgXcAhwOfHLF/B/CXNapJM1DueRlGSZIkSXU0btjLzP8A/iMifisz/2uKatIMEhM3kSRJklQHE03jfGtmfhVYFhHv3f3+zPzkKIcdcPIAHtQ6gF+6JEmSNK1NNI1zVuXv2bUupBHEATzMFY7xSZIkSdPKRNM4/7Xy90enphxJkiRJ0mSo6jp7EbEE+ANg2chjMvP3alOWJEmSJGl/VHtR9e8AtwA3AsXalSNJkiRJmgzVhr2uzPxgTSuRJEmSJE2apirbXR0Rr69pJZIkSZKkSVNt2PsTyoGvLyK2R8SOiNhey8IkSZIkSfuuqmmcmTmn1oVIkiRJkiZPtatxvmCU3duAJzKzMLklSZIkSZL2V7ULtHwWeAFwT2X7ecBdwKKI+KPMvKEWxUmSJEmS9k215+ytAk7NzNMy8zTgFOBe4CzgH2pUmyRJkiRpH1Ub9o7PzPt2bmTm/ZTD32O1KUuSJEmStD+qDXsPRcS/RMRvVP58Fng4ItqBobEOioizI+KhiFgZEReO0+68iMiIOH0v658mst4F1F36HkiSJEnTSrVh7x3ASuBPgT8DHqvsGwJeMdoBEdEMXAq8DjgBeEtEnDBKuznAe4Cf713p08c3frEagCDqXMnUO/BesSRJkjQzVHvphT7gHyt/dtc9xmFnAit3TvWMiMuAc4H7d2v3N5TP+3tfNbVMZ45uSZIkSZouqhrZi4jjIuKKiLg/Ih7b+WeCww4DVo/YXlPZN/JxTwWOyMyr96pqSZIkSdK4qp3G+e/AvwAFytM2/xP4ygTHjDbDb3joKyKagE8Bfz7Rk0fEuyJieUQs37BhQ5UlS5IkSdKBq9qw15mZPwAiM5/IzI8Ar5zgmDXAESO2DweeHrE9BzgJuDkiVgEvAq4abZGWzPx8Zp6emacvWbKkypIlSZIk6cBV7UXV+ysjcY9ExAXAU8DSCY65HTguIo6utD8f+P923pmZ24DFO7cj4mbgfZm5vPryJUmSJEmjqXZk70+BLsqrZp4GvA14+3gHZGYBuAC4HngA+FZm3hcRF0fEOftesiRJkiRpItWuxnl75WY38LvVPnhmXgtcu9u+i8Zo+/JqH1eSJEmSNL5xw15EXDXe/ZnpCJ0kSZIkTUMTjey9mPLlE75B+aLnXkNbuyiWygus/mTlpjpXIkmSJGmkicLewcCrgbdQXlzlGuAbmXlfrQvTzLC5dwiAxzf21LkSSZIkSSONu0BLZhYz87rMfDvlSyOspHyphHdPSXWSJEmSpH0y4QItEdEOvIHy6N4y4J+A/65tWZIkSZKk/THRAi3/QfnC598DPpqZ905JVZIkSZKk/TLRyN7bgB7g2cB7IobXZwkgM3NuDWubccL1ayRJkiRNE+OGvcys9qLrOkAZbyVJkqTpyTAnSZIkSQ3IsKf9cufqrfUuQZIkSdIoDHvaLz96eEO9S5AkSZI0CsOeJEmSJDUgw94UWb25l3vWbKt3GZIkSZIOEBNeVF2T42X/cBMAq/7+DXWuRJIkSdKBwJE9SZIkSWpAhj1JkiRJakCGPUmSJElqQIY9SZIkSWpAhj1JkiRJakCGPUmSJElqQIa9SRRR7wokSZIkqcywJ0mSJEkNyLAn7Yfr71vLrY9srHcZkiRJ0h5a6l2ANJP94VdWALDq799Q50okSZKkXTmyN4l6B4v1LkGSJEmSAMPepBoqlupdgiRJkiQBhr1JlVnvCiRJkiSpzLA3iRLTniRJkqTpwbAnSZIkSQ3IsCdJkiRJDciwJ0mSJEkNyLA3iVygRZIkSdJ0YdiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZk2JMkSZKkBmTYkyRJkqQGZNiTJEmSpAZU07AXEWdHxEMRsTIiLhzl/vdGxP0RcXdE/CAijqplPbWWmfUuQZIkSZKAGoa9iGgGLgVeB5wAvCUiTtit2R3A6Zl5MnAF8A+1qkeSJEmSDiS1HNk7E1iZmY9l5iBwGXDuyAaZeVNm9lY2fwYcXsN6ai4i6l2CJEmSJAG1DXuHAatHbK+p7BvLO4Hv1bCemqtmGmffYHEKKpEkSZJ0oGup4WOPNsw1ahqKiLcCpwO/Mcb97wLeBXDkkUdOVn11MVQq0UlzvcuQJEmS1OBqObK3BjhixPbhwNO7N4qIs4APAedk5sBoD5SZn8/M0zPz9CVLltSkWEmSJElqJLUMe7cDx0XE0RHRBpwPXDWyQUScCvwr5aC3voa1TBue1SdJkiRpKtQs7GVmAbgAuB54APhWZt4XERdHxDmVZh8HZgOXR8SdEXHVGA8nSZIkSdoLtTxnj8y8Frh2t30Xjbh9Vi2ffzpyxU5JkiRJU6GmF1WXJEmSJNWHYW+KOa4nSZIkaSoY9iZRFZfZkyRJkqQpYdibYjPplL3HNnSz7MJruOPJLfUuRZIkSdJeMuxpTDc/tAGA79y5x+URJUmSJE1zhr1JVM0szvCsPUmSJElTwLA3iRrtnL2ZNOVUkiRJ0q4Me5OomnA0EwNUNlqKlSRJkg4Ahr1p5Mo7nuKMj91IsTQ9wtUMzKWSJEmSKgx7k2h/B8A+9O172LBjgL6h4uQUNEmmR/SUJEmStDcMe5PoUzc+PCmPM12mTcZMnHMqSZIkCTDsTSvTNVxNk+wpSZIkaS8Y9qahqchWpVKyfkf/uG2mafaUJEmSVAXD3jSyM1tNxUjap3/wCGd+7Ac8vbWv9k8mSZIkacoZ9qaTKRxJu/mh9QCs3zEwdU8qSZIkacoY9ibZH3/tl/v/IFMwsrc3T5GuxylJkiTNOIa9SXbNPc/s87HD0zinMFyNN5joKXuSJEnSzGXYm0ZcjVOSJEnSZDHsTUNTEa52Pse4+XKahk9JkiRJEzPsTSNTma12ThWNKiZrOrAnSZIkzTyGvQNUNaOHjutJkiRJM5dhbxqaypG08UYTN/cMAp6zJ0mSJM1Ehr1pZCpH0qoJcJ/8/sMA3Pf0thpXI0mSJGmyGfamoZyCobS9eYZC0aE9SZIkaaYx7E0jU3nphZ2BspqndFFOSZIkaeYx7E2x8Qbt+oeK5TZTVAtUtxqnJEmSpJnHsFcD/+/Gh/fpuN7BStibZrMmp1s9kiRJkiZm2KuB/3fjI2PeN1gsTWElk+P+Z7bXuwRJkiRJe8mwVyNrtvSOuv+KFWsmPDa9jLkkSZKk/WTYq5GXXnITxdKeoa0wTUb2Hly7A4D//uXE4VOSJEnSzGPYq6Eb7lu7x76qxuymcGDvvqedoilJkiQ1opZ6F9DIBosl+oeKPLqhe3hfNYudTOUkTqeMSpIkSY3JsFdDV97xFN+/fx1X3/3M8L7pFq5+9tjmepcgSZIkqQYMezV000MbWDKnfZd9A0PT45w9SZIkSY3Nc/ZqbMOOgV22n9nWN+ExXtdOkiRJ0v4y7E0xg5wkSZKkqWDYm2LFKtLeeOf1bewe4Mb7101mSRrDh759D5+4/qF6lyFJkiTtE8PeFPvvXz61X8e//Uu/4Pf/czm9g4Ux22zsHuBff/Qo6TDifvnaz5/kMzetrHcZ01Zmsm57f73LkCRJ0hgMe5PgojeeMKmPN15Ge3xjDwCjXK992Hu/dRd/970HuXvNtkmtS2Pb3j9U7xKm3BdueYwX/u0Pdrm0yO56BwtcsWJNQ/7i4Q3/dAu/eelP6l2GJEnSmAx7k+BVz126V+2XXXgNpXHS2v5+Ld7eVw4ehfES4QQeG+cLvPbUgFlmQrc8shGANVvGXnTow1fey/suv4tfPrllqsqaMvc9vZ07V2+tdxmSJEljMuzVyTF/eS2PrNsx6n03Pbh+zON6B4sAlKpIFxH7VttjG7p55T/+aK+PW7/DKX0Hkq295V8qNI3Tz3YGwcHCAZiGG9jyVZv58JX31LsMSZI0AcNeHb36Uz/m/3xtBT0Du55/d9ujmyY89pmtYwerib5WDxVHv9Zf/1CRY/7iGr546+N73LfswmtYub482rduez//52sr+PHDG3Zpc+bHfjDBMx+Ylq/azD0NOKX2nqfKrykYO+3t7GttLX7UNJLzPncbX/3Zk+POUJAkSfXnRdXr7Np71vKCIxfw+y87ZnjfNfc8w6VjtG+K8vl687tax3zMu8aZWnbvU9t44z/fusu+ZRdew9+9+Xnc+eRWSglf//mTox571id3He279p61e7TpGyzS1ATtLc30DxU5/q+u44JXHMv7XvucMWvaaed5XbGvQ5JTrH+oOHx7xRObeeXxB43a7rzP3QbAqr9/w5TUNdWe3Nw75n3DYa+5scLeQ2tHH5U/0BQzaRon7Deqh9ft4KhFXbS3NNe7FEmSxtVY38BmqP97zQO87tO37LIvM9nUPbBH2/ldbcDYo3MjffamR/n5Y5tYvmozV6xYw62PbNwj6O30F/99D99cvnofqt/Vcy+6jud8+DqWXXgNx//VdQB85qaVrN7cy09Xls/x+tPL7uDUi2/Y5bj+oSJH/8W1fPoHj+x3DdW49KaV+32+1Wl/8/3h25cvX7O/Jc1Yf/ntsafzbe4eBKClefxA8NOVGylOo1GigUJx3Hp6xlkN90BSKE6ff7Opsn5HP6/51I/5qyvvrXcpkiRNyJG9aeKBZ7bvsn3pTSv5xA0P8/0/+3WOO2jO8P6dUz4fWd/N4Qu6uHP1Vp5/+Dz6h0pcdvuT/KQSqABufGAdNz4wPa7J97J/uGmPfS/82xvZ2jvEXX/9muHzFL/801X86VnPrnk9H7/+IT5+/UP7PNr26IZuegZ/NbL3vXv3HOUUPL2tPN34kuse5Mu/e+aoba6++2ku+PodfPDs4/nfL3/WVJY3pud8+DrOPvFgPve200a9v2mGjD7vi5Xrd/CTlZt4+68tm7Dtuu39LFs8q/ZFTSM7+sufwbevarxFhyRJjaemYS8izgY+DTQDX8zMv9/t/nbgP4HTgE3A72TmqlrWNFN84oaHgfJ5faP53X+/fSrLqYl128sjlztHAKG86MdXf/YEb3jeISyY1bbHMRu7B5jf2UrLfkwLLIwxKjpYKDFULPGhb9/DlXc+Pbz/j7/2S37vpcs47aiFbOoeYM2WPs4dZcn9cz9zK9+54KXD232DRb7xi19NiS2VkqYxVjNZv6Of7X1DHLt0zqj3Q3m0t5TQPN6KKNPUzQ9tGPO+C75+BwBf/dkTUxL21m7rZ8mc9jHfx9WVaanX3Td2gP+9L1f3/+/Sm1ZyxrKFnHn0wr0vtAbWbOmlZ6DIcw4eu5+d9cnyZ87/evFRo06pHjniefHV9/Old5wx+YVOY9sqqx3vvAzOgWZLzyC9Q0UOm985Zpud53KO9Xk3kz2xqYcjF3bNmNMNJKlmYS8imoFLgVcDa4DbI+KqzLx/RLN3Alsy89iIOB+4BPidWtVUK+P90NPe+/CV9/Lh3aZIvePXlvHln64a3v7R+1/On33zTt79yuN4aN0O1mzp5YNnH8+jG3qGr302p6Nl+LfwY1l24TUT1nPNPc9wzT3PTNjurjXbxn28Y/7y2lH3L57dxsbKdMfvXvBSjlzUxZz2lj2+KJ380RvY0V/Yr3P/MpPr7l1Lz2CRN558CEPFEt/4xZPcuXor196zll9/9pI9Ft4BOOmwubztRUfxomMW8V8r1vCm5x+6y/2XXPcgv33a4cPXf7xixZo9Lkly/udv4z2vPI47Vm/ltScexJELZ7G1d3D4/qe2jn0Jh8myZksvL73kJk48dC7XvOdlZCYbugdob25mbmcLA4XSLn3vi7c8xouftYgTD51HZrJ2ez8LZ7WxuedXdV9999O84jlLh6eqlkrl1XJbmoOPX/8QAPdf/FqGCsm8Uc61/chV9/HcQ+bwO2ccWeNXDy+9pDzC/ujfvn7UsDtywZUVT2zhpMPm0dG663lpt6/aPHx7tHMXewcLdLQ0T8sv+uP9wgXKtZ9w0fV87q0v4OyTDhm1zZs/+9MJn2db3xDP/+gNvPfVz+Y9rzpun+udbh7f2MMrPnEzUO7TXW2jf4V49ad+xKMbesb9rOoZKDBQKLFwlF/q1UPvYIHP/HAlHzj7+DHbfPeup3n3N+7g0+efwrmnHDZqm1IpuXzFat78gsNpnUHnKvcPFWltbhrzl2ClUvKdu56is7WFs086eMzH+f796+hqa+Ylxy6uVamTrlRKIvZ/vYC+wSI/enj9mJ8d09XK9d08a8ms/X79T23tY2FXG51tM+dc5p6BAq3NTeMuIrduez9L57TP6F/wRK0udhwRLwY+kpmvrWz/BUBm/t2INtdX2twWES3AWmBJjlPU6aefnsuXL69JzfujmtAgSZPhuKWzeWS918I8EH3w7ONZWhmZTpI/++Zd9S5JNfam5x/K4Qs6+ZebH+Xkw+exqXuQ97/2OUTAB664m4HCxOfwz2TPOWgOqzb1MFAo8arjl/KDcS5P1cg+ePbxXHLdg5x94sHc9tgmugcKfOScExvy/OFli7p4Zlv/Ln374nNP5KLv3MebTz2MuZ2tPPugOfQMFPjYtQ9MSU3TcZG9iFiRmadP2K6GYe884OzM/P3K9tuAF2bmBSPa3Ftps6ay/WilzcbdHutdwLsAjjzyyNOeeOKJmtS8P9Zu6+dT33+Y79z1FEFwyLwO3njyITyxuZfvjJgS2Ii62pqHr/8nSZIkNZKZHPZqec7eaOOduyfLatqQmZ8HPg/lkb39L23yHTyvg0vOO5lLzjt5j/s+ff6pu2wXiiWaIoanE+08hyxhl2kfmTk8bJyZFEpJU8Qu0ywGCyWam8r7SqVkR3+Brvbm4ccplZL+QpGmCNqam2hqil0edyxDxRKlzD2WFq/m8giFYoliJs0RRARNI6ZHDBSKtDSVF2svZhKwy/l3mclQMYcv1L2tb4hZ7S20tzSxtXeIrLzm3sEC3QMFFnS10TdUpK25id7KZR9am5tobWpi4ew2mqJ8jlGhmMztbCUz6Rko0tIctLc0MVAosblncPgcrp3n7c3paB0+R66USWtzE4OF0vBQ/0DlPW2vbGeWX8/G7gG6WluY3dFCZtLcFMOvfahYYqBQorU52N5XoJTJnI4WtvUN0dnazKz2FjZ1D7Kjv/yad/QX2NE/RP9QiU09A3S2NlMoJas399La3MTs9hbaW5uG+0DPQIHmpibmdbYSAbPaW1g0q204iDc3letsaW7isPmddA8U6O4vcNTiLkqlZHZ7C9v7C7Q0B9t6h+hobWZORws9AwVKWZ7mtPPfakf/EAu62lg8u51Vm3oolsrv0ay2ZnYMFOgfKrJkTjtzO1ppbgp6B4ps6hmYw8w/AAARc0lEQVTgsAWdbO4ZJCi/1xGwqWeQkw6dy9Nb+1m1qYfugQJNEcxub+GZbX10tbUwUCjSP1SiWCpx+6otLJrdRqFYfn8PX9BJR2szCXS1NvP01j4Wz2mnb7DI4jnt5ecBOlqb6Wht4sRD5/HMtj7mdJTfp47WZgaGinS1tbBgViubewZZNKudNVt6SWBupV3vQLHcf5ubKJaSZYu6GCiUWLe9n9bmJgqlpGegwJyO8sdqECyZ0872/iH6h4o8uHYHyxbNYs2WXg6Z10lTE2zuGeSQeZ30Dha4/+ntlDLZsGOAoWIyWCyxrW+IoxZ2USglT2/tIwIOnd/JS49dzMr13TyxqZeutmYWzmqjuSnY2D1AUwTb+wvM7WhhyZx2SpkMFko8uqGH1554EEcs6GLl+m66BwocNLeDZ7b1cfzBc+keKLBwVhtzOloYLJR4onLuYrHy3LPaWlg6t53HN/Zw7NLZHLN4Nvc/s40tPUMsmt1Gz0CR/qEiS+e2s3BWG3ev2UZXWzPtLc2s3d43/H9oQVcb/UPl93GoWKJ/qMiW3iFam4ND5nUOX7/xqIVdlDL52WObOWbJLA6d38G67QMsmlU+vm+oyFGLZvFfK9YMXwZhVnsLQ8Xy/+FZ7S20NAU7+gts7hnk4HkdnHrkfJoiWLWxh8Wz22ltaWJr7yD9Q0Wef8R8ZrW10NrcxJ2rt9Le0sTczpby50oET23tY9GsNjLhmCWzeGR9N+u397NgVhuDlc+RQ+d3Mr+zlf5CkWcfNIe12/rpGSzSO1CguSnoL5ToGywwv6uN9dv7efGzFnHHk1tpa2li1cZe2lvL/7+WzG6nWJlyWiqVpw+fdNhcOlqaOWR+5/DnWinLUyuborxo1JbeIeZ1trK9b4jugQJ9g0Vmtbdw+ILO4c+ewUKJwxd0Mb+rlSc39zK3o5UFs9oYKpQ4dH4nG7oHKJWSoWKJ1Vv66GxtprU5eGJTLy3NwdknHUxXaws/fHAda7cPsKCrlWOWzGZL7yBtzU0csbCLQ+Z1sKV3kHuf2saGHQMsmNXG9r4CSfLcg+cyWPl3n9vRyobuAfoGi2ztHWJWezOZMFgs0TNQoHewyOz2FuZ1ttJU6d/tLU0cf/BcBgpFBgvlx3liUy8LZ7cxv7ON5x02j829gzy5uZfu/vLn7CHzOjhqURePrOtmoPI5vmR2+TN/Q/cALU1Be0szpxw5n7tXb6WrvYVntvZRyvI5eotnl/8fPe+webS1NHHTQ+uZ19nK1t4hFs1uZ0f/EH2DRV54zEIKxWTR7DYe21A+bqBQYsOOftpbmuloa2ZT9wBPbeljqFjiBUctoFBMkvJUuiVz2glg6dx21m7rp7213P6kQ+fR3Bxs3DHAsUtns277AM9aUl4Y6a4121i7rY/BQoktvUM0NwXHHzyHR9Z3M7+zldkdLRy5sGv4HNP5XW20NgcbdpTf90Pmd3Lw3A5mtTezYccAj23sYdGsNuZ0tLK9f4ihYonZ7S2s2dLHGcsW0jtYYE2l/g07Blg6p52ewSIdrc089+A5zOtq5YFndrC9r/y50D1QoDmC7srick0RLJzVxkCh/O+2bvsAczpaeM0JB/HVnz/BWc89iKvvfobjls6GgFMOn8/SuR2s2dLLwFCJYw+azQPPbGd+Zxvrd/TT1dbM1t4hHlnfzXMOmkNnWzPdAwUOntsx/LNke/8Qczta6R0sMLu9hcFiiR39BQaGipxw6Dxam4OmCLb0DtI9UG5TKCVrt5U/2xd0tbJ6Sy/HLJ7Nhu4BtvUNcfDcDvqGyp958ztbOXheB5lwwqFzueepbQwMlWhqgqFCMlQqf+eb09FCKcuXIvrlk1uY095Cf6HI/M42Fsxq476nt3Hw3A6am4JNPYMcNKe9PFU64Kktfbz6hINYvmozxx8yl3Xb+ymWklOPXMDda7Zy/zPbmdvRyvyuVrb3FVi7rY95na10tbewamMPLzxmIYOFEsVSZVri3Pbhn0fzu9r45ZNbOHrRLHoHi2zpHaSlqYlD5newtXeQ1uYmOlubh2eT7Pw5saVniPld5de+ZksfyxZ1VX5OBY9u6ObFxyxi/Y7y+3Xkwi629g0xu72ZTd2DbOwe5LAFnazb1j/8ni+d286Tlf/Li2e1k5S/9w4MlT/TX3viwdz/zHbaW8rff5YtnkWxlDy+sYf+oSKPb+zh9GULeXjdDn708AaeU1ncsKW5/P03E447aPbwd6Kd36O6Bwps6xtioFBkw44B+odKHLt0Ns1NsGpj+XvA21501JjfeWcCp3FKkiRJ0gxS7cheLc8evh04LiKOjog24Hzgqt3aXAW8vXL7POCH4wU9SZIkSVJ1ajaNMzMLEXEBcD3lSy98KTPvi4iLgeWZeRXwb8BXImIlsJlyIJQkSZIk7aeaXmcvM68Frt1t30UjbvcDv13LGiRJkiTpQDRzLgIjSZIkSaqaYU+SJEmSGpBhT5IkSZIakGFPkiRJkhqQYU+SJEmSGpBhT5IkSZIakGFPkiRJkhpQZGa9a9grEbEBeKLedYxiMbCx3kWoYdm/VEv2L9WafUy1ZP9SLU3X/nVUZi6ZqNGMC3vTVUQsz8zT612HGpP9S7Vk/1Kt2cdUS/Yv1dJM719O45QkSZKkBmTYkyRJkqQGZNibPJ+vdwFqaPYv1ZL9S7VmH1Mt2b9USzO6f3nOniRJkiQ1IEf2JEmSJKkBGfb2UkScHREPRcTKiLhwlPvbI+Kblft/HhHLpr5KzVRV9K9fj4hfRkQhIs6rR42auaroX++NiPsj4u6I+EFEHFWPOjUzVdG//igi7omIOyPi1og4oR51auaaqI+NaHdeRGREzNgVFDX1qvgMe0dEbKh8ht0ZEb9fjzr3lmFvL0REM3Ap8DrgBOAto/yweiewJTOPBT4FXDK1VWqmqrJ/PQm8A/j61Fanma7K/nUHcHpmngxcAfzD1FapmarK/vX1zHxeZp5CuW99corL1AxWZR8jIuYA7wF+PrUVaiartn8B38zMUyp/vjilRe4jw97eORNYmZmPZeYgcBlw7m5tzgX+o3L7CuBVERFTWKNmrgn7V2auysy7gVI9CtSMVk3/uikzeyubPwMOn+IaNXNV07+2j9icBbhogPZGNd/BAP6G8i8T+qeyOM141favGcewt3cOA1aP2F5T2Tdqm8wsANuARVNSnWa6avqXtK/2tn+9E/heTStSI6mqf0XEH0fEo5S/jL9nimpTY5iwj0XEqcARmXn1VBamhlDtz8jfqpzqcEVEHDE1pe0fw97eGW2EbvffTFbTRhqNfUe1VHX/ioi3AqcDH69pRWokVfWvzLw0M58FfBD4cM2rUiMZt49FRBPl02f+fMoqUiOp5jPsu8CyyqkON/KrmXzTmmFv76wBRqb4w4Gnx2oTES3APGDzlFSnma6a/iXtq6r6V0ScBXwIOCczB6aoNs18e/v5dRnwmzWtSI1moj42BzgJuDkiVgEvAq5ykRZVacLPsMzcNOLn4heA06aotv1i2Ns7twPHRcTREdEGnA9ctVubq4C3V26fB/wwvZihqlNN/5L21YT9qzIF6l8pB731dahRM1c1/eu4EZtvAB6Zwvo0843bxzJzW2YuzsxlmbmM8nnH52Tm8vqUqxmmms+wQ0ZsngM8MIX17bOWehcwk2RmISIuAK4HmoEvZeZ9EXExsDwzrwL+DfhKRKykPKJ3fv0q1kxSTf+KiDOAbwMLgDdFxEcz88Q6lq0ZosrPr48Ds4HLK+tKPZmZ59StaM0YVfavCyojx0PAFn71i1FpQlX2MWmfVNm/3hMR5wAFyt/x31G3gvdCOOgkSZIkSY3HaZySJEmS1IAMe5IkSZLUgAx7kiRJktSADHuSJEmS1IAMe5IkSZLUgAx7kqRhEVGMiDtH/FlW75omS0ScHhH/VO86aiEiXh4RvzZi+zcj4oQqjmuNiBURsSwi7q1tlZKkqeZ19iRJI/Vl5ilj3RkRLZlZmMqCJkvl4sqNeoHllwPdwE8r278JXA3cP8FxLx1xzKSbyf1FkhqBI3uSpHFFxDsi4vKI+C5wQ2Xf+yPi9oi4OyI+OqLthyLioYi4MSK+ERHvq+y/OSJOr9xeHBGrKrebI+LjIx7rDyv7X1455oqIeDAivhaVK71HxBkR8dOIuCsifhERcyLilog4ZUQdP4mIk3d7HS+PiKsrtz8SEV+qPMdjEfGeMV57d0RcUhn9ujEizhxxzDkTvIamiPhsRNwXEVdHxLURcV7lvosq7e+NiM9H2dKIWFG5//kRkRFxZGX70Yjoiog3RcTPI+KOSj0HVUZf/wj4s8po7G8A5wAfr2w/KyLeExH3V+q7bMRLPBv4XuV2c0R8oVLvDRHRWXnuUyLiZ5Vjvx0RCyb4N92jv0iS6sOwJ0kaqXPEFM5vj9j/YuDtmfnKiHgNcBxwJnAKcFpE/HpEnAacD5wKvBk4o4rneyewLTPPqLT/g4g4unLfqcCfAicAxwAviYg24JvAn2Tm84GzgD7gi8A7ACLi2UB7Zt49wXMfD7y28jr+OiJaR2kzC7g5M08DdgD/F3g18D+Aiyd4DW8GlgHPA36f8nu402cy84zMPAnoBN6YmeuBjoiYC7yM8ijkyyLiKGB9ZvYCtwIvysxTgcuAD2TmKuBzwKcy85TM/BFwFfD+yvajwIXAqZl5MuVguNMrgJsrt48DLs3ME4GtwG9V9v8n8MHKsfcAfz3B+woj+ksVbSVJNeI0TknSSGNN4/x+Zm6u3H5N5c8dle3ZlIPCHODblVBCRFxVxfO9Bjh554gXMK/yWIPALzJzTeWx7qQcnLYBz2Tm7QCZub1y/+XAX0XE+4HfA75cxXNfk5kDwEBErAcOAtbs1mYQuK5y+x5gIDOHIuKeSj3jvYaXApdnZglYGxE3jXjcV0TEB4AuYCFwH/BdylMqXwL8OvC3lEfeArilctzhwDcj4hCgDXi8itcJcDfwtYi4ErgSICIOBTZnZm9l0PTxzLyz0n4FsCwi5gHzKwES4D+Ay6t4vpH9RZJUJ47sSZKq0TPidgB/Vxk1OiUzj83Mf6vcl2McX+BXP3M6dnusd494rKMzc+fUv4ER7YqUf0EZoz1HJWB+HzgX+J/A16t4TaM9/u6GMnPn85V2HlMJcDvbj/UaYrQnjYgO4LPAeZn5POAL/Oo9uYXyqN5RwHeA51MOjT+u3P/PlEcFnwf8Ibu+l+N5A3ApcBqwIiJagNcB149oU837MdJY/6awa3+RJNWJYU+StLeuB34vImYDRMRhEbGUciD5HxHRGRFzgDeNOGYV5aABcN5uj/W/d06hjIhnR8SscZ77QeDQiDij0n5OJbhAeSrnPwG3T/Go0liv4Vbgtyrn7h1EeREV+FUw2lh5D0e+Hz8G3go8UgmUm4HXAz+p3D8PeKpy++0jjttBeWR1j+2IaAKOyMybgA8A8ymPxo48X29UmbkN2BIRL6vsehuwc5RvFaP/m0qSpgmncUqS9kpm3hARzwVuq0z/6wbempm/jIhvAncCT/CrqYcAnwC+FRFvA344Yv8XKU+H/GWUH2wD5ZUkx3ruwYj4HeCfKwuI9FE+b687M1dExHbg3yfppVZrrNfwX8CrgHuBh4GfUz63b2tEfIHytNBVwO07HygzV1Xe050jebcCh2fmlsr2R4DLI+Ip4GfAzvMbvwtcERHnAu+mfD7fF6K88Mz5wL9VpmQG8CnKYfC4zHywitf3duBzEdEFPAb8bmX/WP+mkqRpIn41O0WSpMkTER+hHMI+MUXPdyjlxUaOr4yK1V1EzM7M7ohYBPwCeElmrp0Gdb2UckD/owkbS5JmLEf2JEkzXkT8L+BjwHunS9CruDoi5lNeTOVvpkPQA8jMWymPGkqSGpgje5IkSZLUgFygRZIkSZIakGFPkiRJkhqQYU+SJEmSGpBhT5IkSZIakGFPkiRJkhqQYU+SJEmSGtD/D16CBF/0XZuaAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 1080x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import fftpack\n",
    "\n",
    "n = len(df)\n",
    "\n",
    "X = fftpack.fft(df['AEP_MW'])\n",
    "freqs = fftpack.fftfreq(n)\n",
    "\n",
    "f, ax = plt.subplots()\n",
    "f.set_size_inches(15, 5)\n",
    "ax.plot(freqs[1:n//2], np.abs(X)[1:n//2])\n",
    "ax.set_xlabel('Frequency in megawatts/hour')\n",
    "ax.set_ylabel('Magnitude')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dominującą częstotliwością jest 0.042 MW/h, czyli około 1 MW/dzień. Możemy więc przypuszczać, że w danych występują wzorce powtarzające się na co dzień."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}