03_Histograms.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "QuUAgQFQlsJU"
      },
      "source": [
        "# Tutorial 03 - Histograms and Thresholding  \n",
        "## Dr. David C. Schedl\n",
        "\n",
        "Note: this tutorial is geared towards students **experienced in general programming** and aims to introduce you to OpenCV.\n",
        "\n",
        "Adapted from: \n",
        "* http://6.869.csail.mit.edu/fa19/schedule.html (written by Julie Ganeshan; @MIT)\n",
        "\n",
        "Useful links:\n",
        "* OpenCV Tutorials: https://docs.opencv.org/master/d9/df8/tutorial_root.html\n",
        "* Image Processing in Pyhton: https://github.com/xn2333/OpenCV/blob/master/Seminar_Image_Processing_in_Python.ipynb\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "weMMc2WZlsJc"
      },
      "source": [
        "# Initilization"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "zkUzZ6lWK_S5"
      },
      "source": [
        "Let's import useful libraries, first."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Gpvh7PPVlsJg"
      },
      "outputs": [],
      "source": [
        "import os\n",
        "import cv2 # openCV\n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "import plotly.graph_objects as go\n",
        "import plotly.express as px\n",
        "import plotly.io as pio"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "DnGyZdM2LH71"
      },
      "source": [
        "Let's download some images to work with. \n",
        "We use the Unix/Windows command `curl`. \n",
        "Images are in the local filesystem after downloading.\n",
        "\n",
        "Image sources:\n",
        "\n",
        "* A picture of a [cat from wikimedia](https://upload.wikimedia.org/wikipedia/commons/thumb/3/3a/Cat03.jpg/481px-Cat03.jpg).\n",
        "* A picture of [Van Gogh from wikimedia](https://upload.wikimedia.org/wikipedia/commons/thumb/b/b2/Vincent_van_Gogh_-_Self-Portrait_-_Google_Art_Project.jpg/842px-Vincent_van_Gogh_-_Self-Portrait_-_Google_Art_Project.jpg) in a decent resolution. \n",
        "* You can use any other image, if you want."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "XTgznuO-LMMH",
        "outputId": "df1f181b-e95e-423f-ddc9-71ae71512ca7",
        "tags": []
      },
      "outputs": [],
      "source": [
        "!curl -o \"cat.jpg\" \"https://raw.githubusercontent.com/Digital-Media/cv_data/main/example_images/cat.jpg\" --silent\n",
        "!curl -o \"gogh.jpg\" \"https://raw.githubusercontent.com/Digital-Media/cv_data/main/example_images/gogh.jpg\" --silent"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# Histogram"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Exercise 1<a id=\"exercise1\" name=\"exercise1\"> </a>📝: Manually compute a histogram!\n",
        "\n",
        "Let's start with a simple exercise. We will compute a histogram of an image manually. \n",
        "For simplicity, we will use an 8-bit grayscale image. \n",
        "You can plot the histogram with plotly as follows:\n",
        "\n",
        "```python\n",
        "fig = px.bar(x=range(256), y=counts, labels={\"x\": \"pixel value\", \"y\": \"count\"})\n",
        "```\n",
        "\n",
        "**(a)** Compute the counts for each pixel value [0 to 255] in the image.\n",
        "\n",
        "Compare your result with the next cell. It features histogram computations with NumPy or Plotly."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "greyscale = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        " \n",
        "\n",
        "# Solution (a)\n",
        "# Todo: Manually generate a histogram of the image\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "As always, NumPy/Python simplifies this a lot. Below you'll find a Histogram computation with NumPy (or Plotly in the comments)."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": []
      },
      "outputs": [],
      "source": [
        "# load the gogh image and convert it to grayscale\n",
        "greyscale = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "\n",
        "\n",
        "# image histogram\n",
        "counts, bins = np.histogram(greyscale.ravel(), bins=range(257))\n",
        "fig = px.bar(x=range(256), y=counts, labels={\"x\": \"pixel value\", \"y\": \"count\"})\n",
        "\n",
        "# optionally use the plotly histogram function directly\n",
        "# fig = px.histogram(x=pixels_grayscale.ravel(), nbins=256, labels={'x':'pixel value', 'y':'count'})\n",
        "\n",
        "fig.show()\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Properties of Histograms\n",
        "\n",
        "Let's look at some properties that can be derived from a histogram. <br>\n",
        "We will plot the histogram and annotate it with the mean, median, and mode. <br>\n",
        "Note: This is very much code, mostly for formatting. You can ignore it, if you want!"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [
          "hide-cell"
        ]
      },
      "outputs": [],
      "source": [
        "from plotly.express.colors import sample_colorscale\n",
        "from plotly.subplots import make_subplots\n",
        "\n",
        "def show_hist_stats(image: np.ndarray, show_stats: bool = True, use_cumulative: bool = False):\n",
        "\n",
        "    x = np.linspace(0, 1, 5)\n",
        "    c = sample_colorscale('HSV', list(x))\n",
        "\n",
        "    # 8-bit (256) image histogram\n",
        "    counts, bins = np.histogram(image.ravel(), bins=range(257))\n",
        "    cumulative = np.cumsum(counts)\n",
        "\n",
        "    fig = px.bar(x=bins[:-1], y=cumulative if use_cumulative else counts, labels={'x':'pixel value', 'y':'count'}, color_discrete_sequence=['black']*256)\n",
        "\n",
        "    fig.update_layout(plot_bgcolor='white', margin=dict(t=0, b=0, r=0, l=0, pad=0))\n",
        "\n",
        "    num_markers = 1000\n",
        "    y_pos = -np.max(cumulative if use_cumulative else counts)*.05\n",
        "    fig.add_traces([\n",
        "        go.Scatter(x=np.linspace(0,255,num_markers), y=[y_pos]*num_markers, mode='markers', marker={'color': np.linspace(0,255,num_markers), 'colorscale': 'gray', 'size': 10, 'symbol': 'square' }),\n",
        "    ])\n",
        "\n",
        "\n",
        "    # show the mean, median, mode and std as vertical lines\n",
        "    mean_value = np.mean(image)\n",
        "    median_value = np.median(image)\n",
        "    std_value = np.std(image)\n",
        "    mode_value = np.argmax(counts)\n",
        "    min_value = np.min(image)\n",
        "    max_value = np.max(image)\n",
        "\n",
        "    if show_stats:\n",
        "        fig.add_trace(\n",
        "            go.Scatter(\n",
        "                x=[mean_value, mean_value],\n",
        "                y=[np.max(cumulative if use_cumulative else counts),0],\n",
        "                mode=\"lines+text\",\n",
        "                line=go.scatter.Line(color=c[0]),\n",
        "                name=\"mean\",\n",
        "                text=[\"Mean\", \"\"],\n",
        "                textposition=\"top center\",\n",
        "                showlegend=True)\n",
        "        )\n",
        "\n",
        "        fig.add_trace(\n",
        "            go.Scatter(\n",
        "                x=[median_value, median_value],\n",
        "                y=[np.max(cumulative if use_cumulative else counts),0],\n",
        "                mode=\"lines+text\",\n",
        "                line=go.scatter.Line(color=c[1]),\n",
        "                name=\"median\",\n",
        "                text=[\"Median\", \"\"],\n",
        "                textposition=\"top center\",\n",
        "                showlegend=True)\n",
        "        )\n",
        "\n",
        "        fig.add_trace(\n",
        "            go.Scatter(\n",
        "                x=[mode_value, mode_value],\n",
        "                y=[np.max(cumulative if use_cumulative else counts),0],\n",
        "                mode=\"lines+text\",\n",
        "                line=go.scatter.Line(color=c[2]),\n",
        "                name=\"mode\",\n",
        "                text=[\"Mode\", \"\"],\n",
        "                textposition=\"top center\",\n",
        "                showlegend=True)\n",
        "        )\n",
        "\n",
        "        fig.add_trace(\n",
        "            go.Scatter(\n",
        "                x=[mean_value - std_value, mean_value + std_value],\n",
        "                y=[np.max(cumulative if use_cumulative else counts)*0.7, np.max(cumulative if use_cumulative else counts)*0.7],\n",
        "                mode=\"text+lines+markers\",\n",
        "                marker_symbol=\"line-ns\",\n",
        "                marker_line_width=1,\n",
        "                marker_line_color=c[0],\n",
        "                marker_size=10,\n",
        "                line=go.scatter.Line(color=c[0]),\n",
        "                name=\"std\",\n",
        "                text=[\"-std\", \"+std\"],\n",
        "                textposition=[\"middle left\",\"middle right\"],\n",
        "                showlegend=True)\n",
        "        )\n",
        "\n",
        "        fig.add_trace(\n",
        "            go.Scatter(\n",
        "                x=[min_value, max_value],\n",
        "                y=[-10, -10],\n",
        "                mode=\"text+lines+markers\",\n",
        "                marker_symbol=\"line-ns\",\n",
        "                marker_line_width=1,\n",
        "                marker_line_color=c[3],\n",
        "                marker_size=10,\n",
        "                line=go.scatter.Line(color=c[3]),\n",
        "                name=\"min/max\",\n",
        "                text=[\"min\", \"max\"],\n",
        "                textposition=[\"middle left\",\"middle right\"],\n",
        "                showlegend=True)\n",
        "        )\n",
        "\n",
        "        if use_cumulative:\n",
        "            fig.add_trace(\n",
        "                go.Scatter(\n",
        "                    x=[0, 255],\n",
        "                    y=[np.max(cumulative if use_cumulative else counts)*0.5, np.max(cumulative if use_cumulative else counts)*0.5],\n",
        "                    mode=\"lines+text\",\n",
        "                    line=go.scatter.Line(color='gray', dash='dot'),\n",
        "                    name=\"fifty_percent\",\n",
        "                    text=[\"50%\", \"\"],\n",
        "                    textposition=\"top center\",\n",
        "                    showlegend=True)\n",
        "            )\n",
        "\n",
        "    return fig\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from plotly.subplots import make_subplots\n",
        "\n",
        "greyscale = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "\n",
        "fig = make_subplots(1, 2)\n",
        "fig.add_trace(go.Image(z=cv2.cvtColor(greyscale, cv2.COLOR_GRAY2BGR), name=\"Image\"), 1, 1)\n",
        "traces = show_hist_stats(greyscale).data\n",
        "for trace in traces:\n",
        "    fig.add_trace(trace, 1, 2)\n",
        "fig.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Exercise 2<a id=\"exercise2\" name=\"exercise2\"> </a>📝: How do point operations change a histogram?\n",
        "\n",
        "What happens if you apply the following point operations to an image? <br>\n",
        "**(a)** multiply the image by a constant $f(a) = 1.5a$ <br>\n",
        "**(b)** add a constant $f(a) = a + 50$ <br>\n",
        "**(c)** invert the image $f(a) = 255 - a$ (for an 8-bit image)<br>\n",
        "\n",
        "Think about what the histogram of the resulting image will look like and then check your answer by plotting the histogram of the resulting image. <br>"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Solutions: Todo"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Color Histograms\n",
        "\n",
        "Just as we can compute a single histogram for a grayscale image, we can compute a histogram for each color channel in a color image."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from plotly.subplots import make_subplots\n",
        "from skimage import data\n",
        "img = cv2.imread('gogh.jpg')\n",
        "img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n",
        "fig = make_subplots(1, 2)\n",
        "# We use go.Image because subplots require traces, whereas px functions return a figure\n",
        "fig.add_trace(go.Image(z=img), 1, 1)\n",
        "for channel, color in enumerate(['red', 'green', 'blue']):\n",
        "    fig.add_trace(go.Histogram(x=img[..., channel].ravel(), opacity=0.5,\n",
        "                               marker_color=color, name=f'{color} channel'), 1, 2)\n",
        "fig.update_layout(height=400)\n",
        "fig.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Cumulative Histograms\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# load the gogh image and convert it to grayscale\n",
        "greyscale = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "\n",
        "\n",
        "# regular histogram\n",
        "counts, bins = np.histogram(greyscale.ravel(), bins=range(257))\n",
        "# cumulative image histogram\n",
        "cumulative = np.cumsum(counts) # <--- this is the computation of the cumulative histogram\n",
        "fig = px.bar(x=range(256), y=cumulative, labels={\"x\": \"pixel value\", \"y\": \"count\"})\n",
        "fig.show()\n",
        "\n",
        "print( f'number of pixels: {np.prod(greyscale.shape[:2])}\\nH[K-1]: {cumulative[-1]}' )\n",
        "print( f'H[median]: {cumulative[np.median(greyscale).astype(np.uint8)]} ~ (number of pixels) / 2: {np.prod(greyscale.shape[:2])/2}' )"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Histogram Equalization\n",
        "\n",
        "Histogram equalization is a technique that can be used to improve the contrast of an image. It is a non-linear point operation that transforms the input image so that the output image has a uniform histogram.\n",
        "Let's try it out on the image of Van Gogh."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "greyscale = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "\n",
        "#OpenCV's histogram equalization\n",
        "equalized = cv2.equalizeHist(greyscale)\n",
        "\n",
        "\n",
        "fig = make_subplots(1, 2)\n",
        "fig.add_trace(go.Image(z=cv2.cvtColor(equalized, cv2.COLOR_GRAY2BGR), name=\"Image\"), 1, 1)\n",
        "traces = show_hist_stats(equalized, use_cumulative=True, show_stats=False).data\n",
        "for trace in traces:\n",
        "    fig.add_trace(trace, 1, 2)\n",
        "fig.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Exercise 3<a id=\"exercise3\" name=\"exercise3\"> </a>📝: Let's manually equalize the histogram of an image!\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "greyscale = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "equalized = greyscale\n",
        "\n",
        "# TODO: Implement histogram equalization\n",
        "\n",
        "fig = make_subplots(1, 2)\n",
        "fig.add_trace(go.Image(z=cv2.cvtColor(equalized, cv2.COLOR_GRAY2BGR), name=\"Image\"), 1, 1)\n",
        "traces = show_hist_stats(equalized, use_cumulative=True, show_stats=False).data\n",
        "for trace in traces:\n",
        "    fig.add_trace(trace, 1, 2)\n",
        "fig.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Histogram Matching\n",
        "\n",
        "Histogram matching is a technique that can be used to transform the histogram of an image to match the histogram of another image. It is a non-linear point operation that transforms the input image so that the output image has the same histogram as the reference image.\n",
        "\n",
        "OpenCV does not have a histogram matching function, but we can use skimage to do it. There is an example [online.](https://scikit-image.org/docs/stable/auto_examples/color_exposure/plot_histogram_matching.html)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "from skimage.exposure import match_histograms\n",
        "\n",
        "reference = cv2.imread(\"gogh.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "image = cv2.imread(\"cat.jpg\", cv2.IMREAD_GRAYSCALE)\n",
        "\n",
        "# using skimage's match_histograms\n",
        "matched = match_histograms(image, reference).astype(np.uint8)\n",
        "\n",
        "# display images\n",
        "imgs = [reference, image, matched]\n",
        "titles = ['Reference', 'Image', 'Matched']\n",
        "\n",
        "fig = make_subplots(2, len(imgs), subplot_titles=titles,\n",
        "    horizontal_spacing = 0.05, vertical_spacing = 0.1)\n",
        "for i, (img, title) in enumerate(zip(imgs, titles)):\n",
        "    fig.add_trace(go.Image(z=cv2.cvtColor(img, cv2.COLOR_GRAY2BGR), name=\"Image\"), 1, i+1)\n",
        "    traces = show_hist_stats(img, use_cumulative=True, show_stats=False).data\n",
        "    for trace in traces:\n",
        "        fig.add_trace(trace, 2, i+1)\n",
        "\n",
        "fig.show()"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "collapsed_sections": [
        "tXV5YxS7lsU2"
      ],
      "name": "03_Histograms.ipynb",
      "provenance": [],
      "toc_visible": true
    },
    "kernelspec": {
      "display_name": "Python 3.8.10 64-bit (microsoft store)",
      "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"
    },
    "varInspector": {
      "cols": {
        "lenName": 16,
        "lenType": 16,
        "lenVar": 40
      },
      "kernels_config": {
        "python": {
          "delete_cmd_postfix": "",
          "delete_cmd_prefix": "del ",
          "library": "var_list.py",
          "varRefreshCmd": "print(var_dic_list())"
        },
        "r": {
          "delete_cmd_postfix": ") ",
          "delete_cmd_prefix": "rm(",
          "library": "var_list.r",
          "varRefreshCmd": "cat(var_dic_list()) "
        }
      },
      "types_to_exclude": [
        "module",
        "function",
        "builtin_function_or_method",
        "instance",
        "_Feature"
      ],
      "window_display": false
    },
    "vscode": {
      "interpreter": {
        "hash": "0bae4c9257a18836fb2e3dc2d0aeb6355625d596c4075009294ab101cd3e0d3c"
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}