{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Laplacian Operator\n",
    "===\n",
    "\n",
    "Here, we show you what the laplacian operator and how it works on the example images"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Filtering in an image means to perform some kind of processing on a pixel value $I(x, y)$\n",
    "using its neighboring pixel values as follows:\n",
    "$$\n",
    "I'(x, y) = \\sum_{i=-1}^{1}\\sum_{j=-1}^{1}K(i, j)I(x + i, y + j),\n",
    "$$\n",
    "where, $I'(x, y)$: performed pixel value, $K$: kernel matrix.\n",
    "\n",
    "Laplacian Operator is a derivative operator which is used to find edges in an image and represented\n",
    "as following kernels:\n",
    "$$\n",
    "\\begin{align}\n",
    "K_4 &=\n",
    "\\begin{pmatrix}\n",
    "0 & 1 & 0\\\\\n",
    "1 & -4 & 1\\\\\n",
    "0 & 1 & 0\n",
    "\\end{pmatrix},\\\\\n",
    "K_8 &=\n",
    "\\begin{pmatrix}\n",
    "1 & 1 & 1\\\\\n",
    "1 & -8 & 1\\\\\n",
    "1 & 1 & 1\n",
    "\\end{pmatrix},\n",
    "\\end{align}\n",
    "$$\n",
    "where, $K_4$: a kernel that considers the contribution of 4 nearest neighbors\n",
    "(top, bottom, left, right) to the pixel of interest, $K_8$: a kernel that considers 8 nearest\n",
    "neighbors (top, bottom, left, right, diagonal) to the pixel of interest.\n",
    "\n",
    "To perform this operator to a image converted 1-D vector array, we generate the laplacian operator\n",
    "matrix.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "from matplotlib.cbook import get_sample_data\n",
    "from matplotlib.colors import CenteredNorm\n",
    "from PIL import Image\n",
    "\n",
    "from cherab.phix.tools import laplacian_matrix"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Visualize simple laplacian matrix\n",
    "---\n",
    "\n",
    "Try to create the simple laplacian matrix (50, 50).\n",
    "Firstly, create the mapping array denoting 2-D image shape and the element of which denotes a index."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# mapping array\n",
    "mapping_array = np.arange(0, 50, dtype=np.int32).reshape(10, 5)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot laplacian matrix as a sparse matrix and compare $K_4$ and $K_8$ laplacian kernel"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axes = plt.subplots(1, 2, dpi=150, tight_layout=True)\n",
    "for ax, neighbors in zip(axes, [4, 8]):\n",
    "\n",
    "    # calculate laplacian matrix\n",
    "    laplacian_mat = laplacian_matrix(mapping_array, dir=neighbors)\n",
    "\n",
    "    # plot sparse matrix\n",
    "    ax.spy(laplacian_mat, markersize=2)\n",
    "    ax.set_title(f\"laplacian matrix K$_{neighbors}$\", pad=25)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "show laplacian matrix $K_8$ in (10, 10) size as a numpy array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "laplacian_mat.toarray()[0:10, 0:10]"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Apply the laplacian matrix to a sample image\n",
    "---\n",
    "Next, let us to apply a laplacian matrix to pixels of a sample image."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Load sample image data from the matplotlib library."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with get_sample_data(\"grace_hopper.jpg\") as file:\n",
    "    arr_image = plt.imread(file)\n",
    "\n",
    "# resize the image deu to the large size.\n",
    "with Image.fromarray(arr_image, mode=\"RGB\") as im:\n",
    "    (width, height) = (im.width // 4, im.height // 4)\n",
    "    arr_image = np.array(im.resize((width, height)))\n",
    "\n",
    "# convert RGB image to monotonic one\n",
    "arr_image = arr_image.mean(axis=2)\n",
    "\n",
    "# show image\n",
    "print(f\"image array shape: {arr_image.shape}\")\n",
    "fig, ax = plt.subplots(dpi=150)\n",
    "ax.imshow(arr_image, cmap=\"gray\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create mapping array\n",
    "image_map = np.arange(0, arr_image.size, dtype=np.int32).reshape(arr_image.shape)\n",
    "\n",
    "# create laplacian matrix with K8\n",
    "neighbors = 8\n",
    "laplacian_mat = laplacian_matrix(image_map, dir=neighbors)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The laplacian filtered image is calculated by multiplying the image vector by the laplacian matrix."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "filterd_image = np.reshape(laplacian_mat @ arr_image.ravel(), arr_image.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "nbsphinx-thumbnail"
    ]
   },
   "outputs": [],
   "source": [
    "# show image\n",
    "fig, ax = plt.subplots(dpi=150)\n",
    "mappable = ax.imshow(filterd_image, cmap=\"seismic\", norm=CenteredNorm())\n",
    "cbar = plt.colorbar(mappable)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here, the edge of the image is emphasized clearly.\n",
    "In the tomography technique, we make use of this operator to smooth reconstructed images."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "cherab-phix-dev",
   "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.9.16"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "2725905a4c02db19e04df9b8fdbbe5ec65a73ea52bebaf9474aa1cc98819834c"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}