[WIP] Overlap-Add math

docs/overlap_add_math.ipynb

@@ -0,0 +1,222 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "e984e902-bf76-43ae-a5d5-55e59357b6bf",
+   "metadata": {},
+   "source": [
+    "# From Linear Convolution to Overlap-Add Method!"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e9725572-a39b-4fb5-992c-5f3d22a80797",
+   "metadata": {},
+   "source": [
+    "This notebook explains the math behind the overlap-add method in convolution. The content is based on the book \"Discrete-Time Signal Processing\" by Alan V. Oppenheim and Ronald W. Schafer. The book is [publicly available in MIT OpenCourseWare](https://ocw.mit.edu/courses/res-6-dtsp-discrete-time-signal-processing/resources/mitres_6-dtsp_s26_thirdedition_pdf/)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "03f9d2f8-1560-4a5d-93a2-ec82fb422f6f",
+   "metadata": {},
+   "source": [
+    "Most textbooks use $x$ and $h$ notation in the convolution. However, in this notebook, we would like to use $Q$ and $T$ as those are the name of arrays in the realm of STUMPY. We would also like to introduce $Q'$, which is simply the reverse of $Q$."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "3dda604a-3a3e-42f4-90c7-7b01ada16337",
+   "metadata": {},
+   "source": [
+    "# Linear Convolution of Two Finite-Length Sequences"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "9ee474bd-1229-4749-80fe-f90c4a135931",
+   "metadata": {},
+   "source": [
+    "Consider two finte-length sequences $Q'$ (of length $m$) and $T$ (of length $n$). Their linear convolution, $C$, can be computed as follows:"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "ad4a63b5-2d91-4875-b4c2-bc7ce11b6af9",
+   "metadata": {},
+   "source": [
+    "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]} $$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "8028ee30-3244-4c79-86ca-22e2b8876694",
+   "metadata": {},
+   "source": [
+    "where $T[.]$ and $Q'[.]$ are zeros for any index that is outside of their range. So, $T[i]$ is zero if $i$ is outside of the range $0 \\le i \\le n-1$. Hence, the equation above is equivalent to:\n",
+    "\n",
+    "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]} $$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "593a8276-49cb-43d2-9300-8d29e67d3104",
+   "metadata": {},
+   "source": [
+    "Note that $Q'[idx-i]$ is not zero when $0 \\le idx-i \\le m-1$, or equivalently $i \\le idx \\le i+m-1$. \n",
+    "\n",
+    "Therefore:\n",
+    "\n",
+    "* $idx \\ge i, \\quad i \\ge 0 \\implies idx \\ge 0$\n",
+    "* $idx \\le i+m-1, \\quad i \\le n-1 \\implies idx \\le n+m-2$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e638d3d7-cf2a-4355-aff3-a9192400defb",
+   "metadata": {},
+   "source": [
+    "This shows that the linear convolution $C$ has values at indices $\\set{0, 1, ..., n + m - 2}$, and it is 0 otherwise. Therefore, the length of output in linear convolution is $n+m-1$. \n",
+    "\n",
+    "Furthermore, let's check out the values for different indices:\n",
+    "* $idx=0 \\implies C[0]=T[0]Q'[0]$\n",
+    "* $idx=1 \\implies C[0]=T[0]Q'[1] + T[1]Q'[0]$\n",
+    "* ...\n",
+    "* $idx=m-1 \\implies C[m-1]=T[0]Q'[m-1] + T[1]Q'[m-2] + ... + T[m-1]Q'[0] = T_{0}.Q$\n",
+    "* $idx=m \\implies C[m]=T[1]Q'[m-1] + T[2]Q'[m-2] + ... + T[m]Q'[0] = T_{1}.Q$\n",
+    "* ...\n",
+    "* $idx=n-1 \\implies C[n-1]=T[n-m]Q'[(n-1)-(n-m)] + T[n-m+1]Q'[(n-1)-(n-m+1)] + ... + T[n-1]Q'[0] = T_{n-m+1}.Q$\n",
+    "* ...\n",
+    "* $idx=n+m-2 \\implies C[n+m-2]=T[n-1]Q'[(n+m-2)-(n-1)] = T[n-1]Q'[m-1]$\n",
+    "\n",
+    "As observed, when $m-1 \\le idx \\le n-1$, $C[idx]$ becomes the dot product between a subsequnce of $T$ and $Q$, which is the reverse of $Q'$. Therefore, the `range(m-1,n)` gives the sliding dot product between $Q$ and $T$."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "48c7f9d8-1b7c-4abc-886b-27a015b78ed6",
+   "metadata": {},
+   "source": [
+    "**How can we leverage this to speed up the computation of sliding dot product?**"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "52bcf566-3ac0-483d-81a3-57fa3b5766c3",
+   "metadata": {},
+   "source": [
+    "# Option I: Convert to Circular Convolution and use FFT-IFFT"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "2c4176e2-60f3-410f-ad2d-6c1d29cca6cf",
+   "metadata": {},
+   "source": [
+    "Circular convolution is defined between two sequences that are both periodic and their period are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$, and it can be computed as follows:\n",
+    "\n",
+    "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{Q'}[i] \\times \\tilde{T}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n",
+    "\n",
+    "where, $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$. $C_{N}$ represents one period of N-Circular convolution. This can also be computed via FFT, i.e. $$C_{N} = IFFT( FFT(\\tilde{Q'}_{N}) \\times FFT(\\tilde{T}_{N}) ) $$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e10c5304-a1a1-4cbd-9b26-266e229bdc57",
+   "metadata": {},
+   "source": [
+    "If there is a way to compute the linear convolution via circular convolution, then we can take advantage of FFT by using eq (4). The good news is that there is a way! The linear convolution between $Q'$ and $T$ can be obtained by performing N-circular convolution between $\\tilde{Q'}$ and $\\tilde{T}$, where:\n",
+    "\n",
+    "* $N \\ge n + m - 2$\n",
+    "* $\\tilde{Q'}_{N}$ is $Q'$ but zero-padded with $N-m$ zeros\n",
+    "* $\\tilde{T}_{N}$ is $T$ but zero-padded with $N-n$ zeros"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "id": "b6529137-fcdc-4d65-aa93-da3746739f94",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# WIP"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "97fcc328-35f9-4b66-9598-1b5d8f1f97ac",
+   "metadata": {},
+   "source": [
+    "# Option II: Overlap-add method\n",
+    "This is a divide and conquer algorithm. It applies divide on \"linear convolution\" and conquer each via circular convolution (See: Option I)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "1e45ccf7-1564-4c41-b32e-4de9b1b201a3",
+   "metadata": {},
+   "source": [
+    "## Linearity in Linear Convolution"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "a084578e-cb20-437c-a70c-3e33f8da8ac3",
+   "metadata": {},
+   "source": [
+    "Suppose the array $T$ can be written as $T_{1} + T_{2}$. In other words: $T[i] = T_{1}[i] + T_{2}[i]$, then:\n",
+    "\n",
+    "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]} $$\n",
+    "\n",
+    "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{(T_{1}[i] \\times Q'[idx-i] + T_{2}[i] \\times Q'[idx-i])} $$\n",
+    "\n",
+    "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T_{1}[i] \\times Q'[idx-i]} + \\sum_{i=-\\infty}^{\\infty}{T_{2}[i] \\times Q'[idx-i]}$$\n",
+    "\n",
+    "$$ C[idx] = C_{1}[idx] + C_{2}[idx] $$"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "1f5be13b-82a7-46d9-9cb1-229c878e67f8",
+   "metadata": {},
+   "source": [
+    "This shows that a linear convolution has the property of \"linearity\". Now, we use this property to show that we can compute the linear convolution between $Q'$ and long $T$ by breaking $T$ into smaller parts."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "fecee5fc-3c10-4c2e-8c66-bce17c3650bd",
+   "metadata": {},
+   "source": []
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "acf5516d-e52d-457c-a8c5-f4d1b7ae8f1c",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "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.14.4"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}