21-10-04 L06.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 06 - Simulation Theorem\n",
    "Lecture 06 - 04 October 2021\n",
    "\n",
    "## Simulation of an analog system by means of a digital one\n",
    "\n",
    "An interesting question is whether an analog system can be simulated by a digital one. The question is not as trivial as it sounds because a digital system deals only with discrete samples while an analog one has a continuou dynamic. A first intuition may be: \"as Nyquist-Shannon allows to perfectly reconstruct a signal from its samples, all the information about the system must be contained in the samples and performing operations on them should be sufficient\". Formally. Given two system as shown in the following image:\n",
    "\n",
    "![Simulation of an analog system by means of a digital one](images/21-10-04_simulation.png)\n",
    "\n",
    "does such $h[n] : y'[n]=\\sum_{k}h[n-k]x_s[k]$ exists? And, if so, how is it possible to find it? To answer this question let's start with a function $y'[n]=\\sum_{k}h[n-k]x_s[k]$ and make some assumptions:\n",
    "\n",
    "### 1a) Both $x_A(t), G(t) \\in \\mathcal{L}^2$.\n",
    "If this is true, it is possible to express the sampled output signal via its inverse Fuorier transform:\n",
    "\n",
    "\\begin{align*}\n",
    "y'[n]=\\sum_k h[n-k]x_A(kT)= \\sum_k h[n-k] \\frac{1}{2\\pi}\\int_{-\\infty}^{+\\infty} \\tilde{x}_A(\\omega)e^{-i\\omega k T} d\\omega\n",
    "\\end{align*}\n",
    "\n",
    "### 2b) The simulator $h[n]$ is BIBO Stable\n",
    "If this is the case, it is possible to swap the summation and the integral in $\\sum_k \\leftrightarrow \\int_{-\\infty}^{+\\infty}$ without dealing with the convergence:\n",
    "\\begin{align*}\n",
    "y'[n] &= \\frac{1}{2\\pi}\\int_{-\\infty}^{+\\infty} \\tilde{x}_A(\\omega) \\sum_k h[n-k] e^{i\\omega (n-k) T} e^{-i\\omega n T} d\\omega \\qquad \\text{define } z=e^{-i\\omega T} \\quad m=n-k\\\\\n",
    "&= \\frac{1}{2\\pi}\\int_{-\\infty}^{+\\infty} \\sum_m h[m] z^{-m} e^{-i\\omega n T} \\tilde{x}_A(\\omega) d\\omega \\\\\n",
    "&= \\frac{1}{2\\pi}\\int_{-\\infty}^{+\\infty} H\\left( z=e^{-i\\omega Tn}\\right) \\tilde{x}_A(\\omega) d\\omega \n",
    "\\end{align*}\n",
    "\n",
    "where the condition $z=e^{-i\\omega T}$ was introduced to obtain the z-transform in conjunction with the condition $z\\in \\Gamma_1 \\subset ROC$ and the BIBO stability of the system.\n",
    "\n",
    "Now, consider the sampling of the output signal $y_A(t)$:\n",
    "\\begin{align*}\n",
    "\\require{cancel}\n",
    "y_A[n]=y_A(nT)=\\frac{1}{2\\pi} \\int_{-\\infty}^{\\infty} e^{-inT\\omega} \\cancelto{\\tilde{y}(\\omega)}{\\tilde{G}(\\omega)\\tilde{x}(\\omega)} d\\omega\n",
    "\\end{align*}\n",
    "\n",
    "In order to satisfy the requirement $y'[n]=y_A[n]$ we need:\n",
    "\\begin{align*}\n",
    "y'[n]-y_A[n]=0 \\iff \\frac{1}{2\\pi} \\int_{-\\infty}^{\\infty} e^{-inT\\omega} \\tilde{x}(\\omega) \\left[ H\\left( z=e^{-i\\omega Tn}\\right) -\\tilde{G}(\\omega)\\right] d\\omega = 0 \n",
    "\\end{align*}\n",
    "\n",
    "The easiest solution to that equation is $H\\left( z=e^{-i\\omega Tn}\\right)=\\tilde{G}(\\omega)$ but it presents many problems, the most important of which is that _it requires $\\tilde{G}(\\omega)$ to be periodic with period T_ as $H\\left( z=e^{-i\\omega Tn}\\right)$ satisfies that property by default.\n",
    "A workaround is to force the whole integral to be zero with an assumption:\n",
    "\n",
    "### Assumption 2: The input signal $\\tilde{x}(\\omega)$ is $\\frac{2\\pi}{T}$ band limited\n",
    "this forces the equality $H\\left( z=e^{-i\\omega Tn}\\right)=\\tilde{G}(\\omega)$ only on $\\left[-\\frac{\\pi}{T};\\frac{\\pi}{T} \\right]$.\n",
    "\n",
    "The discrete transfer function can be found by inverting the z-transform:\n",
    "\n",
    "\\begin{align*}\n",
    "h[n] &=\\frac{1}{2\\pi i} \\oint_{\\Gamma um 0\\\\ \\Gamma\\in ROC} H(z)z^{n-1}dz \\qquad z=e^{-i\\omega T} \\iff \\Gamma=\\Gamma_1\\\\\n",
    "&= -\\frac{1}{2\\pi i} \\int_{-\\frac{\\pi}{T}}^{\\frac{\\pi}{T}} H(z=e^{-i\\omega T})\\left(e^{-i\\omega T}\\right)^{n-1} e^{-i\\omega T} (-iT)d\\omega\\\\\n",
    "&= \\frac{T}{2\\pi} \\int_{-\\frac{\\pi}{T}}^{\\frac{\\pi}{T}} e^{-i\\omega nT}\\tilde{G}(\\omega) d\\omega\n",
    "\\end{align*}\n",
    "\n",
    "#### Properties:\n",
    "1. $G(t)\\in\\mathcal{R} \\iff h[n]\\in\\mathcal{R}$, proof by substitution.\n",
    "2. BIBO STABILITY: must be checked on occurrence but it's generally true\n",
    "3. CAUSALITY IS NOT CONSERVED\n",
    "\n",
    "## Examples\n",
    "### First order low pass filter\n",
    "\\begin{align*}\n",
    "\\require{mathtools}\n",
    "h[n]=\\frac{T}{2\\pi}\\int_{-\\frac{\\pi}{T}}^{\\frac{\\pi}{T}}e^{-i\\omega nT}\\frac{1}{1-i\\omega \\tau} d\\omega \\xrightarrow[\\text{Perfect sampling}]{T\\rightarrow 0} \\frac{T}{2\\pi}\\int_{-\\infty}^{\\infty}\\frac{e^{-i\\omega nT}}{1-i\\omega \\tau} d\\omega = \\frac{\\theta(t)}{\\tau}e^{-t/\\tau}\n",
    "\\end{align*}\n",
    "but if $T\\not \\rightarrow 0$ there is NO WAY to compute the integral.\n",
    "![Not standard integral](images/21-10-04_not_standard_integral.png)\n",
    "\n",
    "### The differentiator\n",
    "The differentiator's equation is $\\tau \\frac{d}{dt}\\rightarrow -i\\omega\\tau=\\tilde{G}(\\omega)$\n",
    "\\begin{align*}\n",
    "\\require{mathtools}\n",
    "h[n]&=\\frac{T}{2\\pi}\\int_{-\\frac{\\pi}{T}}^{\\frac{\\pi}{T}}(-i\\omega \\tau)e^{-i\\omega nT} d\\omega \\qquad x=\\omega T\\\\\n",
    "&= - \\frac{i}{2\\pi T}\\int_{-\\pi}^{\\pi}e^{i nx} x dx = \\left. \\begin{cases}\n",
    "0 \\text{ if } n=0\\\\\n",
    "\\frac{(-1)^n}{nT} \\text{ otherwise}\n",
    "\\end{cases} \\right\\} = (-1)^n \\frac{1-\\delta[n]}{nT}\n",
    "\\end{align*}\n",
    "\n",
    "## Backward interpretation of the simulation theorem\n",
    "The one used above was the straightforward interpretation of the theorem. To overcome some of the problems it is possible to consider a **backward intepretation of the simulation theorem**, that is, the requirement are relaxed to $y'[n]\\simeq y[n]$ so that $\\tilde{G}(\\omega)=H(z=e^{-i\\omega T}) + o(\\omega^n)$ in $\\left[ -\\frac{\\pi}{T}; \\frac{\\pi}{T} \\right]$.\n",
    "\n",
    "### e.g. The differentiator (again)\n",
    "Let's take the approximate transfer function (which is given/found by an intuition): $V[n]=\\delta[n]-\\delta[n-1]$ with z-transform $V(z)=1 - \\frac{1}{z}$. By taking the limit $\\omega T << 1$, it is possible to write:\n",
    "\\begin{align*}\n",
    "\\tilde{G}(\\omega)=V(z=e^{-i\\omega T})=1-e^{i\\omega T} \\xrightarrow[]{\\omega T \\rightarrow0}\n",
    " -i\\omega T + \\mathcal{o}\\left( (\\omega T)^2 \\right)\n",
    "\\end{align*}\n",
    "\n",
    "It is even possible to improve the accuracy by introducing new terms:\n",
    "\\begin{align*}\n",
    "V[n]&= \\alpha \\delta[n]+\\beta \\delta[n-1]+\\gamma\\delta[n-2] \\implies V(z)=\\alpha+\\frac{\\beta}{z}+\\frac{\\gamma}{z^2}\\\\\n",
    "\\tilde{V}(\\omega)&=V(z=e^{-i\\omega T}) = \\alpha + \\beta e^{i\\omega T} + \\gamma e^{2i\\omega T} \\\\\n",
    "&=\\alpha+\\beta+\\gamma + (\\beta+2\\gamma)i\\omega T + (\\beta+4\\gamma)(i\\omega T)^2 + \\mathcal{o}\\left( (\\omega T)^3 \\right)\n",
    "\\end{align*}\n",
    "\n",
    "By comparing the approximated solution with the analitic one $-i\\omega T$ it is possible to write a system of equations and solve it:\n",
    "\\begin{align*}\n",
    "    \\begin{cases}\n",
    "        \\alpha+\\beta+\\gamma=0 \\\\\n",
    "        \\beta+2\\gamma=-1\\\\\n",
    "        \\beta+4\\gamma=0\n",
    "    \\end{cases} \\implies \\begin{cases}\n",
    "        \\alpha = \\frac{3}{2} \\\\\n",
    "        \\beta =-2\\\\\n",
    "        \\gamma=\\frac{1}{2}\n",
    "    \\end{cases}\n",
    "\\end{align*}\n",
    "\n",
    "![](images/21-10-04_differentiator.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "celltoolbar": "Edit 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.8.10"
  },
  "widgets": {
   "application/vnd.jupyter.widget-state+json": {
    "state": {},
    "version_major": 2,
    "version_minor": 0
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}