"## Week 4 - Image deblurring using projections\n",
"[LTS2](https://lts2.epfl.ch)\n",
"\n",
"### Objectives\n",
"In this week's exercises you will need to use the properties of the four fundamental subspaces, as well as the projection operator on those subspaces. You will apply those to deblur an image.\n",
"Compute a basis for each of its four fundamental subspaces."
]
},
{
"cell_type": "markdown",
"id": "c5974723",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"source": [
"_Type your answer here, replacing this text._"
]
},
{
"cell_type": "markdown",
"id": "5504f856-3edd-436f-b14f-4e660a7faf7c",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- END QUESTION -->\n",
"\n",
"<!-- BEGIN QUESTION -->\n",
"\n",
"### 2.\n",
"Let $A \\in \\mathbb{R}^{m \\times n}$ and suppose it is onto. Prove that $AA^T$ is nonsingular."
]
},
{
"cell_type": "markdown",
"id": "658270cd",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"source": [
"_Type your answer here, replacing this text._"
]
},
{
"cell_type": "markdown",
"id": "44bf7882-513e-4c29-9619-96bd357670b0",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- END QUESTION -->\n",
"\n",
"## II. Image deblurring \n",
"### Introduction\n",
"Since we will need to visualize images, just a brief reminder on how to use matplotlib for this task"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "a05d6e95-711e-41ce-bf8a-5b983d8999fd",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import skimage\n",
"import skimage.io\n",
"import os"
]
},
{
"cell_type": "markdown",
"id": "f387c338-7c82-46cf-9950-def341903021",
"metadata": {},
"source": [
"We will use the `imread` function from the `scikit-image` package (pre-installed on noto, make sure to install it if you are using a local environment)"
"Some basic information regarding the contents of `camera`:"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "03e539bd-2d68-4240-9d72-ee4ef02b9081",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"camera.shape"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "3b723f81-ebf9-45e5-83a9-6e3d7048586b",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"camera.dtype"
]
},
{
"cell_type": "markdown",
"id": "f2c5f73f-bfee-47e7-b68d-d364f6478663",
"metadata": {},
"source": [
"The type `uint8` means each value in `camera` is represented using a single byte, and is within $[0, 255]$.\n",
"\n",
"Displaying the image is fairly straightforward :"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "04ce2927-48e4-41d6-961f-f1a5f0791cc3",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"plt.imshow(camera, cmap='gray')"
]
},
{
"cell_type": "markdown",
"id": "18174a60-7b66-4ccf-a5f4-43ec56299c0e",
"metadata": {},
"source": [
"The image here is widely used in many image processing scientific publications.\n",
"\n",
"The `cmap` parameter specifies the colormap used to display the image. Matplotlib default colormap not being grayscale, passing `cmap='gray'` is needed. You can find [several predefined colormaps](https://matplotlib.org/stable/tutorials/colors/colormaps.html) in Matplotlib. \n",
"\n",
"You can also check the images distributed with `scikit-image` on their [github repository](https://github.com/scikit-image/scikit-image/tree/main/skimage/data). Feel free to use another image, or even one of your own."
]
},
{
"cell_type": "markdown",
"id": "30dbc6b0-be1f-43b6-a575-ffdedac904b8",
"metadata": {},
"source": [
"## 1. Blurring operator\n",
"\n",
"In order to efficiently deblur an image, let us design our own blurring operator. Since it will be entirely known, it will be possible to compute its inverse.\n",
"\n",
"Let us consider an input vector $x=\\begin{pmatrix}x_0 \\\\ x_1 \\\\ ... \\\\ x_{N-1}\\end{pmatrix} \\in \\mathbb{R}^N$. \n",
"\n",
"Our blurring operator $B$ will do two things:\n",
"- replace each value $x_k$ by the average of $2p+1$ values surrounding $x_k$, i.e. between $x_{k-p}$ and $x_{k+p}$.\n",
"- compute this average not for every $x_k$ but only every $s$ values, with $s$ being a subsampling factor, $s>0$.\n",
"\n",
"Formally, if we denote by $y$ the output of the blurring operator on $x$, this means that\n",
"- $y\\in\\mathbb{R}^M$, with $M=\\frac{N}{s}$ (NOTE: we will only consider the case where $N$ is a mutiple of $s$)\n",
"As you may have noticed, our summation indices in the definition of $y_m$ can be negative in some cases and bigger than $N-1$ in others. In those cases the missing value of $x_k$ is replaced by 0. This way of handling borders is often referred to as **zero-padding**, as it can be achieved by creating an extended version of $x$ with a number of leading and trailing zeros. There are other ways of handling borders, e.g. by using $x_0$ for negative indices and $x_{N-1}$ for indices greater than $N-1$, or even by mirroring the input vector."
]
},
{
"cell_type": "markdown",
"id": "63fecde4-7ee7-4692-b21e-dfc33e4a4de7",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"1. Write a function that creates the matrix $B$ that performs the operations described above.\n",
"\n",
"Remember that in python the `//` operator performs the integer division.\n",
"\n",
"**!!Warning!!** the automated tests are here to help you, they might however be incomplete and do not fully guarantee that your implementation is correct even if they all pass !"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "e40d4cba-1187-448d-96fb-d0033c129aaf",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"outputs": [],
"source": [
"def blur_matrix(N, p, s):\n",
" \"\"\"\n",
" Computes the blurring matrix \n",
"\n",
" Parameters\n",
" ----------\n",
" N : length of the input signal\n",
" p : half-length of the blurring filter\n",
" s: subsampling factor\n",
"\n",
" Returns\n",
" -------\n",
" The blurring matrix \n",
" \"\"\"\n",
" ..."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "73ed3a03",
"metadata": {
"deletable": false,
"editable": false
},
"outputs": [],
"source": [
"grader.check(\"q2p1\")"
]
},
{
"cell_type": "markdown",
"id": "af393d26-6c86-4ff0-b8fb-2f78fae94c57",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- BEGIN QUESTION -->\n",
"\n",
"2. What is the rank of $B$ ? \n",
"\n",
"Hint: \n",
"- row-rank might be easier to consider\n",
"- you might also need to add constraints on $s$ and $p$"
]
},
{
"cell_type": "markdown",
"id": "d4849ba8",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"source": [
"_Type your answer here, replacing this text._"
]
},
{
"cell_type": "markdown",
"id": "67659a9b-5206-40da-a8ab-0f190884d901",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- END QUESTION -->\n",
"\n",
"## 2. Blurring the image\n",
"\n",
"We now have our blurring matrix ready, let us use it to blur the `camera` image.\n",
"As $B$ is designed to operate on vectors, we will need two steps to blur the image\n",
"- apply $B$ on each column of the input image\n",
"- apply $B$ on each row of the column-blurred image computed above"
]
},
{
"cell_type": "markdown",
"id": "596b49e7-438a-42c6-b5e6-5ea1a7f338f2",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"3. Implement a function that blurs an image using the matrix defined in the previous question"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "0a8714b7-1dfa-4aa9-b09a-a5a9caa317bc",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"outputs": [],
"source": [
"def blur_image(img, B):\n",
" \"\"\"\n",
" Blurs the input image using the blur matrix\n",
"\n",
" Parameters\n",
" ----------\n",
" img : input image\n",
" B : blur matrix\n",
"\n",
" Returns\n",
" -------\n",
" The blurred image \n",
" \"\"\"\n",
" ..."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "6498302c",
"metadata": {
"deletable": false,
"editable": false
},
"outputs": [],
"source": [
"grader.check(\"q2p3\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "c167c207-66da-41a7-b536-5fc4577420e3",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"# build the blur matrix\n",
"B = blur_matrix(512, 3, 2)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "42231cfd-b223-4f90-97c5-dec703aee5af",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"# Blur on rows\n",
"camera_blurred = blur_image(camera, B)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "78aa4a09-ff28-4e21-83f9-2b14c6d6ad83",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"plt.imshow(camera_blurred, cmap='gray') # check the result"
"Now everything is setup, we can proceed with the image restoration. "
]
},
{
"cell_type": "markdown",
"id": "c962b694-45b2-4c2b-b7f5-ee69ffd8e342",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- BEGIN QUESTION -->\n",
"\n",
"4. Let us denote the result of a blur operation by $y=Bx$. Using the properties of the four fundamental subspaces, write $x^*$ as a sum, making use of projections operators those subspaces. Hint: you might need results from question 2.\n"
]
},
{
"cell_type": "markdown",
"id": "c35a1228",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"source": [
"_Type your answer here, replacing this text._"
]
},
{
"cell_type": "markdown",
"id": "8f024ccc-04b6-4437-aa6a-daba627ed4d4",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- END QUESTION -->\n",
"\n",
"<!-- BEGIN QUESTION -->\n",
"\n",
"5. Using this result, find a possible solution for the complete blurring operation (on both rows and columns). "
]
},
{
"cell_type": "markdown",
"id": "ab1c7611",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"source": [
"_Type your answer here, replacing this text._"
]
},
{
"cell_type": "markdown",
"id": "cf20cfff-0c77-4367-8b4e-92f41c41acad",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- END QUESTION -->\n",
"\n",
"6. Implement a function that computes the right inverse of the blur matrix, and uses it preform image restoration. Do not forget to consider row and columns blurring.\n",
"\n",
"You can use `np.linalg.inv` to perform the inverse of a matrix."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "32426413-a8ea-46ac-ae3e-968f480ff150",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"outputs": [],
"source": [
"def inverse_blur(blurred_img, B):\n",
" ..."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "f82aa4d6",
"metadata": {
"deletable": false,
"editable": false
},
"outputs": [],
"source": [
"grader.check(\"q2p6\")"
]
},
{
"cell_type": "markdown",
"id": "812c7443-3492-46d2-ad99-b2c83826163b",
"metadata": {},
"source": [
"Using this function, compute the reconstructed image from the blurred one"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "f270a921-f7c3-40fd-a969-d37010ad2954",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"camera_unblur = inverse_blur(camera_blurred, B)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "75296a85-053c-407c-915f-b562ba9c8acd",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"plt.imshow(camera_unblur) # check the result"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "0dfd3731-b40a-4425-b78d-e0b65590c9f3",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"# Let us compare the original image, the restored image and an upscaled version of the blurred image\n",
"6. Evaluate the performance of the restoration method. What are its limitations for a real-world usage ? Check that the computed solution is valid by blurring the restored image."
]
},
{
"cell_type": "markdown",
"id": "2d3adde9",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"source": [
"_Type your answer here, replacing this text._"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "3658f662-3832-4098-84f4-5d618b95fd02",
"metadata": {
"tags": [
"otter_answer_cell"
]
},
"outputs": [],
"source": [
"..."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "54f26333",
"metadata": {
"deletable": false,
"editable": false
},
"outputs": [],
"source": [
"grader.check(\"q2p6\")"
]
},
{
"cell_type": "markdown",
"id": "650659df",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"<!-- END QUESTION -->\n",
"\n"
]
},
{
"cell_type": "markdown",
"id": "0bd3e058",
"metadata": {
"deletable": false,
"editable": false
},
"source": [
"## Submission\n",
"\n",
"Make sure you have run all cells in your notebook in order before running the cell below, so that all images/graphs appear in the output. The cell below will generate a zip file for you to submit. **Please save before exporting!**"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d10563fe",
"metadata": {
"deletable": false,
"editable": false
},
"outputs": [],
"source": [
"# Save your notebook first, then run this cell to export your submission.\n",
In this week's exercises you will need to use the properties of the four fundamental subspaces, as well as the projection operator on those subspaces. You will apply those to deblur an image.
We will use the `imread` function from the `scikit-image` package (pre-installed on noto, make sure to install it if you are using a local environment)
The image here is widely used in many image processing scientific publications.
The `cmap` parameter specifies the colormap used to display the image. Matplotlib default colormap not being grayscale, passing `cmap='gray'` is needed. You can find [several predefined colormaps](https://matplotlib.org/stable/tutorials/colors/colormaps.html) in Matplotlib.
You can also check the images distributed with `scikit-image` on their [github repository](https://github.com/scikit-image/scikit-image/tree/main/skimage/data). Feel free to use another image, or even one of your own.
In order to efficiently deblur an image, let us design our own blurring operator. Since it will be entirely known, it will be possible to compute its inverse.
Let us consider an input vector $x=\begin{pmatrix}x_0 \\ x_1 \\ ... \\ x_{N-1}\end{pmatrix} \in \mathbb{R}^N$.
Our blurring operator $B$ will do two things:
- replace each value $x_k$ by the average of $2p+1$ values surrounding $x_k$, i.e. between $x_{k-p}$ and $x_{k+p}$.
- compute this average not for every $x_k$ but only every $s$ values, with $s$ being a subsampling factor, $s>0$.
Formally, if we denote by $y$ the output of the blurring operator on $x$, this means that
- $y\in\mathbb{R}^M$, with $M=\frac{N}{s}$ (NOTE: we will only consider the case where $N$ is a mutiple of $s$)
- $y_m = \frac{1}{2p+1}\sum_{k=sm-p}^{sm+p}x_k$
As you may have noticed, our summation indices in the definition of $y_m$ can be negative in some cases and bigger than $N-1$ in others. In those cases the missing value of $x_k$ is replaced by 0. This way of handling borders is often referred to as **zero-padding**, as it can be achieved by creating an extended version of $x$ with a number of leading and trailing zeros. There are other ways of handling borders, e.g. by using $x_0$ for negative indices and $x_{N-1}$ for indices greater than $N-1$, or even by mirroring the input vector.
1. Write a function that creates the matrix $B$ that performs the operations described above.
Remember that in python the `//` operator performs the integer division.
**!!Warning!!** the automated tests are here to help you, they might however be incomplete and do not fully guarantee that your implementation is correct even if they all pass !
4. Let us denote the result of a blur operation by $y=Bx$. Using the properties of the four fundamental subspaces, write $x^*$ as a sum, making use of projections operators those subspaces. Hint: you might need results from question 2.
6. Implement a function that computes the right inverse of the blur matrix, and uses it preform image restoration. Do not forget to consider row and columns blurring.
You can use `np.linalg.inv` to perform the inverse of a matrix.
6. Evaluate the performance of the restoration method. What are its limitations for a real-world usage ? Check that the computed solution is valid by blurring the restored image.
Make sure you have run all cells in your notebook in order before running the cell below, so that all images/graphs appear in the output. The cell below will generate a zip file for you to submit. **Please save before exporting!**
%% Cell type:code id:d10563fe tags:
``` python
# Save your notebook first, then run this cell to export your submission.
