diff --git a/examples_py/test_zpkmodel.ipynb b/examples_py/test_zpkmodel.ipynb new file mode 100644 index 0000000..c5b5668 --- /dev/null +++ b/examples_py/test_zpkmodel.ipynb @@ -0,0 +1,148 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "bcd07c1e-6722-44c7-b162-267f1341c2fe", + "metadata": {}, + "source": [ + "# Compare Butterworth bandpass implementation with Scipy's implementation" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "f0c4d6ef-69b2-4b5c-92fd-d987a11f3cbd", + "metadata": {}, + "outputs": [], + "source": [ + "# Prerequisites, uncomment below in case of errors. Also for ipympl, restart Jupyter Lab if it was not installed\n", + "\n", + "#!cd .. && maturin develop -F python-bindings\n", + "#!pip install ipympl scipy matplotlib" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "23b45cf6-85bc-4eaa-8f27-77c5b354e772", + "metadata": {}, + "outputs": [], + "source": [ + "# If this does not work, install ipympl and reboot Jupyter Lab\n", + "%matplotlib widget" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "96090ac9-2033-411d-8e45-0c6b375a268b", + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "from lasprs.filter import Biquad, SeriesBiquad, BiquadBank, FilterSpec, ZPKModel" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8c137d60-51f3-4a14-892c-925aad662dc5", + "metadata": {}, + "outputs": [], + "source": [ + "#!pip install matplotlib\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np\n", + "from numpy import log10, sqrt, exp, pi\n", + "from scipy.signal import butter, freqs_zpk\n", + "# Cut-\n", + "fc = 100\n", + "fl = fc/1.1\n", + "fu = fc*1.5\n", + "order = 4\n", + "fs = 32e3" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "25d36d28-80e9-40cc-9364-f708d21f9b49", + "metadata": {}, + "outputs": [], + "source": [ + "# Create Scipy filter\n", + "but = butter(order, [2*pi*fl, 2*pi*fu],'bandpass', analog=True, output ='zpk')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "1f877979-55a8-4eee-b1f3-fc8500414664", + "metadata": {}, + "outputs": [], + "source": [ + "\n", + "# fspec = FilterSpec.Lowpass(fc, order)\n", + "# fspec = FilterSpec.Highpass(fc, order)\n", + "fspec = FilterSpec.Bandpass(fl, fu, order)\n", + "m = ZPKModel.butter(fspec)\n", + "freq = np.logspace(log10(fc)-1, log10(fc)+1, 100)\n", + "def level(a):\n", + " return 20*np.log10(np.abs(a))\n", + "#m" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "3ae5384d-2b5e-4688-b20d-148120e66147", + "metadata": {}, + "outputs": [], + "source": [ + "w, hbut = freqs_zpk(*but, worN=2*np.pi*freq)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "ad12c873-a47f-4ed5-a5ad-bb4d45bc5c3a", + "metadata": {}, + "outputs": [], + "source": [ + "plt.figure()\n", + "plt.subplot(211)\n", + "plt.title('Magnitude [dB]')\n", + "plt.semilogx(freq, level(hbut), lw=3)\n", + "plt.semilogx(freq, level(m.tf(1e-16, freq)), '--')\n", + "plt.ylim(-10, 1)\n", + "plt.legend(['Scipy', 'Lasprs'])\n", + "plt.subplot(212)\n", + "plt.title('Phase [deg]')\n", + "plt.semilogx(freq, np.angle(m.tf(1e-16, freq), deg=True))\n", + "plt.semilogx(freq, np.angle(hbut, deg=True))\n", + "plt.xlabel('Freq. [Hz]')" + ] + } + ], + "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.12.4" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}