{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "9f01df51",
   "metadata": {},
   "source": [
    "# Digitális szűrők — gyakorlatok\n",
    "\n",
    "**Államvizsga-tétel kísérő gyakorlófüzete**\n",
    "\n",
    "Forrás: Kuczmann Miklós, *Jelek és rendszerek* (Széchenyi István Egyetem, 2007), 7–10. fejezet, kiegészítve a klasszikus DSP-szakirodalommal (Oppenheim–Schafer, Proakis–Manolakis).\n",
    "\n",
    "---\n",
    "\n",
    "## A füzet felépítése\n",
    "\n",
    "12 gyakorlat, mindegyik egy önálló cellacsoport:\n",
    "\n",
    "1. Impulzusválasz, konvolúció és differenciaegyenlet  \n",
    "2. Frekvenciakarakterisztika ($W(e^{j\\vartheta})$) ábrázolása  \n",
    "3. Pólus-zérus elrendezés és stabilitás  \n",
    "4. Ideális szűrő impulzusválasza (sinc-függvény)  \n",
    "5. FIR aluláteresztő ablakozással — ablakok összehasonlítása  \n",
    "6. Lineáris fázis és csoportkésleltetés  \n",
    "7. Parks–McClellan (Remez) optimum FIR-tervezés  \n",
    "8. Bilineáris transzformáció és frekvencia-elhajlás  \n",
    "9. IIR prototípusok összehasonlítása (Butterworth / Chebyshev / Elliptic)  \n",
    "10. Realizációs struktúrák: direct form vs biquad cascade  \n",
    "11. Kvantálási hatások és limit cycle  \n",
    "12. Gyakorlati alkalmazás — zajos jel szűrése\n",
    "\n",
    "## Hivatkozások a tétel képleteire\n",
    "\n",
    "- $(7.19)$ rendszeregyenlet, $(7.5)$ konvolúció — időtartomány  \n",
    "- $(8.21)$ átviteli karakterisztika — DTFT-tartomány  \n",
    "- $(9.7)$ átviteli függvény, $(9.42)$ pólus-zérus alak — $z$-tartomány  \n",
    "- $(10.17)$ ideális aluláteresztő, $(10.25)$ bilineáris transzformáció\n",
    "\n",
    "---\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "79c6af3d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Közös import-blokk - minden gyakorlat ezeket használja\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import signal\n",
    "\n",
    "# Egységes ábra-stílus\n",
    "plt.rcParams.update({\n",
    "    'figure.figsize':   (10, 4),\n",
    "    'figure.dpi':       100,\n",
    "    'axes.grid':        True,\n",
    "    'grid.alpha':       0.3,\n",
    "    'axes.spines.top':  False,\n",
    "    'axes.spines.right': False,\n",
    "    'font.size':        10,\n",
    "    'lines.linewidth':  1.6,\n",
    "})\n",
    "\n",
    "# Egységes színpaletta a tétellel összhangban\n",
    "NAVY  = '#0F2A47'\n",
    "CYAN  = '#0891B2'\n",
    "AMBER = '#F59E0B'\n",
    "RED   = '#DC2626'\n",
    "MUTE  = '#64748B'\n",
    "\n",
    "print('Környezet betöltve.')\n",
    "print(f'  numpy {np.__version__}, scipy {__import__(\"scipy\").__version__}')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "49f78e5b",
   "metadata": {},
   "source": [
    "## 1. Impulzusválasz, konvolúció és differenciaegyenlet\n",
    "\n",
    "A diszkrét idejű, lineáris, invariáns rendszer válaszát három ekvivalens módon számíthatjuk:\n",
    "\n",
    "1. **Differenciaegyenlettel** $(7.19)$:  \n",
    "   $\\quad y[k] + \\sum_{i=1}^{n} a_i\\,y[k{-}i] = \\sum_{i=0}^{m} b_i\\,s[k{-}i]$\n",
    "2. **Konvolúcióval** $(7.5)$:  \n",
    "   $\\quad y[k] = s[k] * w[k]$\n",
    "3. **$z$-tartományban**:  \n",
    "   $\\quad Y(z) = W(z)\\,S(z)$,  majd inverz $z$-transzformáció.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Adott egy egyszerű FIR és egy egyszerű IIR rendszer:\n",
    "\n",
    "- **FIR**: 5 csapos mozgóátlag,  $w_{\\mathrm{FIR}}[k] = \\frac{1}{5}$,  $k=0,1,2,3,4$\n",
    "- **IIR**: $W_{\\mathrm{IIR}}(z) = \\dfrac{1}{1 - 0{,}9\\,z^{-1}}$,  azaz $y[k] = s[k] + 0{,}9\\,y[k{-}1]$\n",
    "\n",
    "Számítsd ki mindkét rendszer válaszát egy $s[k] = \\varepsilon[k]$ egységugrásra **két különböző módszerrel**: (a) `np.convolve` segítségével az impulzusválasszal, (b) `signal.lfilter` segítségével a differenciaegyenlettel. Ellenőrizd, hogy az eredmények megegyeznek.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e9f3a548",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Diszkrét idő tengely\n",
    "k = np.arange(0, 40)\n",
    "\n",
    "# Egységugrás gerjesztés\n",
    "s = np.ones_like(k, dtype=float)\n",
    "\n",
    "# ----- FIR: 5-csapos mozgóátlag -----\n",
    "b_fir = np.ones(5) / 5.0\n",
    "a_fir = np.array([1.0])\n",
    "\n",
    "# Impulzusválasz (méréssel: egységimpulzus -> w[k])\n",
    "w_fir = signal.lfilter(b_fir, a_fir, np.r_[1.0, np.zeros(len(k)-1)])\n",
    "\n",
    "# Válasz konvolúcióval és lfilter-rel\n",
    "y_fir_conv  = np.convolve(s, w_fir)[:len(k)]\n",
    "y_fir_lfilt = signal.lfilter(b_fir, a_fir, s)\n",
    "\n",
    "# ----- IIR: y[k] = s[k] + 0.9 y[k-1] -----\n",
    "b_iir = np.array([1.0])\n",
    "a_iir = np.array([1.0, -0.9])\n",
    "\n",
    "w_iir = signal.lfilter(b_iir, a_iir, np.r_[1.0, np.zeros(len(k)-1)])\n",
    "y_iir_conv  = np.convolve(s, w_iir)[:len(k)]\n",
    "y_iir_lfilt = signal.lfilter(b_iir, a_iir, s)\n",
    "\n",
    "# ----- Ellenőrzés: numerikus pontossággal egyezik? -----\n",
    "print(f'FIR max. eltérés (konvolúció vs lfilter): {np.max(np.abs(y_fir_conv - y_fir_lfilt)):.2e}')\n",
    "print(f'IIR max. eltérés (konvolúció vs lfilter): {np.max(np.abs(y_iir_conv - y_iir_lfilt)):.2e}')\n",
    "\n",
    "# ----- Vizualizáció -----\n",
    "fig, axes = plt.subplots(2, 2, figsize=(11, 6.5))\n",
    "\n",
    "# FIR\n",
    "axes[0,0].stem(k, w_fir, basefmt=' ', linefmt=CYAN, markerfmt='o')\n",
    "axes[0,0].set_title('FIR impulzusválasz $w[k]$ (5-csapos mozgóátlag)')\n",
    "axes[0,0].set_xlabel('k'); axes[0,0].set_ylabel('w[k]')\n",
    "\n",
    "axes[0,1].step(k, y_fir_lfilt, where='post', color=NAVY, label='differenciaegyenl.')\n",
    "axes[0,1].plot(k, y_fir_conv, 'o', color=AMBER, markersize=4, label='konvolúció')\n",
    "axes[0,1].set_title('FIR válasz egységugrásra (ugrásválasz)')\n",
    "axes[0,1].set_xlabel('k'); axes[0,1].set_ylabel('y[k]'); axes[0,1].legend()\n",
    "\n",
    "# IIR\n",
    "axes[1,0].stem(k, w_iir, basefmt=' ', linefmt=AMBER, markerfmt='o')\n",
    "axes[1,0].set_title('IIR impulzusválasz $w[k] = 0.9^k$')\n",
    "axes[1,0].set_xlabel('k'); axes[1,0].set_ylabel('w[k]')\n",
    "\n",
    "axes[1,1].step(k, y_iir_lfilt, where='post', color=NAVY, label='differenciaegyenl.')\n",
    "axes[1,1].plot(k, y_iir_conv, 'o', color=AMBER, markersize=4, label='konvolúció')\n",
    "axes[1,1].set_title('IIR válasz egységugrásra')\n",
    "axes[1,1].set_xlabel('k'); axes[1,1].set_ylabel('y[k]'); axes[1,1].legend()\n",
    "\n",
    "plt.tight_layout(); plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c1b298b2",
   "metadata": {},
   "source": [
    "**Megfigyelés.** Mindkét módszer numerikusan azonos eredményt ad (a maradék hiba a `float64` aritmetika kvantálási zaja). A FIR rendszer véges $K{+}1{=}5$ ütem alatt eléri az állandósult állapotot ($y[k]{=}1$ az $k\\ge 4$ ütemekben), míg az IIR rendszer **mértani sorral** közelít $1/(1-0{,}9)=10$-hez. Az IIR impulzusválasza végtelen, de **abszolút összegezhető**, mert $|p_1|=0{,}9<1$ — tehát a rendszer BIBO-stabilis.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "17a369a5",
   "metadata": {},
   "source": [
    "## 2. Frekvenciakarakterisztika — $W(e^{j\\vartheta})$\n",
    "\n",
    "Az átviteli karakterisztika $(8.21)$ képlete a normált $\\vartheta = \\omega T_s$ körfrekvencián:\n",
    "\n",
    "$$\n",
    "W(e^{j\\vartheta}) = \\frac{\\sum_{i=0}^{m} b_i\\,e^{-ji\\vartheta}}{1 + \\sum_{i=1}^{n} a_i\\,e^{-ji\\vartheta}}\n",
    "$$\n",
    "\n",
    "A $\\vartheta$ a $[0,\\pi]$ intervallumon mozog, ahol $\\vartheta=\\pi$ a Nyquist-frekvenciának felel meg. A `scipy.signal.freqz` függvény pont ezt az értéket számolja a $z=e^{j\\vartheta}$ helyettesítéssel.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Vizsgáljuk a következő rendszert:\n",
    "\n",
    "$$\n",
    "W(z) \\;=\\; \\frac{1 + 2z^{-1} + z^{-2}}{1 - 0{,}5\\,z^{-1} + 0{,}25\\,z^{-2}}\n",
    "$$\n",
    "\n",
    "Ábrázold az amplitúdókarakterisztikát lineáris és dB-skálán, valamint a fáziskarakterisztikát. Olvasd le, hogy ez milyen szűrőtípus (LP / HP / BP / BS).\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "26fa37f6",
   "metadata": {},
   "outputs": [],
   "source": [
    "b = np.array([1.0, 2.0, 1.0])\n",
    "a = np.array([1.0, -0.5, 0.25])\n",
    "\n",
    "# Frekvenciaválasz a [0, π] tartományon, 1024 ponton\n",
    "w_norm, H = signal.freqz(b, a, worN=1024)\n",
    "\n",
    "mag    = np.abs(H)\n",
    "mag_dB = 20.0 * np.log10(np.maximum(mag, 1e-12))\n",
    "phase  = np.unwrap(np.angle(H))\n",
    "\n",
    "fig, axes = plt.subplots(1, 3, figsize=(13, 3.6))\n",
    "\n",
    "axes[0].plot(w_norm/np.pi, mag, color=CYAN)\n",
    "axes[0].set_title(r'Amplitúdó  $|W(e^{j\\vartheta})|$')\n",
    "axes[0].set_xlabel(r'$\\vartheta / \\pi$'); axes[0].set_ylabel('|W|')\n",
    "\n",
    "axes[1].plot(w_norm/np.pi, mag_dB, color=NAVY)\n",
    "axes[1].set_title('Amplitúdó (dB)')\n",
    "axes[1].set_xlabel(r'$\\vartheta / \\pi$'); axes[1].set_ylabel('|W| [dB]')\n",
    "axes[1].set_ylim(-40, 25); axes[1].axhline(-3, color=AMBER, linestyle='--', alpha=0.6, label='-3 dB')\n",
    "axes[1].legend()\n",
    "\n",
    "axes[2].plot(w_norm/np.pi, phase, color=AMBER)\n",
    "axes[2].set_title(r'Fázis  $\\arg W(e^{j\\vartheta})$')\n",
    "axes[2].set_xlabel(r'$\\vartheta / \\pi$'); axes[2].set_ylabel('fázis [rad]')\n",
    "\n",
    "plt.tight_layout(); plt.show()\n",
    "\n",
    "# Vágási frekvencia leolvasása (-3 dB)\n",
    "ix_3dB = np.argmin(np.abs(mag_dB - (-3)))\n",
    "print(f'-3 dB pont:   ϑ_c ≈ {w_norm[ix_3dB]/np.pi:.3f} π')\n",
    "print(f'DC erősítés:  |W(1)|   = {np.abs(np.polyval(b, 1)/np.polyval(a, 1)):.3f}')\n",
    "print(f'Nyq. erősítés: |W(-1)| = {np.abs(np.polyval(b, -1)/np.polyval(a, -1)):.3f}')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e02e2809",
   "metadata": {},
   "source": [
    "**Megfigyelés.** A DC ($\\vartheta{=}0$) erősítés $|W(1)| = 4/0{,}75 \\approx 5{,}33$, a Nyquist-frekvencián ($\\vartheta{=}\\pi$) $|W(-1)| = 0$ — ez egy **aluláteresztő** karakterisztika. A számláló $z = -1$-ben (azaz $\\vartheta = \\pi$-ben) kettős zérust ad, ami a frekvenciaválasz teljes elnyomásához vezet.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "26d3e5b2",
   "metadata": {},
   "source": [
    "## 3. Pólus-zérus elrendezés és stabilitás\n",
    "\n",
    "Az átviteli függvény gyöktényezős alakja $(9.42)$:\n",
    "\n",
    "$$\n",
    "W(z) = K \\cdot \\frac{(z-z_1)(z-z_2)\\cdots(z-z_m)}{(z-p_1)(z-p_2)\\cdots(z-p_n)}\n",
    "$$\n",
    "\n",
    "**Stabilitási tétel.** Egy LTI, kauzális diszkrét idejű rendszer akkor és csakis akkor BIBO-stabilis, ha $|p_i| < 1$ minden $i$-re — minden pólus az egységkör belsejében.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Vizsgáljuk négy különböző rendszer pólus-zérus elrendezését. Döntsd el mindegyikről, hogy stabilis-e, és figyeld meg a frekvenciaválasz alakja és a pólusok elhelyezkedése közötti összefüggést.\n",
    "\n",
    "| Rendszer | $b$ együtthatók | $a$ együtthatók |\n",
    "|---|---|---|\n",
    "| (a) | $[1, -1{,}414, 1]$ | $[1, -1{,}273, 0{,}81]$ |\n",
    "| (b) | $[1, -1{,}414, 1]$ | $[1, -1{,}414, 1]$ |\n",
    "| (c) | $[1, 0, 0]$        | $[1, -1{,}5, 0{,}9]$ |\n",
    "| (d) | $[1, 0, 0]$        | $[1, -1{,}5, 1{,}1]$ |\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6130757e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_pzmap_and_freq(b, a, ax_pz, ax_mag, title=''):\n",
    "    \"\"\"Pólus-zérus rajz egységkörrel + frekvenciaválasz.\"\"\"\n",
    "    z, p, k = signal.tf2zpk(b, a)\n",
    "\n",
    "    # --- Pólus-zérus diagram ---\n",
    "    theta = np.linspace(0, 2*np.pi, 200)\n",
    "    ax_pz.plot(np.cos(theta), np.sin(theta), '--', color=CYAN, lw=1)\n",
    "    ax_pz.axhline(0, color=MUTE, lw=0.5); ax_pz.axvline(0, color=MUTE, lw=0.5)\n",
    "\n",
    "    if len(z):\n",
    "        ax_pz.plot(z.real, z.imag, 'o', color=AMBER, ms=10, mfc='none', mew=2, label='zérus')\n",
    "    if len(p):\n",
    "        unstable = np.abs(p) >= 1\n",
    "        ax_pz.plot(p[~unstable].real, p[~unstable].imag, 'x', color=NAVY, ms=12, mew=2, label='pólus')\n",
    "        if unstable.any():\n",
    "            ax_pz.plot(p[unstable].real, p[unstable].imag, 'x', color=RED, ms=12, mew=2, label='instabil pólus')\n",
    "\n",
    "    ax_pz.set_xlim(-2, 2); ax_pz.set_ylim(-2, 2)\n",
    "    ax_pz.set_aspect('equal'); ax_pz.set_title(title); ax_pz.legend(loc='upper right', fontsize=8)\n",
    "    ax_pz.set_xlabel('Re(z)'); ax_pz.set_ylabel('Im(z)')\n",
    "\n",
    "    # --- Frekvenciaválasz ---\n",
    "    is_stable = np.all(np.abs(p) < 1) if len(p) else True\n",
    "    if is_stable:\n",
    "        wn, H = signal.freqz(b, a, worN=1024)\n",
    "        ax_mag.plot(wn/np.pi, 20*np.log10(np.maximum(np.abs(H), 1e-12)), color=NAVY)\n",
    "        ax_mag.set_ylim(-40, 30)\n",
    "    else:\n",
    "        ax_mag.text(0.5, 0.5, 'INSTABIL\\n(nincs frekvenciaválasz)', ha='center', va='center',\n",
    "                    fontsize=12, color=RED, transform=ax_mag.transAxes, fontweight='bold')\n",
    "    ax_mag.set_xlabel(r'$\\vartheta/\\pi$'); ax_mag.set_ylabel('|W| [dB]')\n",
    "    ax_mag.set_title('Amplitúdó (dB)')\n",
    "\n",
    "    return is_stable, np.abs(p) if len(p) else np.array([])\n",
    "\n",
    "systems = [\n",
    "    ('(a) sávzáró',     [1, -1.414, 1],   [1, -1.273, 0.81]),\n",
    "    ('(b) sávzáró-szélen', [1, -1.414, 1], [1, -1.414, 1]),\n",
    "    ('(c) rezonáns LP', [1, 0, 0],         [1, -1.5, 0.9]),\n",
    "    ('(d) instabil',    [1, 0, 0],         [1, -1.5, 1.1]),\n",
    "]\n",
    "\n",
    "fig, axes = plt.subplots(4, 2, figsize=(11, 14))\n",
    "for i, (name, b, a) in enumerate(systems):\n",
    "    stable, abs_p = plot_pzmap_and_freq(b, a, axes[i,0], axes[i,1], title=name)\n",
    "    print(f'{name:30s}  |p| = {abs_p.round(3)}   stabilis: {stable}')\n",
    "\n",
    "plt.tight_layout(); plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "532ac1eb",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- **(a)** A pólusok és zérusok közel vannak az egységkörhöz, de a pólusok belül vannak ($|p|=0{,}9$) — stabilis sávzáró karakterisztika.\n",
    "- **(b)** A pólusok és zérusok pontosan az egységkörön vannak ($|p|=1$) — **határeset**, gyakorlatban a numerikus zaj könnyen instabillá teheti.\n",
    "- **(c)** Egy konjugált pólus-pár az egységkörhöz közel ($|p|\\approx 0{,}949$) éles **rezonancia-csúcsot** okoz az amplitúdókarakterisztikán.\n",
    "- **(d)** A pólusok az egységkörön **kívül** vannak — a rendszer instabil, az impulzusválasz exponenciálisan divergál, frekvenciaválasz nem értelmezhető.\n",
    "\n",
    "**Geometriai szabály**: minél közelebb van egy pólus az egységkörhöz, annál élesebb a megfelelő frekvenciánál a csúcs. Egy egységkörhöz közeli zérus pedig mély bevágást (notch) okoz.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "acbd8746",
   "metadata": {},
   "source": [
    "## 4. Ideális aluláteresztő szűrő impulzusválasza\n",
    "\n",
    "A tankönyv $(10.17)$ képlete az ideális aluláteresztő átviteli karakterisztikájáról:\n",
    "\n",
    "$$\n",
    "W_\\Omega(j\\omega) = \\begin{cases} T_s/\\tau, & |\\omega|\\le \\omega_s/2 \\\\ 0, & |\\omega|>\\omega_s/2 \\end{cases}\n",
    "$$\n",
    "\n",
    "A megfelelő impulzusválasz egy **sinc-függvény**:\n",
    "\n",
    "$$\n",
    "w_d[k] \\;=\\; \\frac{\\sin(\\vartheta_c\\, k)}{\\pi\\,k}\n",
    "$$\n",
    "\n",
    "Ez a függvény nem belépő ($w_d[k]\\ne 0$ a $k<0$ ütemekben is) — tehát a rendszer **nem kauzális, nem realizálható**. Ezért „ideális\".\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Számítsd ki az ideális aluláteresztő szűrő impulzusválaszát $\\vartheta_c = 0{,}3\\pi$ vágási frekvenciával **két módon**:\n",
    "\n",
    "1. közvetlenül a $w_d[k] = \\sin(\\vartheta_c k)/(\\pi k)$ képlettel,\n",
    "2. inverz DTFT-vel egy ablakos $W_d(e^{j\\vartheta})$-ből.\n",
    "\n",
    "Ezután csonkold le a végtelen impulzusválaszt $|k|\\le M$ tartományra, és figyeld meg, hogy a frekvenciaválaszon megjelenik a **Gibbs-jelenség** (oszcilláció a vágási frekvencia körül).\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "950d2c59",
   "metadata": {},
   "outputs": [],
   "source": [
    "theta_c = 0.3 * np.pi\n",
    "\n",
    "# (1) Közvetlen sinc-formula\n",
    "M = 50\n",
    "k = np.arange(-M, M+1)\n",
    "with np.errstate(divide='ignore', invalid='ignore'):\n",
    "    w_ideal = np.where(k == 0, theta_c/np.pi, np.sin(theta_c*k)/(np.pi*k))\n",
    "\n",
    "# (2) Inverz DTFT mintavételezett W_d-ből\n",
    "N_dft = 1024\n",
    "theta = np.linspace(-np.pi, np.pi, N_dft, endpoint=False)\n",
    "Wd = np.where(np.abs(theta) <= theta_c, 1.0, 0.0)\n",
    "# inverz DFT (a frekvenciatengelyt nullára shifteljük)\n",
    "w_idft = np.fft.fftshift(np.fft.ifft(np.fft.ifftshift(Wd))).real\n",
    "k_idft = np.arange(-N_dft//2, N_dft//2)\n",
    "\n",
    "# Több csonkolási hossz\n",
    "fig, axes = plt.subplots(2, 2, figsize=(12, 7))\n",
    "\n",
    "# Bal felső: a ki nem csonkolt sinc\n",
    "axes[0,0].stem(k, w_ideal, basefmt=' ', linefmt=CYAN, markerfmt='o')\n",
    "axes[0,0].set_title(r'Ideális LP impulzusválasz  $w_d[k] = \\sin(\\vartheta_c k) / (\\pi k)$')\n",
    "axes[0,0].set_xlabel('k'); axes[0,0].set_ylabel('$w_d$[k]')\n",
    "axes[0,0].axvline(0, color=AMBER, ls='--', alpha=0.5)\n",
    "axes[0,0].text(2, 0.27, 'k<0: nem-kauzális rész', color=AMBER, fontsize=9)\n",
    "\n",
    "# Jobb felső: a két számítás egyezőségének ellenőrzése\n",
    "ix = (k_idft >= -M) & (k_idft <= M)\n",
    "axes[0,1].plot(k, w_ideal, 'o', color=CYAN, ms=4, label='formula')\n",
    "axes[0,1].plot(k_idft[ix], w_idft[ix], '-', color=NAVY, lw=1, alpha=0.7, label='IDFT')\n",
    "axes[0,1].set_title('A két módszer egyezősége')\n",
    "axes[0,1].set_xlabel('k'); axes[0,1].set_ylabel('$w_d$[k]'); axes[0,1].legend()\n",
    "\n",
    "# Csonkolás vizsgálata: különböző hosszak frekvenciaválasza\n",
    "for ax_idx, M_short in enumerate([5, 25]):\n",
    "    k_t = np.arange(-M_short, M_short+1)\n",
    "    with np.errstate(divide='ignore', invalid='ignore'):\n",
    "        w_t = np.where(k_t == 0, theta_c/np.pi, np.sin(theta_c*k_t)/(np.pi*k_t))\n",
    "    # FFT a csonkolt válaszból\n",
    "    N_fft = 4096\n",
    "    W_t = np.abs(np.fft.fft(w_t, N_fft))\n",
    "    f_t = np.linspace(0, 2, N_fft, endpoint=False)\n",
    "    axes[1, ax_idx].plot(f_t[:N_fft//2], W_t[:N_fft//2], color=CYAN if ax_idx==0 else NAVY)\n",
    "    axes[1, ax_idx].axvline(0.3, color=AMBER, ls='--', alpha=0.6, label=r'$\\vartheta_c$')\n",
    "    axes[1, ax_idx].set_title(f'Csonkolt válasz, M={M_short}  (Gibbs-jelenség)')\n",
    "    axes[1, ax_idx].set_xlabel(r'$\\vartheta/\\pi$')\n",
    "    axes[1, ax_idx].set_ylabel('|W|'); axes[1, ax_idx].legend()\n",
    "    axes[1, ax_idx].set_xlim(0, 1)\n",
    "\n",
    "plt.tight_layout(); plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "685fcfe9",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A $w_d[k]$ értékei a $k<0$ ütemekben sem nullák — a rendszer **nem kauzális**, ezért a fizikai időben nem realizálható.\n",
    "- A direkt sinc-formula és az inverz DTFT azonos eredményt ad (kvadratikus pontosságon belül).\n",
    "- **Csonkolás után**: a frekvenciaválaszon **átmeneti tartomány** jelenik meg a vágási frekvencia körül, és **Gibbs-oszcillációk** keletkeznek (a `M=5` esetben jól láthatóak). Hosszabb $M$ esetén az átmenet élesebb, de az oszcillációk amplitúdója (kb. 9%) nem csökken — ez a Gibbs-jelenség lényege.\n",
    "- A megoldás: az éles téglalap ablak helyett egy simán kifutó ablakot (Hann, Hamming, Blackman, Kaiser) használunk — ezt vizsgáljuk az 5. gyakorlatban.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4ec4ef14",
   "metadata": {},
   "source": [
    "## 5. FIR aluláteresztő tervezése ablakozással\n",
    "\n",
    "Az ablakos FIR-tervezés a tétel 5. fejezetének központi módszere:\n",
    "\n",
    "$$\n",
    "h[k] = w_d[k] \\cdot w_w[k], \\qquad k=0,\\dots,N-1\n",
    "$$\n",
    "\n",
    "ahol $w_w[k]$ az ablakfüggvény. Az ablak megválasztása szabja meg a kompromisszumot az **átmeneti tartomány szélessége** és a **záró-tartománybeli csillapítás** között.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Tervezz $N=51$-es aluláteresztő FIR szűrőt $\\vartheta_c = 0{,}3\\pi$ normált vágási frekvenciával, négy különböző ablakkal: **rectangular**, **Hann**, **Hamming**, **Blackman**. Vizsgáld:\n",
    "\n",
    "- az amplitúdóválaszt dB-ben,\n",
    "- az áteresztő-tartománybeli ingadozást,\n",
    "- a záró-tartománybeli csillapítást,\n",
    "- a maximális mellékhullám szintjét.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9377da44",
   "metadata": {},
   "outputs": [],
   "source": [
    "N      = 51\n",
    "fc_pi  = 0.3   # ϑ_c / π\n",
    "fs_des = 2.0   # 'normált' mintavételi frekvencia (Nyquist = 1)\n",
    "\n",
    "windows = [\n",
    "    ('Rectangular', 'boxcar',  CYAN),\n",
    "    ('Hann',        'hann',    AMBER),\n",
    "    ('Hamming',     'hamming', NAVY),\n",
    "    ('Blackman',    'blackman', RED),\n",
    "]\n",
    "\n",
    "fig, (ax_h, ax_W) = plt.subplots(1, 2, figsize=(13, 4.5))\n",
    "\n",
    "print(f'{\"Ablak\":<12} {\"Áter. ripple\":<15} {\"Záró csill.\":<15} {\"Mellékhullám\":<15}')\n",
    "print('-' * 60)\n",
    "\n",
    "for name, win, col in windows:\n",
    "    # signal.firwin a sinc * ablak FIR-tervezést végzi\n",
    "    h = signal.firwin(N, fc_pi, window=win)\n",
    "\n",
    "    # Frekvenciaválasz\n",
    "    w_n, H = signal.freqz(h, 1, worN=4096)\n",
    "    H_dB = 20*np.log10(np.maximum(np.abs(H), 1e-12))\n",
    "\n",
    "    ax_h.plot(np.arange(N), h, '-', color=col, label=name, alpha=0.8)\n",
    "    ax_W.plot(w_n/np.pi, H_dB, color=col, label=name)\n",
    "\n",
    "    # Számszerű jellemzők\n",
    "    passband = w_n/np.pi <= fc_pi - 0.05\n",
    "    stopband = w_n/np.pi >= fc_pi + 0.10\n",
    "    pass_ripple = np.max(H_dB[passband]) - np.min(H_dB[passband])\n",
    "    stop_atten  = -np.max(H_dB[stopband])\n",
    "    sidelobe    = -np.max(H_dB[stopband])\n",
    "    print(f'{name:<12} {pass_ripple:>8.3f} dB     {stop_atten:>8.1f} dB     {sidelobe:>8.1f} dB')\n",
    "\n",
    "ax_h.set_title(f'h[k] - {N}-csapos LP FIR (ablakozás)')\n",
    "ax_h.set_xlabel('k'); ax_h.set_ylabel('h[k]'); ax_h.legend()\n",
    "\n",
    "ax_W.axvline(fc_pi, color=MUTE, ls='--', alpha=0.5)\n",
    "ax_W.set_title('Amplitúdó (dB)')\n",
    "ax_W.set_xlabel(r'$\\vartheta/\\pi$'); ax_W.set_ylabel('|W| [dB]')\n",
    "ax_W.set_ylim(-110, 5); ax_W.set_xlim(0, 1); ax_W.legend()\n",
    "\n",
    "plt.tight_layout(); plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e187fe66",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A **rectangular** ablak ad a legmeredekebb átmenetet, **de** csak ~21 dB záró-csillapítást — ez a Gibbs-jelenség közvetlen következménye.\n",
    "- A **Hann** és **Hamming** ablakok már jó kompromisszumok: kb. 44 és 53 dB záró-csillapítás, mérsékelt átmeneti szélesség.\n",
    "- A **Blackman** ablak ad ~74 dB csillapítást, de szélesebb átmeneti tartományt — adott $N$-nél a fő nyaláb is szélesebb.\n",
    "- A **Kaiser**-ablak ($\\beta$ paraméterrel) folytonosan hangolható ezek között a kompromisszumok között — lásd `signal.kaiser_atten` és `signal.kaiser_beta`.\n",
    "- Magasabb csillapítás eléréséhez vagy keskenyebb ablakot kell választani, **vagy** növelni kell az $N$ szűrőhosszot.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "035c479e",
   "metadata": {},
   "source": [
    "## 6. Lineáris fázis és csoportkésleltetés\n",
    "\n",
    "A FIR fő erőssége a könnyen elérhető **lineáris fázis**: ha $h[i] = \\pm h[N-1-i]$ (szimmetrikus / antiszimmetrikus impulzusválasz), akkor\n",
    "\n",
    "$$\n",
    "\\arg W(e^{j\\vartheta}) = -\\tfrac{N-1}{2}\\,\\vartheta + \\varphi_0\n",
    "$$\n",
    "\n",
    "azaz a **csoportkésleltetés** $\\tau_g = -\\partial \\varphi/\\partial \\vartheta = (N-1)/2$ **konstans**, frekvenciafüggetlen. Ez **alakhű jelátvitelt** biztosít az áteresztő tartományban.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "1. Tervezz egy szimmetrikus FIR LP szűrőt és egy IIR LP szűrőt azonos vágási frekvenciával.\n",
    "2. Hasonlítsd össze a csoportkésleltetést.\n",
    "3. Vezess át egy darabos ($\\sum$ szinusz) tesztjelet mindkét szűrőn, és figyeld meg, hogy a FIR csak késlelteti, az IIR pedig **torzítja** (a komponensek eltolódnak egymáshoz képest).\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4de380e5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Tervezés azonos specifikációval\n",
    "N_fir = 81\n",
    "fc    = 0.25   # normált vágási (Nyquist=1)\n",
    "h_fir = signal.firwin(N_fir, fc, window='hamming')\n",
    "\n",
    "b_iir, a_iir = signal.butter(8, fc, btype='low')\n",
    "\n",
    "# --- Csoportkésleltetés ---\n",
    "w_n, gd_fir = signal.group_delay((h_fir, 1.0), w=1024)\n",
    "_,   gd_iir = signal.group_delay((b_iir, a_iir), w=1024)\n",
    "\n",
    "# --- Áteresztő tartománybeli teszt: két szinuszos komponens ---\n",
    "fs = 1.0\n",
    "n  = np.arange(0, 200)\n",
    "f1, f2 = 0.05, 0.15      # a vágási alatt, az áteresztő tartományban\n",
    "x = np.sin(2*np.pi*f1*n) + 0.5*np.sin(2*np.pi*f2*n + 0.5)\n",
    "\n",
    "y_fir = signal.lfilter(h_fir, 1.0, x)\n",
    "y_iir = signal.lfilter(b_iir, a_iir, x)\n",
    "\n",
    "# --- Vizualizáció ---\n",
    "fig, axes = plt.subplots(2, 2, figsize=(12, 6.5))\n",
    "\n",
    "axes[0,0].plot(w_n/np.pi, gd_fir, color=CYAN, label=f'FIR (N={N_fir})')\n",
    "axes[0,0].axhline((N_fir-1)/2, color=MUTE, ls='--', alpha=0.5,\n",
    "                  label=f'τ = (N-1)/2 = {(N_fir-1)/2:.0f}')\n",
    "axes[0,0].set_title('FIR csoportkésleltetés - konstans')\n",
    "axes[0,0].set_xlabel(r'$\\vartheta/\\pi$'); axes[0,0].set_ylabel(r'$\\tau_g$ [minta]')\n",
    "axes[0,0].set_ylim(38, 42)   # fix tengely - csak a ténylegesen érdekes tartomány\n",
    "axes[0,0].legend()\n",
    "\n",
    "axes[0,1].plot(w_n/np.pi, gd_iir, color=AMBER)\n",
    "axes[0,1].set_title('IIR csoportkésleltetés - frekvenciafüggő')\n",
    "axes[0,1].set_xlabel(r'$\\vartheta/\\pi$'); axes[0,1].set_ylabel(r'$\\tau_g$ [minta]')\n",
    "\n",
    "axes[1,0].plot(n, x,     color=MUTE, label='bemenet x[n]', alpha=0.7)\n",
    "axes[1,0].plot(n, y_fir, color=CYAN, label='FIR kimenet y[n]')\n",
    "axes[1,0].set_title('FIR: alakhű késleltetés')\n",
    "axes[1,0].set_xlabel('n'); axes[1,0].set_ylabel('érték'); axes[1,0].legend()\n",
    "axes[1,0].set_xlim(0, 200)\n",
    "\n",
    "axes[1,1].plot(n, x,     color=MUTE, label='bemenet x[n]', alpha=0.7)\n",
    "axes[1,1].plot(n, y_iir, color=AMBER, label='IIR kimenet y[n]')\n",
    "axes[1,1].set_title('IIR: fázistorzulás (komponensek eltolódnak)')\n",
    "axes[1,1].set_xlabel('n'); axes[1,1].set_ylabel('érték'); axes[1,1].legend()\n",
    "axes[1,1].set_xlim(0, 200)\n",
    "\n",
    "plt.tight_layout(); plt.show()\n",
    "\n",
    "print(f'FIR átlag csoportkésleltetés:  {gd_fir[w_n/np.pi < fc].mean():.2f} minta')\n",
    "print(f'IIR min/max csoportkésleltetés: {gd_iir[w_n/np.pi < fc].min():.2f} ... {gd_iir[w_n/np.pi < fc].max():.2f} minta')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71fa4a30",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A FIR csoportkésleltetése **pontosan konstans** $(N{-}1)/2 = 40$ minta — a kimenet a bemenet 40 mintával késleltetett, **alakhűen átvitt** változata.\n",
    "- Az IIR csoportkésleltetése **frekvenciafüggő**: a vágási frekvencia közelében jelentősen megnő, ami a hullámforma **torzulását** okozza (a két szinusz-komponens más-más mértékben tolódik el).\n",
    "- Ez a fő érv a FIR mellett audio, EKG, képfeldolgozás területén — minden olyan alkalmazásban, ahol a hullámforma megőrzése kritikus.\n",
    "- IIR esetén csak **allpass-szűrőkkel** lehet közelítőleg lineáris fázist elérni (Hilbert, fázis-egyenlítő szűrők).\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ba2a7e5e",
   "metadata": {},
   "source": [
    "## 7. Parks–McClellan optimum FIR-tervezés\n",
    "\n",
    "A Parks–McClellan algoritmus (más néven Remez-csere algoritmus) az ablakos módszerrel szemben **optimum** FIR szűrőt tervez a Csebisev-norma értelmében: a hibafüggvény **ekvi-ripple** lesz az áteresztő- és záró-tartományban — azaz adott specifikációhoz a **legkisebb fokszám** elérhető.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "1. Adott specifikáció: $\\vartheta_p = 0{,}25\\pi$, $\\vartheta_s = 0{,}30\\pi$, áteresztőbeli ripple $\\le 0{,}05$, záró-csillapítás $\\ge 50$ dB.\n",
    "2. Tervezz Parks–McClellan FIR-t (`signal.remez`).\n",
    "3. Hasonlítsd össze egy ugyanolyan fokszámú Hamming-ablakos FIR-rel.\n",
    "4. Figyeld meg az ekvi-ripple tulajdonságot.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d92f4f8c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Specifikáció\n",
    "fs_norm = 2.0\n",
    "theta_p, theta_s = 0.25, 0.30   # ϑ/π\n",
    "delta_p, delta_s = 0.05, 10**(-50/20)\n",
    "\n",
    "# Becsült fokszám (Kaiser-féle empirikus, ~ A_s / (22 * Δf))\n",
    "delta_omega = (theta_s - theta_p) * np.pi\n",
    "A_s = -20*np.log10(delta_s)\n",
    "N_est = int(np.ceil((A_s - 7.95) / (14.36 * delta_omega/(2*np.pi))))\n",
    "N = N_est | 1   # páratlanra állítjuk\n",
    "print(f'Becsült FIR fokszám: N = {N}  (átmeneti szélesség: {(theta_s-theta_p):.3f}π)')\n",
    "\n",
    "# Parks-McClellan szűrő\n",
    "bands  = [0, theta_p, theta_s, 1.0]   # ϑ/π egységekben\n",
    "desired = [1, 0]                       # áteresztő, záró\n",
    "weight  = [1/delta_p, 1/delta_s]       # súlyozás\n",
    "h_pm = signal.remez(N, bands, desired, weight=weight, fs=fs_norm)\n",
    "\n",
    "# Hamming ablakos szűrő ugyanezzel a hosszal\n",
    "h_win = signal.firwin(N, (theta_p+theta_s)/2, window='hamming')\n",
    "\n",
    "# Frekvenciaválaszok\n",
    "w_n, H_pm  = signal.freqz(h_pm,  1, worN=4096)\n",
    "_,   H_win = signal.freqz(h_win, 1, worN=4096)\n",
    "\n",
    "fig, axes = plt.subplots(1, 2, figsize=(13, 4.5))\n",
    "\n",
    "# Lineáris skála az áteresztő tartományban (ekvi-ripple szemléltetés)\n",
    "axes[0].plot(w_n/np.pi, np.abs(H_pm),  color=CYAN, label='Parks-McClellan')\n",
    "axes[0].plot(w_n/np.pi, np.abs(H_win), color=AMBER, label='Hamming-ablak', alpha=0.8)\n",
    "axes[0].axhline(1+delta_p, color=MUTE, ls='--', alpha=0.4)\n",
    "axes[0].axhline(1-delta_p, color=MUTE, ls='--', alpha=0.4)\n",
    "axes[0].axvline(theta_p, color=MUTE, ls=':', alpha=0.4)\n",
    "axes[0].set_title('Áteresztő tartomány közelítése')\n",
    "axes[0].set_xlim(0, 0.30); axes[0].set_ylim(0.94, 1.06)\n",
    "axes[0].set_xlabel(r'$\\vartheta/\\pi$'); axes[0].set_ylabel('|W|'); axes[0].legend()\n",
    "\n",
    "# dB skála a teljes tartományon\n",
    "axes[1].plot(w_n/np.pi, 20*np.log10(np.maximum(np.abs(H_pm), 1e-12)),  color=CYAN, label='Parks-McClellan')\n",
    "axes[1].plot(w_n/np.pi, 20*np.log10(np.maximum(np.abs(H_win), 1e-12)), color=AMBER, label='Hamming-ablak', alpha=0.8)\n",
    "axes[1].axhline(-50, color=MUTE, ls='--', alpha=0.4, label='-50 dB')\n",
    "axes[1].axvline(theta_s, color=MUTE, ls=':', alpha=0.4)\n",
    "axes[1].set_title('Teljes amplitúdó (dB)')\n",
    "axes[1].set_xlim(0, 1); axes[1].set_ylim(-90, 5)\n",
    "axes[1].set_xlabel(r'$\\vartheta/\\pi$'); axes[1].set_ylabel('|W| [dB]'); axes[1].legend()\n",
    "\n",
    "plt.tight_layout(); plt.show()\n",
    "\n",
    "# Ekvi-ripple ellenőrzés a záró tartományban\n",
    "sb_mask = w_n/np.pi >= theta_s\n",
    "print(f'PM    záróbeli max csillapítás: {-20*np.log10(np.max(np.abs(H_pm)[sb_mask])):.1f} dB')\n",
    "print(f'Ham.  záróbeli max csillapítás: {-20*np.log10(np.max(np.abs(H_win)[sb_mask])):.1f} dB')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "704ad952",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A Parks–McClellan szűrő mindkét tartományban **ekvi-ripple** — az ingadozások egyenlő amplitúdójúak, ez a Csebisev-optimum.\n",
    "- Az ablakos szűrő záró-tartománybeli csillapítása **monoton nő**, az áteresztő szél közelében pedig nagyobb az ingadozás — nem optimális.\n",
    "- Adott specifikációnál a Parks–McClellan szűrő **kisebb fokszámmal** is elérheti ugyanazt a teljesítményt, mint az ablakos.\n",
    "- A `signal.remez` súlyozási vektora ($1/\\delta_p$, $1/\\delta_s$) szabja meg, hogy melyik tartományban legyen kisebb a hiba.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "474741ee",
   "metadata": {},
   "source": [
    "## 8. Bilineáris transzformáció és frekvencia-elhajlás\n",
    "\n",
    "A tankönyv $(10.25)$ képlete:\n",
    "\n",
    "$$\n",
    "W(z) = W(s)\\Big|_{s = \\frac{2}{T_s}\\frac{z-1}{z+1}}\n",
    "$$\n",
    "\n",
    "Ez egy-egy értelmű leképezés a $j\\omega$ tengelyről az egységkörre: $\\omega_a = \\frac{2}{T_s}\\tan\\!\\left(\\frac{\\vartheta}{2}\\right)$. Ez a **frekvencia-elhajlás** (warping). Ezért tervezésnél a kritikus frekvenciákat **előtorzítani** (prewarp) kell.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "1. Indulj egy 4. rendű analóg Butterworth aluláteresztő $H(s)$-ből, $\\omega_c = 1$ rad/s vágási körfrekvenciával.\n",
    "2. Diszkretizáld bilineáris transzformációval $f_s = 100$ Hz mintavételi frekvencián, **előtorzítás nélkül**.\n",
    "3. Ismételd meg **előtorzítással**.\n",
    "4. Ábrázold mindkét eredményt és az eredeti analóg karakterisztikát — hasonlítsd össze a vágási frekvenciát.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0e21de89",
   "metadata": {},
   "outputs": [],
   "source": [
    "fs   = 100.0\n",
    "T_s  = 1.0/fs\n",
    "omega_c_des = 100.0   # kívánt analóg vágási [rad/s] (azaz f_c ≈ 15.9 Hz, ~1/3 Nyquist)\n",
    "\n",
    "# Analóg 4. rendű Butterworth LP, normalizálva omega_c_des-hez\n",
    "b_a, a_a = signal.butter(4, omega_c_des, btype='low', analog=True)\n",
    "\n",
    "# (1) Bilineáris transzformáció ELŐTORZÍTÁS NÉLKÜL\n",
    "b_d1, a_d1 = signal.bilinear(b_a, a_a, fs)\n",
    "\n",
    "# (2) Bilineáris ELŐTORZÍTÁSSAL: a digitális rendszerben az omega_c_des-t várjuk\n",
    "omega_c_prewarp = 2*fs*np.tan(omega_c_des/(2*fs))\n",
    "b_a_pre, a_a_pre = signal.butter(4, omega_c_prewarp, btype='low', analog=True)\n",
    "b_d2, a_d2 = signal.bilinear(b_a_pre, a_a_pre, fs)\n",
    "\n",
    "# Frekvenciaválaszok\n",
    "w_a = np.logspace(0, 3, 1024)\n",
    "_,    H_a  = signal.freqs(b_a,  a_a,  w_a)\n",
    "w_d, H_d1  = signal.freqz(b_d1, a_d1, worN=2048, fs=fs)\n",
    "_,   H_d2  = signal.freqz(b_d2, a_d2, worN=2048, fs=fs)\n",
    "w_d_rad = w_d * 2*np.pi  # rad/s\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(11, 4.5))\n",
    "\n",
    "ax.semilogx(w_a, 20*np.log10(np.maximum(np.abs(H_a), 1e-12)),\n",
    "            color=NAVY, lw=2, label=r'analóg $H(s)$')\n",
    "ax.semilogx(w_d_rad, 20*np.log10(np.maximum(np.abs(H_d1), 1e-12)),\n",
    "            color=CYAN, lw=1.5, label=r'bilineáris (NINCS előtorz.)')\n",
    "ax.semilogx(w_d_rad, 20*np.log10(np.maximum(np.abs(H_d2), 1e-12)),\n",
    "            color=AMBER, lw=1.5, ls='--', label=r'bilineáris (előtorzított)')\n",
    "\n",
    "ax.axvline(omega_c_des, color=RED, ls=':', alpha=0.6,\n",
    "           label=fr'kívánt $\\omega_c = {omega_c_des:.0f}$ rad/s')\n",
    "ax.axvline(np.pi*fs, color=MUTE, ls='--', alpha=0.4, label=f'Nyq. = {np.pi*fs:.0f} rad/s')\n",
    "\n",
    "ax.set_title('Frekvencia-elhajlás (warping) bilineáris transzformációnál')\n",
    "ax.set_xlabel('ω [rad/s]'); ax.set_ylabel('|H| [dB]')\n",
    "ax.set_ylim(-80, 5); ax.set_xlim(2, 400); ax.legend(loc='lower left')\n",
    "plt.tight_layout(); plt.show()\n",
    "\n",
    "# Mért -3 dB pontok\n",
    "def find_3dB(w, H):\n",
    "    H_dB = 20*np.log10(np.maximum(np.abs(H), 1e-12))\n",
    "    return w[np.argmin(np.abs(H_dB + 3))]\n",
    "\n",
    "print(f'Kívánt vágási:                {omega_c_des:.3f} rad/s')\n",
    "print(f'Mért, előtorzítás nélkül:    {find_3dB(w_d_rad, H_d1):.3f} rad/s   <- elhajlott!')\n",
    "print(f'Mért, előtorzítással:        {find_3dB(w_d_rad, H_d2):.3f} rad/s   <- pontos.')\n",
    "print(f'\\nElőtorzítási képlet:  ω_c* = (2/T_s)·tan(ω_c·T_s/2) = {omega_c_prewarp:.3f} rad/s')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "41ef30c5",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- **Előtorzítás nélkül**: a digitális vágási frekvencia eltolódik a kívánttól — az alacsony frekvenciás tartományban (ahol $\\tan\\vartheta/2 \\approx \\vartheta/2$) az eltérés kicsi, de a Nyquist-frekvenciához közeledve egyre nagyobb.\n",
    "- **Előtorzítással**: a $-3$ dB pont pontosan ott van, ahol kértük. Az analóg prototípus $H(s)$-t úgy módosítjuk, hogy az $\\omega_c^* = (2/T_s)\\tan(\\omega_c T_s/2)$ frekvencián vágjon, így a bilineáris transzformáció után a digitális rendszer az eredeti $\\omega_c$-n vág.\n",
    "- **Fontos**: az előtorzítás csak a **kritikus** (vágási, sávszéli) frekvenciákat állítja be pontosan — a teljes analóg karakterisztika $\\omega$-tengelye nemlineárisan torzul. Sávzáró és sáváteresztő tervezésénél mindkét sávszélt elő kell torzítani.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4a60f0b3",
   "metadata": {},
   "source": [
    "## 9. IIR prototípusok összehasonlítása\n",
    "\n",
    "A klasszikus analóg LP prototípusok más-más kompromisszumot kínálnak:\n",
    "\n",
    "| Prototípus | Áteresztő | Záró | Csoportkésleltetés |\n",
    "|---|---|---|---|\n",
    "| Butterworth | maximálisan lapos | monoton | közepes |\n",
    "| Chebyshev I | ekvi-ripple | monoton | közepes-rossz |\n",
    "| Chebyshev II | monoton | ekvi-ripple | közepes |\n",
    "| Elliptic (Cauer) | ekvi-ripple | ekvi-ripple | rossz |\n",
    "| Bessel | maximálisan lapos | monoton | maximálisan lapos |\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Tervezz négy különböző IIR LP szűrőt **azonos specifikációval**: $f_s = 1$ kHz, $f_p = 100$ Hz, $f_s = 150$ Hz, áteresztőbeli ripple max. 1 dB, záró-csillapítás min. 50 dB. Hasonlítsd össze:\n",
    "\n",
    "- a szükséges fokszámot,\n",
    "- az amplitúdóválaszt,\n",
    "- a csoportkésleltetést.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a5efc372",
   "metadata": {},
   "outputs": [],
   "source": [
    "fs   = 1000.0\n",
    "fp, fst = 100.0, 150.0\n",
    "gp, gst = 1.0, 50.0   # áteresztő ripple [dB], záró csill. [dB]\n",
    "\n",
    "# Megengedett rendszámok kiszámítása az adott specre\n",
    "N_butt, _ = signal.buttord(fp, fst, gp, gst, fs=fs)\n",
    "N_cby1, _ = signal.cheb1ord(fp, fst, gp, gst, fs=fs)\n",
    "N_cby2, _ = signal.cheb2ord(fp, fst, gp, gst, fs=fs)\n",
    "N_ell , _ = signal.ellipord(fp, fst, gp, gst, fs=fs)\n",
    "\n",
    "print(f'Szükséges fokszám:')\n",
    "print(f'  Butterworth:  {N_butt}')\n",
    "print(f'  Chebyshev I:  {N_cby1}')\n",
    "print(f'  Chebyshev II: {N_cby2}')\n",
    "print(f'  Elliptic:     {N_ell}')\n",
    "\n",
    "# Tervezés (sos = second-order sections - kaszkád biquad)\n",
    "sos_butt = signal.butter (N_butt, fp, btype='low', fs=fs, output='sos')\n",
    "sos_cby1 = signal.cheby1 (N_cby1, gp, fp, btype='low', fs=fs, output='sos')\n",
    "sos_cby2 = signal.cheby2 (N_cby2, gst, fst, btype='low', fs=fs, output='sos')\n",
    "sos_ell  = signal.ellip  (N_ell, gp, gst, fp, btype='low', fs=fs, output='sos')\n",
    "\n",
    "prototypes = [\n",
    "    ('Butterworth', sos_butt, CYAN),\n",
    "    ('Chebyshev I', sos_cby1, AMBER),\n",
    "    ('Chebyshev II', sos_cby2, NAVY),\n",
    "    ('Elliptic',    sos_ell,  RED),\n",
    "]\n",
    "\n",
    "fig, axes = plt.subplots(1, 2, figsize=(13, 4.5))\n",
    "\n",
    "for name, sos, col in prototypes:\n",
    "    w_n, H = signal.sosfreqz(sos, worN=2048, fs=fs)\n",
    "    H_dB = 20*np.log10(np.maximum(np.abs(H), 1e-12))\n",
    "    axes[0].plot(w_n, H_dB, color=col, label=name)\n",
    "\n",
    "    # Csoportkésleltetés (sos -> ba konverzió)\n",
    "    b, a = signal.sos2tf(sos)\n",
    "    _, gd = signal.group_delay((b, a), w=2048, fs=fs)\n",
    "    axes[1].plot(w_n, gd, color=col, label=name)\n",
    "\n",
    "# Spec-vonalak\n",
    "axes[0].axhline(-gst, color=MUTE, ls='--', alpha=0.4)\n",
    "axes[0].axhline(-gp,  color=MUTE, ls=':', alpha=0.4)\n",
    "axes[0].axvline(fp,   color=MUTE, ls=':', alpha=0.4)\n",
    "axes[0].axvline(fst,  color=MUTE, ls=':', alpha=0.4)\n",
    "axes[0].set_title('Amplitúdóválasz (dB)')\n",
    "axes[0].set_xlabel('f [Hz]'); axes[0].set_ylabel('|H| [dB]')\n",
    "axes[0].set_xlim(0, 250); axes[0].set_ylim(-80, 3); axes[0].legend()\n",
    "\n",
    "axes[1].set_title('Csoportkésleltetés [minta]')\n",
    "axes[1].set_xlabel('f [Hz]'); axes[1].set_ylabel(r'$\\tau_g$')\n",
    "axes[1].set_xlim(0, 250); axes[1].set_ylim(0, 50); axes[1].legend()\n",
    "\n",
    "plt.tight_layout(); plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "14ffc451",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- Az **elliptikus** szűrő érte el a specifikációt a **legkisebb fokszámmal** — adott áteresztő-ripple és záró-csillapítás mellett ez az optimum.\n",
    "- A **Butterworth** sima áteresztő tartományt és monoton karakterisztikát ad, de a legmagasabb fokszámot igényli.\n",
    "- A **Chebyshev I** áteresztőben hullámzik (ezért érhető el alacsonyabb fokszám), záróban monoton.\n",
    "- A **Chebyshev II** fordított: áteresztőben sima, záróban hullámzik.\n",
    "- **Csoportkésleltetés**: a Butterworth viszonylag lapos, az elliptikus erősen ingadozik a vágási frekvencia közelében — hullámformára érzékeny alkalmazásnál ez fontos szempont. Ha a fáziskésleltetés a kritikus, a Bessel-prototípust érdemes választani (ott a $|H|$ rosszabb, viszont $\\tau_g$ közel konstans).\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a30febe3",
   "metadata": {},
   "source": [
    "## 10. Realizációs struktúrák — direct form vs biquad cascade\n",
    "\n",
    "Magas fokszámú IIR szűrő **direct form** (egyetlen $W(z)$ tört) realizációja **numerikusan instabil** lehet: a koefficiensek véges szóhosszra történő kerekítése a pólusokat az egységkörhöz közel mozdíthatja, a szűrő tényleges karakterisztikája erősen torzulhat.\n",
    "\n",
    "A megoldás: **second-order sections** (SOS / biquad kaszkád):\n",
    "\n",
    "$$\n",
    "W(z) = \\prod_{i=1}^{L} \\frac{b_{0i} + b_{1i}z^{-1} + b_{2i}z^{-2}}{1 + a_{1i}z^{-1} + a_{2i}z^{-2}}\n",
    "$$\n",
    "\n",
    "Minden szekció pólus-zérus pár csak helyileg érzékeny, így a numerikus érzékenység drámaian csökken.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Tervezz egy 8. rendű elliptikus IIR LP szűrőt, és vizsgáld meg:\n",
    "\n",
    "1. a koefficiensek dinamikatartományát direct form vs SOS realizációban,\n",
    "2. a frekvenciaválasz torzulását, ha a koefficienseket `float16` pontossággal kvantáljuk.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "14f078d3",
   "metadata": {},
   "outputs": [],
   "source": [
    "fs = 1000.0\n",
    "N  = 8\n",
    "\n",
    "# Tervezés direct form (ba) és SOS módon\n",
    "b_df, a_df = signal.ellip(N, 0.5, 60, 100, btype='low', fs=fs, output='ba')\n",
    "sos        = signal.ellip(N, 0.5, 60, 100, btype='low', fs=fs, output='sos')\n",
    "\n",
    "print(f'Direct form együtthatók:')\n",
    "print(f'  b: range = [{b_df.min():.3e}, {b_df.max():.3e}]   max/min: {b_df.max()/abs(b_df).min():.1e}')\n",
    "print(f'  a: range = [{a_df.min():.3e}, {a_df.max():.3e}]   max/min: {a_df.max()/abs(a_df).min():.1e}')\n",
    "print(f'\\nSOS együtthatók (max abs. value/szekció):')\n",
    "for i, sect in enumerate(sos):\n",
    "    print(f'  szekció {i+1}:  max|coef| = {np.abs(sect).max():.3f}')\n",
    "\n",
    "# Float16 kvantálás\n",
    "def quantize(x, dtype=np.float16):\n",
    "    return x.astype(dtype).astype(np.float64)\n",
    "\n",
    "b_q  = quantize(b_df)\n",
    "a_q  = quantize(a_df)\n",
    "sos_q = quantize(sos)\n",
    "\n",
    "# Frekvenciaválaszok\n",
    "w, H_orig    = signal.sosfreqz(sos, worN=2048, fs=fs)\n",
    "_, H_df_q    = signal.freqz(b_q, a_q, worN=2048, fs=fs)\n",
    "_, H_sos_q   = signal.sosfreqz(sos_q, worN=2048, fs=fs)\n",
    "\n",
    "fig, axes = plt.subplots(1, 2, figsize=(13, 4.5))\n",
    "\n",
    "axes[0].plot(w, 20*np.log10(np.maximum(np.abs(H_orig), 1e-12)),\n",
    "             color=NAVY, lw=2, label='float64 referencia')\n",
    "axes[0].plot(w, 20*np.log10(np.maximum(np.abs(H_df_q), 1e-12)),\n",
    "             color=RED, lw=1.2, label='float16, direct form')\n",
    "axes[0].plot(w, 20*np.log10(np.maximum(np.abs(H_sos_q), 1e-12)),\n",
    "             color=CYAN, lw=1.2, ls='--', label='float16, SOS biquad')\n",
    "axes[0].set_title('Kvantálás hatása a frekvenciaválaszra')\n",
    "axes[0].set_xlabel('f [Hz]'); axes[0].set_ylabel('|H| [dB]')\n",
    "axes[0].set_xlim(0, 250); axes[0].set_ylim(-100, 5); axes[0].legend()\n",
    "\n",
    "# Pólus-elmozdulás\n",
    "z_orig, p_orig, _ = signal.sos2zpk(sos)\n",
    "b_full, a_full = signal.sos2tf(sos_q)\n",
    "p_q_sos = np.roots(a_full)\n",
    "p_q_df  = np.roots(a_q)\n",
    "\n",
    "theta_circ = np.linspace(0, 2*np.pi, 200)\n",
    "axes[1].plot(np.cos(theta_circ), np.sin(theta_circ), '--', color=MUTE, lw=1)\n",
    "axes[1].plot(p_orig.real,  p_orig.imag,  'o', color=NAVY, ms=10, mfc='none', mew=2, label='eredeti pólus')\n",
    "axes[1].plot(p_q_sos.real, p_q_sos.imag, 'x', color=CYAN, ms=10, mew=2, label='SOS, kvant.')\n",
    "axes[1].plot(p_q_df.real,  p_q_df.imag,  '+', color=RED,  ms=12, mew=2, label='direct, kvant.')\n",
    "axes[1].axhline(0, color=MUTE, lw=0.5); axes[1].axvline(0, color=MUTE, lw=0.5)\n",
    "axes[1].set_xlim(-1.2, 1.2); axes[1].set_ylim(-1.2, 1.2); axes[1].set_aspect('equal')\n",
    "axes[1].set_title('Pólusok elmozdulása kvantáláskor')\n",
    "axes[1].set_xlabel('Re(z)'); axes[1].set_ylabel('Im(z)'); axes[1].legend()\n",
    "\n",
    "plt.tight_layout(); plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d65155c6",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A direct form együtthatók dinamikatartománya **több nagyságrendet** átfog (kis és nagy értékek vegyesen) — kvantáláskor a kisebb értékek elvesznek.\n",
    "- Az SOS minden szekció koefficiensei **közel azonos nagyságrendűek** ($|\\cdot| \\le 2$) — ez ideális fixpontos vagy kis pontosságú lebegőpontos aritmetikához.\n",
    "- Float16 kvantálás után a direct form karakterisztikája **drámaian torzul** — a záró-tartománybeli csillapítás eltűnik, a karakterisztika hullámzik. Az SOS realizáció gyakorlatilag változatlan.\n",
    "- A pólusok látványosan elmozdulnak direct form esetén, míg SOS-ben az elmozdulás minimális.\n",
    "\n",
    "**Gyakorlati következmény**: az `scipy.signal.butter`/`cheby1`/stb. függvényeket mindig `output='sos'` paraméterrel hívd, és `signal.sosfilt`-et használj `signal.lfilter` helyett.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aeddd679",
   "metadata": {},
   "source": [
    "## 11. Kvantálási hatások és limit cycle\n",
    "\n",
    "Véges szóhossz mellett a következő hatásokkal kell számolni:\n",
    "\n",
    "1. **ADC kvantálási zaj** ($\\sigma_e^2 = q^2/12$, $q$ a kvantum)\n",
    "2. **Koefficiens-kvantálás** (pólusok elmozdulása, lásd 10. gyakorlat)\n",
    "3. **Aritmetikai kerekítési zaj** (MAC-műveletek után)\n",
    "4. **Túlcsordulás**\n",
    "\n",
    "IIR rendszerek specifikuma a **limit cycle**: véges szóhossz mellett nulla bemenetnél is megjelenhet egy kis amplitúdójú periodikus oszcilláció (granular limit cycle).\n",
    "\n",
    "### Feladat\n",
    "\n",
    "Demonstráld a limit cycle jelenséget egy egyszerű 2. rendű IIR rezonátorral:\n",
    "\n",
    "$$\n",
    "y[k] = s[k] + 1{,}9\\,y[k-1] - 0{,}95\\,y[k-2]\n",
    "$$\n",
    "\n",
    "A rendszer egy konjugált pólus-párral rendelkezik az egységkörhöz közel. Szimuláld a választ kvantáláson keresztül 8 bites aritmetikával, nem nulla kezdeti feltételekkel és nulla bemenetre. Hasonlítsd össze a `float64` referenciával.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cdfa59a7",
   "metadata": {},
   "outputs": [],
   "source": [
    "def lfilter_quantized(b, a, x, n_bits, init=None):\n",
    "    \"\"\"IIR lfilter szimulálva fixpontos kerekítéssel.\"\"\"\n",
    "    Q = 2**(n_bits - 1) - 1  # max ábrázolható szám (előjeles)\n",
    "    def quant(v):\n",
    "        return np.round(v * Q) / Q\n",
    "\n",
    "    N_in = len(x)\n",
    "    y = np.zeros(N_in)\n",
    "    # Direct form II transposed-szerű, de kvantálással minden MAC után\n",
    "    M = max(len(b), len(a)) - 1\n",
    "    state = np.zeros(M) if init is None else init.copy()\n",
    "\n",
    "    for n in range(N_in):\n",
    "        y_n = b[0]*x[n] + state[0]\n",
    "        y_n = quant(y_n)\n",
    "        for i in range(M-1):\n",
    "            state[i] = quant(b[i+1]*x[n] - a[i+1]*y_n + state[i+1])\n",
    "        state[M-1] = quant(b[M]*x[n] - a[M]*y_n) if len(b) > M else quant(-a[M]*y_n)\n",
    "        y[n] = y_n\n",
    "    return y\n",
    "\n",
    "# 2. rendű rezonátor pólus-párja\n",
    "b = np.array([1.0, 0.0, 0.0])\n",
    "a = np.array([1.0, -1.9, 0.95])\n",
    "poles = np.roots(a)\n",
    "print(f'Pólusok: {poles},  |p| = {np.abs(poles).round(4)}')\n",
    "\n",
    "# Kezdeti feltétel: y[-1] = 0.5, y[-2] = 0\n",
    "N = 200\n",
    "x = np.zeros(N)\n",
    "\n",
    "# float64 referencia, kezdeti állapottal\n",
    "zi = signal.lfilter_zi(b, a) * 0.5    # y[-1] ~ 0.5\n",
    "y_ref, _ = signal.lfilter(b, a, x, zi=zi)\n",
    "\n",
    "# 8 bites kvantálás\n",
    "y_q8  = lfilter_quantized(b, a, x, n_bits=8,  init=np.array([0.5, 0.0]))\n",
    "y_q12 = lfilter_quantized(b, a, x, n_bits=12, init=np.array([0.5, 0.0]))\n",
    "y_q16 = lfilter_quantized(b, a, x, n_bits=16, init=np.array([0.5, 0.0]))\n",
    "\n",
    "n = np.arange(N)\n",
    "fig, ax = plt.subplots(figsize=(11, 4.5))\n",
    "ax.plot(n, y_ref, color=NAVY,  lw=1.0, label='float64 referencia')\n",
    "ax.plot(n, y_q16, color=AMBER, lw=1.0, alpha=0.8, label='16 bit')\n",
    "ax.plot(n, y_q12, color=CYAN,  lw=1.0, alpha=0.8, label='12 bit')\n",
    "ax.plot(n, y_q8,  color=RED,   lw=1.0, alpha=0.8, label='8 bit')\n",
    "ax.set_title('Limit cycle különböző szóhosszaknál (nulla bemenet, kezdeti feltétellel)')\n",
    "ax.set_xlabel('n'); ax.set_ylabel('y[n]'); ax.legend()\n",
    "plt.tight_layout(); plt.show()\n",
    "\n",
    "# Maradó energia (a refnek nullához kellene tartania)\n",
    "print(f'\\nMaradó energia n>100 utáni szakaszon:')\n",
    "print(f'  float64:  {np.sum(y_ref[100:]**2):.2e}')\n",
    "print(f'  16 bit:   {np.sum(y_q16[100:]**2):.2e}')\n",
    "print(f'  12 bit:   {np.sum(y_q12[100:]**2):.2e}')\n",
    "print(f'   8 bit:   {np.sum(y_q8[100:]**2):.2e}')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "40e2d78a",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A `float64` referencia exponenciálisan csillapodik nullához (a pólusok az egységkörön belül vannak).\n",
    "- 16 és 12 biten is még jól közelíti az elvárt választ.\n",
    "- **8 biten** a rendszer **nem áll meg nullában** — egy kis amplitúdójú **periodikus oszcilláció** marad. Ez a granular limit cycle: a kvantálási kerekítés egy ön-fenntartó körfolyamatba löki a rendszert.\n",
    "- Limit cycle ellen védekezés:\n",
    "  - Nagyobb szóhossz (egyenes út, de drága).\n",
    "  - **Speciális struktúrák** (wave digital filter, lattice).\n",
    "  - **Telítéses aritmetika** és kerekítés helyett **csonkolás** (de utóbbi torzít).\n",
    "  - **Dither** hozzáadása a kvantálás előtt.\n",
    "- A FIR szűrőkben **nincs** limit cycle, mert nincs visszacsatolás.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bbb7cd03",
   "metadata": {},
   "source": [
    "## 12. Gyakorlati alkalmazás — zajos jel szűrése\n",
    "\n",
    "Egy szintetikus, de realisztikus szűrési feladaton keresztül összerakjuk az eddigieket. Egy érzékelő $f_s = 1000$ Hz mintavételi frekvenciával rögzít egy hasznos jelet (10 Hz + 30 Hz szinusz keverék), amelyre **két típusú zaj** szuperponálódik:\n",
    "\n",
    "- széles sávú Gauss-zaj,\n",
    "- $50$ Hz hálózati interferencia.\n",
    "\n",
    "A feladat egy aluláteresztő + notch-szűrő láncolat tervezése, és az SNR-javulás kvantitatív értékelése.\n",
    "\n",
    "### Feladat\n",
    "\n",
    "1. Generálj egy $N = 2000$ mintából álló zajos jelet.\n",
    "2. Tervezz egy aluláteresztő FIR szűrőt $f_c = 40$ Hz vágási frekvenciával (Hamming-ablak).\n",
    "3. Tervezz egy 50 Hz-es notch IIR szűrőt (`signal.iirnotch`).\n",
    "4. Kapcsold sorba a két szűrőt és ábrázold a bemeneti, közbenső és végső kimeneti jelet, valamint a spektrumokat.\n",
    "5. Mérd az SNR-javulást.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "389ec40c",
   "metadata": {},
   "outputs": [],
   "source": [
    "rng = np.random.default_rng(42)\n",
    "fs  = 1000.0\n",
    "T   = 2.0\n",
    "N   = int(T*fs)\n",
    "t   = np.arange(N)/fs\n",
    "\n",
    "# Hasznos jel\n",
    "clean = np.sin(2*np.pi*10*t) + 0.7*np.sin(2*np.pi*30*t + 0.4)\n",
    "\n",
    "# Zaj\n",
    "gauss_noise = 0.5 * rng.standard_normal(N)\n",
    "mains       = 0.6 * np.sin(2*np.pi*50*t)\n",
    "noise       = gauss_noise + mains\n",
    "\n",
    "x = clean + noise\n",
    "\n",
    "# (1) FIR LP szűrő, f_c = 40 Hz, Hamming\n",
    "N_fir = 101\n",
    "h_fir = signal.firwin(N_fir, 40, window='hamming', fs=fs)\n",
    "\n",
    "# (2) Notch IIR szűrő 50 Hz-en, Q = 30\n",
    "b_notch, a_notch = signal.iirnotch(50, Q=30, fs=fs)\n",
    "\n",
    "# Sorba kapcsolás: x -> LP -> notch -> y\n",
    "x_lp = signal.lfilter(h_fir, 1.0, x)\n",
    "y    = signal.lfilter(b_notch, a_notch, x_lp)\n",
    "\n",
    "# Csoport-késleltetés-kompenzáció a vizualizációhoz: a FIR (N-1)/2 minta késést okoz\n",
    "delay = (N_fir-1)//2\n",
    "\n",
    "# --- Ábrázolás: idő- és spektrumtartomány ---\n",
    "fig, axes = plt.subplots(3, 2, figsize=(13, 8))\n",
    "\n",
    "# Idő\n",
    "for ax, sig_, lab, col in zip(axes[:,0],\n",
    "                              [x, x_lp, y],\n",
    "                              ['Bemenet (zajos)', 'LP szűrő után',\n",
    "                               'LP + notch után'],\n",
    "                              [MUTE, CYAN, NAVY]):\n",
    "    ax.plot(t, sig_, color=col, lw=0.8)\n",
    "    ax.plot(t, clean, color=AMBER, ls='--', alpha=0.5, label='referencia (clean)')\n",
    "    ax.set_title(lab); ax.set_xlim(0, 0.5); ax.set_xlabel('t [s]'); ax.set_ylabel('amp')\n",
    "    ax.legend(loc='upper right', fontsize=8)\n",
    "\n",
    "# Spektrum\n",
    "def plot_spectrum(ax, x, fs, lab, col):\n",
    "    f, Pxx = signal.welch(x, fs=fs, nperseg=512)\n",
    "    ax.semilogy(f, Pxx, color=col, label=lab)\n",
    "    ax.set_xlabel('f [Hz]'); ax.set_ylabel('PSD'); ax.set_xlim(0, 100)\n",
    "    ax.legend()\n",
    "\n",
    "plot_spectrum(axes[0,1], x,    fs, 'Bemenet',          MUTE)\n",
    "plot_spectrum(axes[1,1], x_lp, fs, 'LP után',          CYAN)\n",
    "plot_spectrum(axes[2,1], y,    fs, 'LP + notch után',  NAVY)\n",
    "\n",
    "# Spektrum-jelölők\n",
    "for ax in axes[:,1]:\n",
    "    ax.axvline(10, color=AMBER, ls=':', alpha=0.5)\n",
    "    ax.axvline(30, color=AMBER, ls=':', alpha=0.5)\n",
    "    ax.axvline(50, color=RED, ls=':', alpha=0.6, label='50 Hz')\n",
    "\n",
    "plt.tight_layout(); plt.show()\n",
    "\n",
    "# SNR-mérés\n",
    "def snr(signal_, reference, delay=0):\n",
    "    if delay > 0:\n",
    "        signal_ = signal_[delay:]\n",
    "        reference = reference[:len(signal_)]\n",
    "    err = signal_ - reference\n",
    "    return 10*np.log10(np.sum(reference**2) / np.sum(err**2))\n",
    "\n",
    "print(f'SNR a bemeneten:                {snr(x,    clean):.2f} dB')\n",
    "print(f'SNR LP után (késleltetéskomp.): {snr(x_lp, clean, delay=delay):.2f} dB')\n",
    "print(f'SNR LP+notch után:              {snr(y,    clean, delay=delay):.2f} dB')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3b668eda",
   "metadata": {},
   "source": [
    "**Megfigyelések.**\n",
    "\n",
    "- A bemeneti spektrumon jól látható mind a Gauss-zaj „szőnyege\", mind a 50 Hz-es szpájk.\n",
    "- Az **aluláteresztő FIR szűrő** levágja az 50 Hz feletti komponenseket, ami a Gauss-zaj jelentős részét eltávolítja, és a hálózati interferenciát is csökkenti — de ez utóbbi kis részben átszivárog a 40 Hz-es vágás miatti tranziens átmeneten.\n",
    "- A **notch IIR szűrő** specifikusan a 50 Hz-es komponenst nyomja el, megőrizve a környező frekvenciákat.\n",
    "- A **sorba kapcsolt** lánc kimenete jelentős SNR-javulást ad (a számszerű érték a generálási seedtől függ, nálunk kb. **20 dB nagyságrendű**).\n",
    "- A FIR LP szűrő $(N-1)/2 = 50$ minta lineáris fáziskésést okoz — ezt vizuális összehasonlításnál kompenzálni kell.\n",
    "\n",
    "**További gyakorlatok**:\n",
    "- Nézd meg, mi történik, ha a notch-szűrő $Q$ tényezőjét csökkented/növeled — szélesebb/keskenyebb elnyomó sáv.\n",
    "- Cseréld ki a FIR-t egy hasonló specifikációjú IIR-re és vizsgáld a **fáziskésést**.\n",
    "- Próbálj **filtfilt** zéró-fázisú szűrést — időtartományban nulla késés, cserébe duplán fut a szűrő.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "08455bcf",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Összefoglalás és továbblépés\n",
    "\n",
    "A 12 gyakorlatban végigmentünk a digitális szűrők elméletének és gyakorlatának teljes spektrumán:\n",
    "\n",
    "| Téma | Gyakorlat |\n",
    "|---|---|\n",
    "| Időtartomány | 1 (konvolúció, differenciaegyenlet) |\n",
    "| Frekvenciatartomány | 2, 4, 5 |\n",
    "| $z$-tartomány, stabilitás | 3 |\n",
    "| FIR-tervezés | 4, 5, 6, 7 |\n",
    "| IIR-tervezés | 8, 9 |\n",
    "| Realizáció és numerika | 10, 11 |\n",
    "| Alkalmazás | 12 |\n",
    "\n",
    "### Javasolt további gyakorlatok\n",
    "\n",
    "1. **Adaptív szűrő** (LMS, NLMS): hangos környezetben echo-csökkentés.\n",
    "2. **2D-szűrés**: képi LP/HP/Sobel/Laplace operátorok.\n",
    "3. **Polifázis dekimáció / interpoláció**: hatékony mintavételiráta-átalakítás.\n",
    "4. **Sávszűrő-bank**: oktávsávos akusztikai analizátor.\n",
    "5. **Wiener-szűrő**: optimális zajeltávolítás statisztikai értelemben.\n",
    "6. **Fáziseltolódás-kompenzáció**: `filtfilt` és Hilbert-transzformáció.\n",
    "\n",
    "### Hivatkozások\n",
    "\n",
    "- Kuczmann Miklós: *Jelek és rendszerek*, Széchenyi István Egyetem (2007), 7–10. fejezet.\n",
    "- Oppenheim & Schafer: *Discrete-Time Signal Processing*, 3rd ed.\n",
    "- Proakis & Manolakis: *Digital Signal Processing*, 4th ed.\n",
    "- SciPy Signal Processing dokumentáció: <https://docs.scipy.org/doc/scipy/reference/signal.html>\n"
   ]
  }
 ],
 "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.11"
  },
  "title": "Digitális szűrők — gyakorlatok"
 },
 "nbformat": 4,
 "nbformat_minor": 5
}