{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.10","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":49349,"databundleVersionId":5447706,"sourceType":"competition"}],"dockerImageVersionId":30474,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Image matching Challenge\n\nAs a robotics engineer 3D reconstruction is an interesting problem.\nEspecially one of my interest, as the perception side is my focus of work","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19"}},{"cell_type":"code","source":"mode = \"development\"\nmode = \"submission\"","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:29:19.325350Z","iopub.execute_input":"2024-03-05T09:29:19.325870Z","iopub.status.idle":"2024-03-05T09:29:19.337337Z","shell.execute_reply.started":"2024-03-05T09:29:19.325830Z","shell.execute_reply":"2024-03-05T09:29:19.335803Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"if mode != \"submission\":\n    !pip install open3d","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:29:19.846739Z","iopub.execute_input":"2024-03-05T09:29:19.847817Z","iopub.status.idle":"2024-03-05T09:29:19.854299Z","shell.execute_reply.started":"2024-03-05T09:29:19.847766Z","shell.execute_reply":"2024-03-05T09:29:19.852629Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import seaborn as sns\nimport matplotlib.pyplot as plt\nimport matplotlib as mpl\nimport cv2\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport os\nimport pycolmap\nif mode != \"submission\":\n    import open3d as o3d\nfrom tqdm import tqdm","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:29:21.163263Z","iopub.execute_input":"2024-03-05T09:29:21.163747Z","iopub.status.idle":"2024-03-05T09:29:21.173403Z","shell.execute_reply.started":"2024-03-05T09:29:21.163712Z","shell.execute_reply":"2024-03-05T09:29:21.170834Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Discover a single object","metadata":{}},{"cell_type":"code","source":"category = \"phototourism\"\nobject_name = \"st_peters_square\"\npath_object = f\"/kaggle/input/image-matching-challenge-2023/train/{category}/{object_name}/images/\"\nimage_files = os.listdir(path_object)\nlen(image_files)","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:29:23.868411Z","iopub.execute_input":"2024-03-05T09:29:23.868922Z","iopub.status.idle":"2024-03-05T09:29:24.173123Z","shell.execute_reply.started":"2024-03-05T09:29:23.868885Z","shell.execute_reply":"2024-03-05T09:29:24.171785Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 查看图片","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(16,12))\nf, axarr = plt.subplots(3, 5) \n\ni = 0\nj = 0\n\nfor file in image_files[:15]:\n    im = plt.imread(path_object+file, format='jpeg')\n    dim = (im.shape[1] // 2, im.shape[0] // 2)\n    im = cv2.resize(im, dim, interpolation = cv2.INTER_AREA)\n    axarr[i,j].imshow(im)\n    \n    i += 1\n    if i == 3:\n        j += 1\n        i = 0\n","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:30:10.484432Z","iopub.execute_input":"2024-03-05T09:30:10.485032Z","iopub.status.idle":"2024-03-05T09:30:13.563781Z","shell.execute_reply.started":"2024-03-05T09:30:10.484990Z","shell.execute_reply":"2024-03-05T09:30:13.561987Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"## Select two similair images\nimages_compare = [image_files[3], image_files[12]]","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:30:20.385365Z","iopub.execute_input":"2024-03-05T09:30:20.385879Z","iopub.status.idle":"2024-03-05T09:30:20.391601Z","shell.execute_reply.started":"2024-03-05T09:30:20.385845Z","shell.execute_reply":"2024-03-05T09:30:20.390718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# 计算SIFT特征\n","metadata":{}},{"cell_type":"code","source":"from PIL import Image, ImageOps\n\ndescriptor_list = []\nkeypoint_list = []\nimgs = []\nmax_features = 512\n\nfor file in images_compare:\n    im = Image.open(path_object+file).convert('RGB')\n    img = ImageOps.grayscale(im)\n    img = np.array(img).astype(np.float) / 255.\n\n    # Optional parameters:\n    # - options: dict or pycolmap.SiftExtractionOptions\n    # - device: default pycolmap.Device.auto uses the GPU if available\n    sift = pycolmap.Sift()\n\n    # Parameters:\n    # - image: HxW float array\n    keypoints, scores, descriptors = sift.extract(img)\n    scores_idx = np.argsort(scores)\n    descriptor_list.append(descriptors[scores_idx[-(max_features):]])\n    keypoint_list.append(keypoints[scores_idx[-(max_features):]])\n\n    im = np.array(im)\n    im_append = im.copy()\n    imgs.append(im_append)\n","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:34:58.367369Z","iopub.execute_input":"2024-03-05T09:34:58.367987Z","iopub.status.idle":"2024-03-05T09:34:59.698899Z","shell.execute_reply.started":"2024-03-05T09:34:58.367941Z","shell.execute_reply":"2024-03-05T09:34:59.697744Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"keypoints","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:35:22.535408Z","iopub.execute_input":"2024-03-05T09:35:22.536037Z","iopub.status.idle":"2024-03-05T09:35:22.545921Z","shell.execute_reply.started":"2024-03-05T09:35:22.535990Z","shell.execute_reply":"2024-03-05T09:35:22.544454Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"scores","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:35:59.054225Z","iopub.execute_input":"2024-03-05T09:35:59.055586Z","iopub.status.idle":"2024-03-05T09:35:59.063591Z","shell.execute_reply.started":"2024-03-05T09:35:59.055533Z","shell.execute_reply":"2024-03-05T09:35:59.062420Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"descriptors","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:38:16.680392Z","iopub.execute_input":"2024-03-05T09:38:16.680856Z","iopub.status.idle":"2024-03-05T09:38:16.689907Z","shell.execute_reply.started":"2024-03-05T09:38:16.680822Z","shell.execute_reply":"2024-03-05T09:38:16.688699Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 计算图像之间的loss","metadata":{}},{"cell_type":"code","source":"print(descriptor_list[0].shape) # 图0描述子\nprint(descriptor_list[1].shape) # 图1描述子\n\n\nscore_matrix = np.zeros((descriptor_list[0].shape[0], descriptor_list[1].shape[0]))\n\nmatching_keypoints = []\nerrors = []\nerror_threshold = 0.1\n\nfrom tqdm import tqdm\nfor i in tqdm(range(score_matrix.shape[0])):\n    for j in range(score_matrix.shape[1]):\n        if abs(np.sum((descriptor_list[0][i] - descriptor_list[1][j])**2)) < error_threshold:  # 描述子欧几里得距离小于<0.1，接近\n            matching_keypoints.append([i,j])\n\nprint(len(matching_keypoints), score_matrix.shape[0]*score_matrix.shape[1])\n# 262144对图像中，40个相似","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:41:51.800558Z","iopub.execute_input":"2024-03-05T09:41:51.802124Z","iopub.status.idle":"2024-03-05T09:41:55.625034Z","shell.execute_reply.started":"2024-03-05T09:41:51.802072Z","shell.execute_reply":"2024-03-05T09:41:55.623877Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"score_matrix","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:41:55.627419Z","iopub.execute_input":"2024-03-05T09:41:55.633492Z","iopub.status.idle":"2024-03-05T09:41:55.643013Z","shell.execute_reply.started":"2024-03-05T09:41:55.633448Z","shell.execute_reply":"2024-03-05T09:41:55.641872Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## plotting the similair pixels\n","metadata":{}},{"cell_type":"code","source":"# Radius of circle\nradius = 5\n# Line thickness of 2 px\nthickness = 2\n\ncolormaps = mpl.colormaps['tab20c']\nnorm = mpl.colors.Normalize(vmin=0, vmax=(len(matching_keypoints)-1))\n\nif imgs[0].shape[0] > imgs[1].shape[0]:\n    add = (imgs[0].shape[0] - imgs[1].shape[0])\n    pad = np.zeros((add, imgs[1].shape[1], 3))*255\n    imgs[1] = np.append(imgs[1], pad, axis=0)\nelif imgs[1].shape[0] > imgs[0].shape[0]:\n    add = (imgs[1].shape[0] - imgs[0].shape[0])\n    pad = np.zeros((add, imgs[0].shape[1], 3))*255\n    imgs[0] = np.append(imgs[0], pad, axis=0)\n\n\ncombined_image = np.append(imgs[0],imgs[1], axis=1).astype(np.int32)\nwidth = imgs[0].shape[1]\n\ni=0\nfor point in matching_keypoints:\n    im1_point_id = descriptor_list[0][point[0]]\n    im2_point_id = descriptor_list[1][point[1]]\n    im1_point = np.array([int(keypoint_list[0][point[0]][0]), int(keypoint_list[0][point[0]][1])])\n    im2_point = np.array([int(keypoint_list[1][point[1]][0]), int(keypoint_list[1][point[1]][1])])\n    im2_point[0] += width\n\n    color = colormaps(norm(i))\n    color = (int(color[0]*255), int(color[1]*255), int(color[2]*255))\n\n    cv2.circle(combined_image, (im1_point[0], im1_point[1]), radius, color, thickness)\n    cv2.circle(combined_image, (im2_point[0], im2_point[1]), radius, color, thickness)\n    cv2.line(combined_image, (im1_point[0], im1_point[1]), (int(im2_point[0]), int(im2_point[1])), color, thickness)\n    i += 1\n\nplt.figure(figsize=(16,8))\nplt.imshow(combined_image)","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:44:48.820279Z","iopub.execute_input":"2024-03-05T09:44:48.820820Z","iopub.status.idle":"2024-03-05T09:44:49.950617Z","shell.execute_reply.started":"2024-03-05T09:44:48.820783Z","shell.execute_reply":"2024-03-05T09:44:49.949122Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"images_compare = [image_files[3], image_files[13]]\n\nfrom PIL import Image, ImageOps\n\ndescriptor_list = []\nkeypoint_list = []\nimgs = []\nmax_features = 512\n\nfor file in images_compare:\n    im = Image.open(path_object+file).convert('RGB')\n    img = ImageOps.grayscale(im)\n    img = np.array(img).astype(np.float) / 255.\n\n    # Optional parameters:\n    # - options: dict or pycolmap.SiftExtractionOptions\n    # - device: default pycolmap.Device.auto uses the GPU if available\n    sift = pycolmap.Sift()\n\n    # Parameters:\n    # - image: HxW float array\n    keypoints, scores, descriptors = sift.extract(img)\n    scores_idx = np.argsort(scores)\n    descriptor_list.append(descriptors[scores_idx[-(max_features):]])\n    keypoint_list.append(keypoints[scores_idx[-(max_features):]])\n\n    im = np.array(im)\n    im_append = im.copy()\n    imgs.append(im_append)\n\n    \nprint(descriptor_list[0].shape) # 图0描述子\nprint(descriptor_list[1].shape) # 图1描述子\n\n\nscore_matrix = np.zeros((descriptor_list[0].shape[0], descriptor_list[1].shape[0]))\n\nmatching_keypoints = []\nerrors = []\nerror_threshold = 0.1\n\nfrom tqdm import tqdm\nfor i in tqdm(range(score_matrix.shape[0])):\n    for j in range(score_matrix.shape[1]):\n        if abs(np.sum((descriptor_list[0][i] - descriptor_list[1][j])**2)) < error_threshold:  # 描述子欧几里得距离小于<0.1，接近\n            matching_keypoints.append([i,j])\n\nprint(len(matching_keypoints), score_matrix.shape[0]*score_matrix.shape[1])\n# 262144对图像中，40个相似\n\n\n# Radius of circle\nradius = 5\n# Line thickness of 2 px\nthickness = 2\n\ncolormaps = mpl.colormaps['tab20c']\nnorm = mpl.colors.Normalize(vmin=0, vmax=(len(matching_keypoints)-1))\n\nif imgs[0].shape[0] > imgs[1].shape[0]:\n    add = (imgs[0].shape[0] - imgs[1].shape[0])\n    pad = np.zeros((add, imgs[1].shape[1], 3))*255\n    imgs[1] = np.append(imgs[1], pad, axis=0)\nelif imgs[1].shape[0] > imgs[0].shape[0]:\n    add = (imgs[1].shape[0] - imgs[0].shape[0])\n    pad = np.zeros((add, imgs[0].shape[1], 3))*255\n    imgs[0] = np.append(imgs[0], pad, axis=0)\n\n\ncombined_image = np.append(imgs[0],imgs[1], axis=1).astype(np.int32)\nwidth = imgs[0].shape[1]\n\ni=0\nfor point in matching_keypoints:\n    im1_point_id = descriptor_list[0][point[0]]\n    im2_point_id = descriptor_list[1][point[1]]\n    im1_point = np.array([int(keypoint_list[0][point[0]][0]), int(keypoint_list[0][point[0]][1])])\n    im2_point = np.array([int(keypoint_list[1][point[1]][0]), int(keypoint_list[1][point[1]][1])])\n    im2_point[0] += width\n\n    color = colormaps(norm(i))\n    color = (int(color[0]*255), int(color[1]*255), int(color[2]*255))\n\n    cv2.circle(combined_image, (im1_point[0], im1_point[1]), radius, color, thickness)\n    cv2.circle(combined_image, (im2_point[0], im2_point[1]), radius, color, thickness)\n    cv2.line(combined_image, (im1_point[0], im1_point[1]), (int(im2_point[0]), int(im2_point[1])), color, thickness)\n    i += 1\n\nplt.figure(figsize=(16,8))\nplt.imshow(combined_image)","metadata":{"execution":{"iopub.status.busy":"2024-03-05T09:49:21.990291Z","iopub.execute_input":"2024-03-05T09:49:21.990869Z","iopub.status.idle":"2024-03-05T09:49:28.396844Z","shell.execute_reply.started":"2024-03-05T09:49:21.990827Z","shell.execute_reply":"2024-03-05T09:49:28.395409Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}