commit e340d6b99e6e79c8191e80a771280d0504190b23 Author: Johannes Kattnig Date: Mon Dec 2 20:08:23 2024 +0100 Initial Commit diff --git a/Kattnig_report_A0.pdf b/Kattnig_report_A0.pdf new file mode 100755 index 0000000..e983104 Binary files /dev/null and b/Kattnig_report_A0.pdf differ diff --git a/Tutorial/0-numpy.ipynb b/Tutorial/0-numpy.ipynb new file mode 100755 index 0000000..a917a5d --- /dev/null +++ b/Tutorial/0-numpy.ipynb @@ -0,0 +1,490 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# NumPy - Numeric Python\n", + "NumPy is the core library for scientifc computing Python. The module provides high-performance implementations for dealing with n-dimensional arrays and tools for working with these arrays." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## How to work with NumPy\n", + "The fundamental data structure in NumPy is called _array_. An array is a n-dimensional container. It can be used for \n", + "\n", + "* vectors\n", + "* matrices\n", + "* n-dimensional data\n", + "\n", + "Internally a NumPy array is stored row major (C order) by default. Note that this is different to Matlab where arrays are stored column major." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Import NumPy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib notebook\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Creating NumPy arrays" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# integer vector, 1D array\n", + "v = np.array([1, 2, 3])\n", + "v" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# matrix, 2D array\n", + "M = np.array([[1, 2], \n", + " [3, 4]], dtype=float)\n", + "M" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# volume, 3D array\n", + "V = np.array([[[1, 2, 3], [4, 5, 6]], [[7, 8 , 9], [10, 11, 12]]])\n", + "V" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Initial placeholders" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create array of zeros with 3 rows and 4 columns\n", + "np.zeros((3, 4))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create array of ones \n", + "np.ones((2, 3, 4), dtype=np.int16)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create array of evenly spaced values\n", + "# args: start, stop, step: inclusive start, exclusive stop\n", + "np.arange(3, 11, 2)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create array with evenly spaced values using number of samples\n", + "# start, stop, num_samples\n", + "x = np.linspace(-1, 4, 20)\n", + "\n", + "_, ax = plt.subplots(figsize=(5, 3.5))\n", + "ax.plot(x, np.sin(x))\n", + "\n", + "_, ax = plt.subplots(figsize=(5, 3.5))\n", + "ax.plot(np.sin(x))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create a constant array \n", + "np.full((2, 3), 7)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create a n x n identity matrix\n", + "np.eye(3)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create an array with random values\n", + "np.random.random((2, 3))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# create an empty array (attention: memory is uninitialized!)\n", + "np.empty((3, 2))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "base = np.array([1, 2, 3])\n", + "\n", + "a = []\n", + "for i in range(10):\n", + " a.append(base.copy())\n", + "a = np.array(a)\n", + " \n", + "b = np.stack([base.copy()] * 10)\n", + "c = np.repeat(base[None], 10, axis=0)\n", + "d = np.tile(base, (10, 1))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "print(a)\n", + "assert np.allclose(a, b) and np.allclose(b, c) and np.allclose(c, d)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Inspecting array" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# array dimensions\n", + "V.shape" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# number of dimensions\n", + "V.ndim" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# number of elements\n", + "V.size" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# data type of array elements\n", + "V.dtype" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Array operations" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# sum\n", + "V.sum()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# min\n", + "V.min()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# max\n", + "V.max()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# mean\n", + "V.mean()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# standard deviation\n", + "V.std()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Arithmetic Operations" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "a = np.arange(3)\n", + "a" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# addition\n", + "a + a" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "a + 1" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# exp\n", + "np.exp(a)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# * is elementwise multiplication\n", + "# @ is __matmul__\n", + "# a is one-dimensional\n", + "for inner in [(a * a).sum(), a @ a, a.dot(a), a.T.dot(a)]:\n", + " print(inner)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "for outer in [np.outer(a, a), a[None] * a[:, None]]:\n", + " print(outer)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Slicing, Indexing\n", + "In python indexing starts at 0!" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# get first element\n", + "a\n", + "a[0]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "M" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# get 1st row\n", + "M[0] # = M[0,:]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# get last column\n", + "M[:,-1]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "v = np.arange(10)\n", + "v" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# get every second element\n", + "# [start:stop:step]\n", + "v[1::2]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# boolean indexing\n", + "v[v > 3]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "v > 3" + ] + } + ], + "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.10.7" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/Tutorial/1-images.ipynb b/Tutorial/1-images.ipynb new file mode 100755 index 0000000..0d258c0 --- /dev/null +++ b/Tutorial/1-images.ipynb @@ -0,0 +1,339 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Working with images" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib notebook\n", + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "import imageio.v2 as imageio" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Image loading and display" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "image = imageio.imread('taj_mahal_gray.jpg')\n", + "\n", + "start, step = 300, 30\n", + "print(f'The image \\n{image}\\nhas {image.shape[0]} rows and {image.shape[1]} columns.')\n", + "print(f'The {step} pixels in the zero-th row, starting at column {start} contain the values \\n{image[0, start:start + step]}.')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "scrolled": false + }, + "outputs": [], + "source": [ + "_, ax = plt.subplots(1, 2, figsize=(9, 3))\n", + "ax[0].imshow(image)\n", + "ax[1].imshow(image, cmap='gray')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Deep and shallow copy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "image = imageio.imread('taj_mahal_gray.jpg')\n", + "\n", + "# Assign another name to the same underlying storage (\"Reference\")\n", + "# copy = image\n", + "\n", + "# Copy the underlying data to a new storage location (\"Copy\")\n", + "copy = image.copy()\n", + "\n", + "copy[200:350, 300:450] = 0\n", + "\n", + "fig, ax = plt.subplots(1, 2, figsize=(9, 3))\n", + "ax[0].imshow(image)\n", + "ax[1].imshow(copy, cmap='gray')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Image statistics & manipulations" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "im_integral = imageio.imread('taj_mahal_gray.jpg')\n", + "im_float = im_integral / 255. # normalize grayvalues to [0,1]\n", + "\n", + "fig, axs = plt.subplots(1, 2, figsize=(9, 3))\n", + "for ax, im in zip(axs, [im_integral, im_float]):\n", + " artist = ax.imshow(im, cmap='gray')\n", + " ax.set_title(f'min={im.min()}, max={im.max()}')\n", + " fig.colorbar(artist, ax=ax)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "scrolled": false + }, + "outputs": [], + "source": [ + "im = imageio.imread('taj_mahal_small.jpg')\n", + "print(im.shape)\n", + "_, ax = plt.subplots()\n", + "ax.imshow(im)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Conversion to grayscale image" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "im = im / 255.0\n", + "\n", + "coeffs = np.array([0.299, 0.587, 0.114])\n", + "gray_explicit = im[:,:,0] * 0.299 + im[:,:,1] * 0.587 + im[:,:,2] * 0.114\n", + "gray_broadcasting = (im * coeffs).sum(axis=2)\n", + "assert np.allclose(gray_broadcasting, gray_explicit)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "gray_mean = im.mean(-1)\n", + "\n", + "_, ax = plt.subplots(1, 2, figsize=(9, 3))\n", + "ax[0].imshow(gray_broadcasting, cmap='gray')\n", + "ax[1].imshow(gray_mean, cmap='gray')" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "collapsed": true + }, + "source": [ + "## Image filtering" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import scipy.ndimage as sn\n", + "from skimage.transform import rescale" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# load image\n", + "im = imageio.imread('taj_mahal_gray.jpg') / 255.0\n", + "im_aa = rescale(im, 0.2, anti_aliasing=True)\n", + "im_noaa = rescale(im, 0.2, anti_aliasing=False)\n", + "\n", + "_, ax = plt.subplots(1, 2, figsize=(9,3))\n", + "ax[0].imshow(im_aa, cmap='gray')\n", + "ax[0].set_title('AA')\n", + "ax[1].imshow(im_noaa, cmap='gray')\n", + "ax[1].set_title('No AA')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Gaussian filter" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "sigma = 2.3\n", + "kernel_size = 15\n", + "\n", + "space = np.arange(-kernel_size // 2 + 1, kernel_size // 2 + 1)\n", + "x, y = np.meshgrid(space, space)\n", + "\n", + "kernel = np.exp(-(x ** 2 + y ** 2) / (2 * sigma ** 2))\n", + "kernel = kernel / np.sum(kernel)\n", + "_, ax = plt.subplots(subplot_kw={'projection': '3d'})\n", + "ax.plot_surface(x, y, kernel, cmap='viridis')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Apply gaussian filter to image" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# convolution with the gaussian filter\n", + "filtered = sn.convolve(im, kernel)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "_, ax = plt.subplots(1, 2, figsize=(9,3));\n", + "ax[0].imshow(im, cmap='gray')\n", + "ax[0].set_title('Original')\n", + "ax[1].imshow(filtered, cmap='gray')\n", + "ax[1].set_title('2D filter')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Edge filter" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "dx = np.array([\n", + " [-1, 0, 1],\n", + " [-1, 0, 1],\n", + " [-1, 0, 1]\n", + "])\n", + "dy = dx.T\n", + "\n", + "_, ax = plt.subplots(1, 2, figsize=(5,2))\n", + "ax[0].imshow(dx, cmap='gray')\n", + "ax[1].imshow(dy, cmap='gray')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "nabla_x = sn.convolve(im, dx)\n", + "nabla_y = sn.convolve(im, dy)\n", + "grad_magnitude = np.sqrt(nabla_x ** 2 + nabla_y ** 2)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fi, ax = plt.subplots(1, 3, figsize=(9,3))\n", + "ax[0].imshow(nabla_x, cmap='gray')\n", + "ax[1].imshow(nabla_y, cmap='gray')\n", + "ax[2].imshow(grad_magnitude, cmap='gray')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "nabla_x = sn.convolve(filtered, dx)\n", + "nabla_y = sn.convolve(filtered, dy)\n", + "grad_magnitude = np.sqrt(nabla_x ** 2 + nabla_x ** 2)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fi, ax = plt.subplots(1, 3, figsize=(9,3))\n", + "ax[0].imshow(nabla_x, cmap='gray')\n", + "ax[1].imshow(nabla_y, cmap='gray')\n", + "ax[2].imshow(grad_magnitude, cmap='gray')" + ] + } + ], + "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.10.7" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/Tutorial/2-implementation.ipynb b/Tutorial/2-implementation.ipynb new file mode 100755 index 0000000..f0d0c11 --- /dev/null +++ b/Tutorial/2-implementation.ipynb @@ -0,0 +1,361 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Explicit for-loop vs. built-ins" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "def index_sum_for(arr, idx):\n", + " sum_ = 0\n", + " for a, i in zip(arr, idx):\n", + " if i:\n", + " sum_ += a\n", + " return sum_\n", + "\n", + "\n", + "def index_sum_np(arr, idx):\n", + " return arr[idx].sum()" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [], + "source": [ + "num_points = int(5e5)\n", + "array = np.random.randn((num_points))\n", + "indices = np.random.choice([False, True], size=num_points)" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "44.7 ms ± 192 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n" + ] + } + ], + "source": [ + "%timeit sum_for = index_sum_for(array, indices)\n", + "sum_for = index_sum_for(array, indices)" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "2.17 ms ± 2.25 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n" + ] + } + ], + "source": [ + "%timeit sum_np = index_sum_np(array, indices)\n", + "sum_np = index_sum_np(array, indices)" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [], + "source": [ + "assert np.allclose(sum_for, sum_np)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Pure allocation vs. allocation + initialization" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [], + "source": [ + "num_points = int(1e4)" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "6.5 µs ± 33.6 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)\n" + ] + } + ], + "source": [ + "%timeit out = np.empty((num_points, num_points))" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "72.1 ms ± 626 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n" + ] + } + ], + "source": [ + "%timeit out = np.full((num_points, num_points), 1)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Explicit for-loop vs. \"batched\" operations" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "metadata": {}, + "outputs": [], + "source": [ + "rng = np.random.default_rng()\n", + "A = rng.standard_normal((3, 4))\n", + "batched_x = rng.standard_normal((1_000_000, 4))" + ] + }, + { + "cell_type": "code", + "execution_count": 21, + "metadata": {}, + "outputs": [], + "source": [ + "def apply_matrix_for(A, xs):\n", + " b = np.empty((xs.shape[0], A.shape[0]))\n", + " for i, x in enumerate(xs):\n", + " b[i] = A @ x\n", + " return b" + ] + }, + { + "cell_type": "code", + "execution_count": 22, + "metadata": {}, + "outputs": [], + "source": [ + "def apply_matrix_batched(A, xs):\n", + " return xs @ A.T" + ] + }, + { + "cell_type": "code", + "execution_count": 23, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "1.37 s ± 20.9 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n" + ] + } + ], + "source": [ + "%timeit b_for = apply_matrix_for(A, batched_x)\n", + "b_for = apply_matrix_for(A, batched_x)" + ] + }, + { + "cell_type": "code", + "execution_count": 24, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "2.15 ms ± 270 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n" + ] + } + ], + "source": [ + "%timeit b_batched = apply_matrix_batched(A, batched_x)\n", + "b_batched = apply_matrix_batched(A, batched_x)" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "metadata": {}, + "outputs": [], + "source": [ + "assert np.allclose(b_for, b_batched)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Hints for Assignment 0" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([[ 6., 2., 3., 4.],\n", + " [ 5., 11., 7., 8.],\n", + " [ 9., 10., 16., 12.],\n", + " [13., 14., 15., 21.]])" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "np.set_printoptions(precision=1)\n", + "X = np.arange(1, 17, dtype=float).reshape(4, 4) + np.eye(4) * 5\n", + "kernel = np.array([\n", + " [1, 2, 1],\n", + " [2, 4, 2],\n", + " [1, 2, 1],\n", + "]) / 16\n", + "X" + ] + }, + { + "cell_type": "code", + "execution_count": 49, + "metadata": {}, + "outputs": [], + "source": [ + "def variance_naive(X):\n", + " M, N = X.shape\n", + " var_im = np.zeros_like(X)\n", + " for i in range(1, M - 1):\n", + " for j in range(1, N - 1):\n", + " x = X[i - 1:i + 2, j - 1:j + 2]\n", + " mu_x = (x * kernel).sum()\n", + " var_im[i, j] = (kernel * (x - mu_x) ** 2).sum()\n", + "\n", + " return var_im" + ] + }, + { + "cell_type": "code", + "execution_count": 50, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([[0. , 0. , 0. , 0. ],\n", + " [0. , 8.7, 9.4, 0. ],\n", + " [0. , 7.9, 8.7, 0. ],\n", + " [0. , 0. , 0. , 0. ]])" + ] + }, + "execution_count": 50, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "naive = variance_naive(X)\n", + "naive" + ] + }, + { + "cell_type": "code", + "execution_count": 54, + "metadata": {}, + "outputs": [], + "source": [ + "import scipy.signal as ss\n", + "def variance(X):\n", + " # in your code this corresponds to gaussian_filter calls\n", + " mu_x = ss.convolve2d(X, kernel, mode='valid')\n", + " mu_xx = ss.convolve2d(X * X, kernel, mode='valid')\n", + " return np.pad(mu_xx - mu_x ** 2, 1)" + ] + }, + { + "cell_type": "code", + "execution_count": 55, + "metadata": {}, + "outputs": [], + "source": [ + "smart = variance(X)\n", + "assert np.allclose(smart, naive)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "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.10.8" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/Tutorial/taj_mahal_gray.jpg b/Tutorial/taj_mahal_gray.jpg new file mode 100755 index 0000000..3a25146 Binary files /dev/null and b/Tutorial/taj_mahal_gray.jpg differ diff --git a/Tutorial/taj_mahal_small.jpg b/Tutorial/taj_mahal_small.jpg new file mode 100755 index 0000000..7598907 Binary files /dev/null and b/Tutorial/taj_mahal_small.jpg differ diff --git a/ass0.pdf b/ass0.pdf new file mode 100755 index 0000000..72c207a Binary files /dev/null and b/ass0.pdf differ diff --git a/girl.png b/girl.png new file mode 100755 index 0000000..7ba4fc7 Binary files /dev/null and b/girl.png differ diff --git a/girl_cubic.png b/girl_cubic.png new file mode 100755 index 0000000..eea4456 Binary files /dev/null and b/girl_cubic.png differ diff --git a/girl_linear.png b/girl_linear.png new file mode 100755 index 0000000..5a3f2b8 Binary files /dev/null and b/girl_linear.png differ diff --git a/girl_nearest.png b/girl_nearest.png new file mode 100755 index 0000000..8d23806 Binary files /dev/null and b/girl_nearest.png differ diff --git a/interpolation_error.py b/interpolation_error.py new file mode 100755 index 0000000..a93e736 --- /dev/null +++ b/interpolation_error.py @@ -0,0 +1,109 @@ +import numpy as np +import imageio +from skimage.transform import resize +from scipy.ndimage import gaussian_filter + + +def mssim( + x: np.ndarray, + y: np.ndarray, +) -> float: + # Standard choice for the parameters + K1 = 0.01 + K2 = 0.03 + sigma = 1.5 + truncate = 3.5 + m = 1 + C1 = (K1 * m) ** 2 + C2 = (K2 * m) ** 2 + + # radius size of the local window (needed for + # normalizing the standard deviation) + r = int(truncate * sigma + 0.5) + win_size = 2 * r + 1 + # use these arguments for the gaussian filtering + # e.g. filtered = gaussian_filter(x, **filter_args) + filter_args = { + 'sigma': sigma, + 'truncate': truncate + } + + # Implement Eq. (9) from assignment sheet + # S should be an "image" of the SSIM evaluated for a window + # centered around the corresponding pixel in the original input image + S = np.ones_like(x) + + mu_x = gaussian_filter(x, **filter_args) + mu_y = gaussian_filter(y, **filter_args) + + # Compute local variances and covariance + n = win_size ** 2 # number of pixels in the window + # n~ = n/(n-1) + # E[x^2] −E[x]^2). = (gaussian_filter(x * x, **filter_args) - mu_x * mu_x) + sigma_x_sq = (gaussian_filter(x * x, **filter_args) - mu_x * mu_x) * (n / (n - 1)) + sigma_y_sq = (gaussian_filter(y * y, **filter_args) - mu_y * mu_y) * (n / (n - 1)) + sigma_xy = (gaussian_filter(x * y, **filter_args) - mu_x * mu_y) * (n / (n - 1)) + + # Compute SSIM + S = ((2 * mu_x * mu_y + C1) * (2 * sigma_xy + C2)) / \ + ((mu_x ** 2 + mu_y ** 2 + C1) * (sigma_x_sq + sigma_y_sq + C2)) + + # crop to remove boundary artifacts, return MSSIM + pad = (win_size - 1) // 2 + return S[pad:-pad, pad:-pad].mean() + + +def psnr( + x: np.ndarray, + y: np.ndarray, +) -> float: + # Implement Eq. (2) without for loops + mse = np.mean((x - y) ** 2) + m = 1 + psnr_value = 10 * np.log10(m ** 2 / mse) + return psnr_value + + +def psnr_for( + x: np.ndarray, + y: np.ndarray, +) -> float: + # Implement Eq. (2) using for loops + m, n = x.shape + mse = 0 + # Eq. (1) from the assignment sheet for mse + for i in range(m): + for j in range(n): + mse += (x[i, j] - y[i, j]) ** 2 + mse /= (m * n) + m = 1 + psnr_value = 10 * np.log10(m ** 2 / mse) + return psnr_value + + +def interpolation_error(): + x = imageio.imread('./girl.png') / 255. + shape_lower = (x.shape[0] // 2, x.shape[1] // 2) + # downsample image to half the resolution + # and successively upsample to the original resolution + # using no nearest neighbor, linear and cubic interpolation + nearest, linear, cubic = [ + resize(resize( + x, shape_lower, order=order, anti_aliasing=False + ), x.shape, order=order, anti_aliasing=False) + for order in [0, 1, 3] + ] + + for label, rescaled in zip( + ['nearest', 'linear', 'cubic'], + [nearest, linear, cubic] + ): + print(label) + print(mssim(x, rescaled)) + print(psnr(x, rescaled)) + print(psnr_for(x, rescaled)) + imageio.imwrite('girl_' + label + '.png', (rescaled * 255).astype(np.uint8)) + + +if __name__ == '__main__': + interpolation_error() diff --git a/reference_output.txt b/reference_output.txt new file mode 100755 index 0000000..f8ab901 --- /dev/null +++ b/reference_output.txt @@ -0,0 +1,12 @@ +nearest +0.8031736958539067 +27.16729887950422 +27.167298879503985 +linear +0.8890060634234201 +31.092116935553634 +31.092116935553662 +cubic +0.9139147893771699 +31.97622728671044 +31.976227286710447 diff --git a/test_compare_output.py b/test_compare_output.py new file mode 100755 index 0000000..d69b7cc --- /dev/null +++ b/test_compare_output.py @@ -0,0 +1,29 @@ +import unittest +import subprocess + +class TestOutputComparison(unittest.TestCase): + + def test_compare_output(self): + # Define the Python script to run and the expected output file + script_to_run = 'interpolation_error.py' # Replace with your script's name + expected_file = 'reference_output.txt' + + # Run the script and capture the output + result = subprocess.run(['python3', script_to_run], capture_output=True, text=True) + output_lines = result.stdout.splitlines() + + # Open and read the expected output file + with open(expected_file, 'r') as expected: + expected_lines = expected.read().splitlines() + + # Ensure the number of lines is the same + self.assertEqual(len(output_lines), len(expected_lines), + "The number of lines in the output does not match the expected output.") + + # Compare each line + for i, (output_line, expected_line) in enumerate(zip(output_lines, expected_lines)): + self.assertEqual(output_line.strip(), expected_line.strip(), + f"Line {i+1} does not match:\nOutput: {output_line}\nExpected: {expected_line}") + +if __name__ == '__main__': + unittest.main()