1""" 2======================== 3Sliding window histogram 4======================== 5 6Histogram matching can be used for object detection in images [1]_. This 7example extracts a single coin from the ``skimage.data.coins`` image and uses 8histogram matching to attempt to locate it within the original image. 9 10First, a box-shaped region of the image containing the target coin is 11extracted and a histogram of its grayscale values is computed. 12 13Next, for each pixel in the test image, a histogram of the grayscale values in 14a region of the image surrounding the pixel is computed. 15``skimage.filters.rank.windowed_histogram`` is used for this task, as it employs 16an efficient sliding window based algorithm that is able to compute these 17histograms quickly [2]_. The local histogram for the region surrounding each 18pixel in the image is compared to that of the single coin, with a similarity 19measure being computed and displayed. 20 21The histogram of the single coin is computed using ``numpy.histogram`` on a box 22shaped region surrounding the coin, while the sliding window histograms are 23computed using a disc shaped structural element of a slightly different size. 24This is done in aid of demonstrating that the technique still finds similarity 25in spite of these differences. 26 27To demonstrate the rotational invariance of the technique, the same test is 28performed on a version of the coins image rotated by 45 degrees. 29 30References 31---------- 32.. [1] Porikli, F. "Integral Histogram: A Fast Way to Extract Histograms 33 in Cartesian Spaces" CVPR, 2005. Vol. 1. IEEE, 2005 34.. [2] S.Perreault and P.Hebert. Median filtering in constant time. 35 Trans. Image Processing, 16(9):2389-2394, 2007. 36""" 37 38import numpy as np 39import matplotlib 40import matplotlib.pyplot as plt 41 42from skimage import data, transform 43from skimage.util import img_as_ubyte 44from skimage.morphology import disk 45from skimage.filters import rank 46 47 48matplotlib.rcParams['font.size'] = 9 49 50 51def windowed_histogram_similarity(image, footprint, reference_hist, n_bins): 52 # Compute normalized windowed histogram feature vector for each pixel 53 px_histograms = rank.windowed_histogram(image, footprint, n_bins=n_bins) 54 55 # Reshape coin histogram to (1,1,N) for broadcast when we want to use it in 56 # arithmetic operations with the windowed histograms from the image 57 reference_hist = reference_hist.reshape((1, 1) + reference_hist.shape) 58 59 # Compute Chi squared distance metric: sum((X-Y)^2 / (X+Y)); 60 # a measure of distance between histograms 61 X = px_histograms 62 Y = reference_hist 63 64 num = (X - Y) ** 2 65 denom = X + Y 66 denom[denom == 0] = np.infty 67 frac = num / denom 68 69 chi_sqr = 0.5 * np.sum(frac, axis=2) 70 71 # Generate a similarity measure. It needs to be low when distance is high 72 # and high when distance is low; taking the reciprocal will do this. 73 # Chi squared will always be >= 0, add small value to prevent divide by 0. 74 similarity = 1 / (chi_sqr + 1.0e-4) 75 76 return similarity 77 78 79# Load the `skimage.data.coins` image 80img = img_as_ubyte(data.coins()) 81 82# Quantize to 16 levels of grayscale; this way the output image will have a 83# 16-dimensional feature vector per pixel 84quantized_img = img // 16 85 86# Select the coin from the 4th column, second row. 87# Co-ordinate ordering: [x1,y1,x2,y2] 88coin_coords = [184, 100, 228, 148] # 44 x 44 region 89coin = quantized_img[coin_coords[1]:coin_coords[3], 90 coin_coords[0]:coin_coords[2]] 91 92# Compute coin histogram and normalize 93coin_hist, _ = np.histogram(coin.flatten(), bins=16, range=(0, 16)) 94coin_hist = coin_hist.astype(float) / np.sum(coin_hist) 95 96# Compute a disk shaped mask that will define the shape of our sliding window 97# Example coin is ~44px across, so make a disk 61px wide (2 * rad + 1) to be 98# big enough for other coins too. 99footprint = disk(30) 100 101# Compute the similarity across the complete image 102similarity = windowed_histogram_similarity(quantized_img, footprint, coin_hist, 103 coin_hist.shape[0]) 104 105# Now try a rotated image 106rotated_img = img_as_ubyte(transform.rotate(img, 45.0, resize=True)) 107# Quantize to 16 levels as before 108quantized_rotated_image = rotated_img // 16 109# Similarity on rotated image 110rotated_similarity = windowed_histogram_similarity(quantized_rotated_image, 111 footprint, coin_hist, 112 coin_hist.shape[0]) 113 114fig, axes = plt.subplots(nrows=2, ncols=2, figsize=(10, 10)) 115 116axes[0, 0].imshow(quantized_img, cmap='gray') 117axes[0, 0].set_title('Quantized image') 118axes[0, 0].axis('off') 119 120axes[0, 1].imshow(coin, cmap='gray') 121axes[0, 1].set_title('Coin from 2nd row, 4th column') 122axes[0, 1].axis('off') 123 124axes[1, 0].imshow(img, cmap='gray') 125axes[1, 0].imshow(similarity, cmap='hot', alpha=0.5) 126axes[1, 0].set_title('Original image with overlaid similarity') 127axes[1, 0].axis('off') 128 129axes[1, 1].imshow(rotated_img, cmap='gray') 130axes[1, 1].imshow(rotated_similarity, cmap='hot', alpha=0.5) 131axes[1, 1].set_title('Rotated image with overlaid similarity') 132axes[1, 1].axis('off') 133 134plt.tight_layout() 135plt.show() 136