master-thesis/lab-windows/grid_frequency_spectra.ipynb
2020-02-20 19:21:09 +00:00

354 lines
13 KiB
Text

{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import csv\n",
"\n",
"import numpy as np\n",
"from matplotlib import pyplot as plt\n",
"import scipy.fftpack"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib widget"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"data = np.genfromtxt('data/Netzfrequenz_Sekundenwerte_2012_KW37.csv', delimiter=',')[1:,1:]"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "36f1f4d7970e41afa4737b6f63b7c449",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"[<matplotlib.lines.Line2D at 0x7f952a6e4580>]"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"fig, ax = plt.subplots()\n",
"ax.plot(data[:3600*24, 0])"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.02051102806199375"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"np.std(data[:,0])"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "d7fe0512f4254efeb15235a5617ef064",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"(1e-06, 0.5)"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Number of samplepoints\n",
"N = len(data[:,0])\n",
"# sample spacing\n",
"T = 1.0\n",
"x = np.linspace(0.0, N*T, N)\n",
"yf = scipy.fftpack.fft(data[:,0])\n",
"xf = np.linspace(0.0, 1.0/(2.0*T), N//2)\n",
"\n",
"yf = 2.0/N * np.abs(yf[:N//2])\n",
"\n",
"#yf = sum(yf[s::10] for s in range(10)) / 10\n",
"#xf = sum(xf[s::10] for s in range(10)) / 10\n",
"\n",
"fig, ax = plt.subplots()\n",
"ax.loglog(xf, yf)\n",
"ax.xaxis.set_major_formatter(plt.FuncFormatter(lambda x, _pos: f'{1/x:.1f}'))\n",
"ax.set_xlabel('T in s')\n",
"ax.set_ylabel('Amplitude Δf')\n",
"ax.grid()\n",
"ax.set_xlim([1/1000000, 0.5])"
]
},
{
"cell_type": "code",
"execution_count": 72,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"<ipython-input-72-51d3a7cc1678>:20: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).\n",
" fig, ax = plt.subplots()\n"
]
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "2991d932b113496a9135d569f9577abe",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"(5e-07, 0.02)"
]
},
"execution_count": 72,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Number of samplepoints\n",
"N = len(data[:,0])\n",
"# sample spacing\n",
"T = 1.0\n",
"x = np.linspace(0.0, N*T, N)\n",
"yf = scipy.fftpack.fft(data[:,0])\n",
"xf = np.linspace(0.0, 1.0/(2.0*T), N//2)\n",
"\n",
"yf = 2.0/N * np.abs(yf[:N//2])\n",
"\n",
"average_from = lambda val, start, average_width: np.hstack([val[:start], [ np.mean(val[i:i+average_width]) for i in range(start, len(val), average_width) ]])\n",
"\n",
"average_width = 20\n",
"average_start = 100\n",
"yf = average_from(yf, average_start, average_width)\n",
"xf = average_from(xf, average_start, average_width)\n",
"yf = average_from(yf, 300, average_width)\n",
"xf = average_from(xf, 300, average_width)\n",
"\n",
"fig, ax = plt.subplots()\n",
"ax.loglog(xf, yf)\n",
"ax.xaxis.set_major_formatter(plt.FuncFormatter(lambda x, _pos: f'{1/x:.1f}'))\n",
"ax.set_xlabel('T in s')\n",
"ax.set_ylabel('Amplitude Δf')\n",
"\n",
"for i, t in enumerate([45, 60, 600, 1200, 1800, 3600]):\n",
" ax.axvline(1/t, color='red', alpha=0.5)\n",
" ax.annotate(f'{t} s', xy=(1/t, 3e-3), xytext=(-15, 0), xycoords='data', textcoords='offset pixels', rotation=90)\n",
"#ax.text(1/60, 10,'60 s', ha='left')\n",
"ax.grid()\n",
"ax.set_xlim([1/60000, 0.5])\n",
"ax.set_ylim([5e-7, 2e-2])"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"ename": "TypeError",
"evalue": "object of type <class 'float'> cannot be safely interpreted as an integer.",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mTypeError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m~/safety-reset/lab-windows/env/lib/python3.8/site-packages/numpy/core/function_base.py\u001b[0m in \u001b[0;36mlinspace\u001b[0;34m(start, stop, num, endpoint, retstep, dtype, axis)\u001b[0m\n\u001b[1;32m 116\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 117\u001b[0;31m \u001b[0mnum\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0moperator\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mindex\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mnum\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 118\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;31mTypeError\u001b[0m: 'float' object cannot be interpreted as an integer",
"\nDuring handling of the above exception, another exception occurred:\n",
"\u001b[0;31mTypeError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m<ipython-input-7-75728c9461c4>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m 5\u001b[0m \u001b[0mys\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mconvolve\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mys\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mones\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m/\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmode\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'valid'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 6\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 7\u001b[0;31m \u001b[0mxs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mlinspace\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m1.0\u001b[0m\u001b[0;34m/\u001b[0m\u001b[0;36m2.0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdata\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m/\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 8\u001b[0m \u001b[0;31m#xs = np.linspace(len(data)/2, 1, len(data)/2)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 9\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m<__array_function__ internals>\u001b[0m in \u001b[0;36mlinspace\u001b[0;34m(*args, **kwargs)\u001b[0m\n",
"\u001b[0;32m~/safety-reset/lab-windows/env/lib/python3.8/site-packages/numpy/core/function_base.py\u001b[0m in \u001b[0;36mlinspace\u001b[0;34m(start, stop, num, endpoint, retstep, dtype, axis)\u001b[0m\n\u001b[1;32m 117\u001b[0m \u001b[0mnum\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0moperator\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mindex\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mnum\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 118\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 119\u001b[0;31m raise TypeError(\n\u001b[0m\u001b[1;32m 120\u001b[0m \u001b[0;34m\"object of type {} cannot be safely interpreted as an integer.\"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 121\u001b[0m .format(type(num)))\n",
"\u001b[0;31mTypeError\u001b[0m: object of type <class 'float'> cannot be safely interpreted as an integer."
]
}
],
"source": [
"ys = scipy.fftpack.fft(data[:,0])\n",
"ys = 2.0/len(data) * np.abs(ys[:len(data)//2])\n",
"s = 60\n",
"\n",
"ys = np.convolve(ys, np.ones((s,))/s, mode='valid')\n",
"\n",
"xs = np.linspace(0, 1.0/2.0, len(data)/2)\n",
"#xs = np.linspace(len(data)/2, 1, len(data)/2)\n",
"\n",
"fig, ax = plt.subplots()\n",
"ax.loglog(xs[s//2:-s//2+1], ys)\n",
"ax.xaxis.set_major_formatter(plt.FuncFormatter(lambda x, _pos: 1/x))\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"ys = scipy.fftpack.fft(data[:,0])\n",
"ys = 2.0/len(data) * np.abs(ys[:len(data)//2])\n",
"s = 1\n",
"\n",
"ys = np.convolve(ys, np.ones((s,))/s, mode='valid')\n",
"\n",
"xs = np.linspace(0, 1.0/2.0, len(data)/2)\n",
"#xs = np.linspace(len(data)/2, 1, len(data)/2)\n",
"\n",
"fig, ax = plt.subplots()\n",
"ax.loglog(xs[s//2:-s//2+1 if s > 1 else None], ys)\n",
"ax.xaxis.set_major_formatter(plt.FuncFormatter(lambda x, _pos: 1/x))\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"ys = scipy.fftpack.fft(data[:,0])\n",
"ys = 2.0/len(data) * np.abs(ys[:len(data)//2])\n",
"s = 1\n",
"\n",
"ys = np.convolve(ys, np.ones((s,))/s, mode='valid')\n",
"\n",
"xs = np.linspace(0, 1.0/2.0, len(data)/2)\n",
"\n",
"ys *= 2*np.pi*xs\n",
"#xs = np.linspace(len(data)/2, 1, len(data)/2)\n",
"\n",
"fig, ax = plt.subplots()\n",
"ax.loglog(xs[s//2:-s//2+1 if s > 1 else None], ys)\n",
"ax.xaxis.set_major_formatter(plt.FuncFormatter(lambda x, _pos: 1/x))\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"ys = scipy.fftpack.fft(data[:,0])\n",
"ys = 2.0/len(data) * np.abs(ys[:len(data)//2])\n",
"s = 30\n",
"\n",
"ys = np.convolve(ys, np.ones((s,))/s, mode='valid')\n",
"\n",
"xs = np.linspace(0, 1.0/2.0, len(data)/2)\n",
"\n",
"ys *= 2*np.pi*xs[s//2:-s//2+1]\n",
"\n",
"#xs = np.linspace(len(data)/2, 1, len(data)/2)\n",
"\n",
"fig, ax = plt.subplots(figsize=(9,5))\n",
"ax.loglog(xs[s//2:-s//2+1], ys)\n",
"ax.xaxis.set_major_formatter(plt.FuncFormatter(lambda x, _pos: 1/x))\n",
"ax.grid()\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"1/0.0628"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "labenv",
"language": "python",
"name": "labenv"
},
"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.8.1"
}
},
"nbformat": 4,
"nbformat_minor": 4
}