{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 代码实现"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-04T08:46:10.663241Z",
     "iopub.status.busy": "2022-07-04T08:46:10.662713Z",
     "iopub.status.idle": "2022-07-04T08:46:11.722391Z",
     "shell.execute_reply": "2022-07-04T08:46:11.721158Z",
     "shell.execute_reply.started": "2022-07-04T08:46:10.663147Z"
    }
   },
   "outputs": [],
   "source": [
    "from utils import *\n",
    "\n",
    "\n",
    "# 计算角点响应函数\n",
    "def responseFunc(M):\n",
    "    # 设置超参数\n",
    "    k = 0.04\n",
    "    \n",
    "    # 计算M的行列式\n",
    "    det = np.linalg.det(M)\n",
    "    # 计算M的迹\n",
    "    trace = np.trace(M)\n",
    "    \n",
    "    # 计算角点响应函数\n",
    "    R = det - k * trace ** 2\n",
    "    \n",
    "    return R\n",
    "\n",
    "\n",
    "# Harris角点检测算法\n",
    "def harris_corners(src, NMS=False):\n",
    "    \n",
    "    # 获得输入图像的长宽\n",
    "    h, w = src.shape[:2]\n",
    "    \n",
    "    # 将图像转为灰度图像\n",
    "    gray_image = cv2.cvtColor(src, cv2.COLOR_BGR2GRAY)\n",
    "    \n",
    "    # 初始化角点矩阵\n",
    "    cornerPoint = np.zeros_like(gray_image,dtype=np.float32)\n",
    "    \n",
    "    # 计算图像沿横轴和纵轴的梯度\n",
    "    grad = np.zeros((h, w, 2), dtype=np.float32)\n",
    "    # x轴梯度\n",
    "    grad[:,:,0] = cv2.Sobel(gray_image, cv2.CV_16S, 1, 0)\n",
    "    # y轴梯度\n",
    "    grad[:,:,1] = cv2.Sobel(gray_image, cv2.CV_16S, 0, 1)\n",
    "    \n",
    "    # 计算黑塞矩阵内元素的值，此时的值是求和的结果，即M矩阵的元素\n",
    "    Ixx = grad[:,:,0] ** 2\n",
    "    Iyy = grad[:,:,1] ** 2\n",
    "    Ixy = grad[:,:,0] * grad[:,:,1]\n",
    "     \n",
    "    # 计算窗口内黑塞矩阵元素的值，窗函数使用高斯函数\n",
    "    Ixx = cv2.GaussianBlur(Ixx, (3, 3), sigmaX=2)\n",
    "    Iyy = cv2.GaussianBlur(Iyy, (3, 3), sigmaX=2)\n",
    "    Ixy = cv2.GaussianBlur(Ixy, (3, 3), sigmaX=2)\n",
    "    \n",
    "    for i in range(gray_image.shape[0]):\n",
    "        for j in range(gray_image.shape[1]):\n",
    "            \n",
    "            # 构建M矩阵\n",
    "            struture_matrix = [[Ixx[i][j], Ixy[i][j]], \n",
    "                               [Ixy[i][j], Iyy[i][j]]]\n",
    "\n",
    "            # 计算角点响应函数\n",
    "            R = responseFunc(struture_matrix)\n",
    "            cornerPoint[i][j] = R\n",
    "\n",
    "    # 非最大抑制\n",
    "    corners = np.zeros_like(gray_image, dtype=np.float32)\n",
    "    threshold = 0.01\n",
    "    \n",
    "    # 返回所有角点响应的最大值\n",
    "    maxValue = np.max(cornerPoint)\n",
    "    \n",
    "    # 我们将角点响应函数的阈值设定为 threshold * maxValue\n",
    "    \n",
    "    for i in range(cornerPoint.shape[0]):\n",
    "        for j in range(cornerPoint.shape[1]):\n",
    "            \n",
    "            # 如果进行NMS操作\n",
    "            if NMS:\n",
    "                # 当前角点响应值大于阈值，同时也是邻居周围的最大值\n",
    "                if cornerPoint[i][j] > threshold * maxValue and \\\n",
    "                    cornerPoint[i][j] == np.max(\n",
    "                    cornerPoint[max(0, i - 1):min(i + 1, h - 1), \n",
    "                    max(0, j - 1):min(j + 1, w - 1)]):\n",
    "                    \n",
    "                    corners[i][j] = 255\n",
    "            else:\n",
    "                # 当前角点响应值大于阈值\n",
    "                if cornerPoint[i][j] > threshold * maxValue:\n",
    "                    corners[i][j] = 255\n",
    "                    \n",
    "    # 返回检测到的角点\n",
    "    return corners\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2022-07-04T08:46:11.724725Z",
     "iopub.status.busy": "2022-07-04T08:46:11.724326Z",
     "iopub.status.idle": "2022-07-04T08:46:28.125897Z",
     "shell.execute_reply": "2022-07-04T08:46:28.124575Z",
     "shell.execute_reply.started": "2022-07-04T08:46:11.724692Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 173684 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 193280 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# 图像素材来自opencvpython网站\n",
    "img = cv2.imread('sudoku.png')\n",
    "\n",
    "plt.imshow(img[:, :, ::-1])\n",
    "plt.axis('off')\n",
    "plt.title('Original Input')\n",
    "plt.show()\n",
    "\n",
    "response = harris_corners(img)\n",
    "img[response == 255] = [0, 0, 255]\n",
    "plt.imshow(img[:, :, ::-1])\n",
    "plt.axis('off')\n",
    "plt.title('Harris Corner Response')\n",
    "plt.show()"
   ]
  }
 ],
 "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.8.16"
  },
  "vscode": {
   "interpreter": {
    "hash": "31f2aee4e71d21fbe5cf8b01ff0e069b9275f58929596ceb00d14d90e3e16cd6"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
