_hog.py 13 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341
  1. import numpy as np
  2. from . import _hoghistogram
  3. from .._shared import utils
  4. def _hog_normalize_block(block, method, eps=1e-5):
  5. if method == 'L1':
  6. out = block / (np.sum(np.abs(block)) + eps)
  7. elif method == 'L1-sqrt':
  8. out = np.sqrt(block / (np.sum(np.abs(block)) + eps))
  9. elif method == 'L2':
  10. out = block / np.sqrt(np.sum(block**2) + eps**2)
  11. elif method == 'L2-Hys':
  12. out = block / np.sqrt(np.sum(block**2) + eps**2)
  13. out = np.minimum(out, 0.2)
  14. out = out / np.sqrt(np.sum(out**2) + eps**2)
  15. else:
  16. raise ValueError('Selected block normalization method is invalid.')
  17. return out
  18. def _hog_channel_gradient(channel):
  19. """Compute unnormalized gradient image along `row` and `col` axes.
  20. Parameters
  21. ----------
  22. channel : (M, N) ndarray
  23. Grayscale image or one of image channel.
  24. Returns
  25. -------
  26. g_row, g_col : channel gradient along `row` and `col` axes correspondingly.
  27. """
  28. g_row = np.empty(channel.shape, dtype=channel.dtype)
  29. g_row[0, :] = 0
  30. g_row[-1, :] = 0
  31. g_row[1:-1, :] = channel[2:, :] - channel[:-2, :]
  32. g_col = np.empty(channel.shape, dtype=channel.dtype)
  33. g_col[:, 0] = 0
  34. g_col[:, -1] = 0
  35. g_col[:, 1:-1] = channel[:, 2:] - channel[:, :-2]
  36. return g_row, g_col
  37. @utils.channel_as_last_axis(multichannel_output=False)
  38. def hog(
  39. image,
  40. orientations=9,
  41. pixels_per_cell=(8, 8),
  42. cells_per_block=(3, 3),
  43. block_norm='L2-Hys',
  44. visualize=False,
  45. transform_sqrt=False,
  46. feature_vector=True,
  47. *,
  48. channel_axis=None,
  49. ):
  50. """Extract Histogram of Oriented Gradients (HOG) for a given image.
  51. Compute a Histogram of Oriented Gradients (HOG) by
  52. 1. (optional) global image normalization
  53. 2. computing the gradient image in `row` and `col`
  54. 3. computing gradient histograms
  55. 4. normalizing across blocks
  56. 5. flattening into a feature vector
  57. Parameters
  58. ----------
  59. image : (M, N[, C]) ndarray
  60. Input image.
  61. orientations : int, optional
  62. Number of orientation bins.
  63. pixels_per_cell : 2-tuple (int, int), optional
  64. Size (in pixels) of a cell.
  65. cells_per_block : 2-tuple (int, int), optional
  66. Number of cells in each block.
  67. block_norm : str {'L1', 'L1-sqrt', 'L2', 'L2-Hys'}, optional
  68. Block normalization method:
  69. ``L1``
  70. Normalization using L1-norm.
  71. ``L1-sqrt``
  72. Normalization using L1-norm, followed by square root.
  73. ``L2``
  74. Normalization using L2-norm.
  75. ``L2-Hys``
  76. Normalization using L2-norm, followed by limiting the
  77. maximum values to 0.2 (`Hys` stands for `hysteresis`) and
  78. renormalization using L2-norm. (default)
  79. For details, see [3]_, [4]_.
  80. visualize : bool, optional
  81. Also return an image of the HOG. For each cell and orientation bin,
  82. the image contains a line segment that is centered at the cell center,
  83. is perpendicular to the midpoint of the range of angles spanned by the
  84. orientation bin, and has intensity proportional to the corresponding
  85. histogram value.
  86. transform_sqrt : bool, optional
  87. Apply power law compression to normalize the image before
  88. processing. DO NOT use this if the image contains negative
  89. values. Also see `notes` section below.
  90. feature_vector : bool, optional
  91. Return the data as a feature vector by calling .ravel() on the result
  92. just before returning.
  93. channel_axis : int or None, optional
  94. If None, the image is assumed to be a grayscale (single channel) image.
  95. Otherwise, this parameter indicates which axis of the array corresponds
  96. to channels.
  97. .. versionadded:: 0.19
  98. `channel_axis` was added in 0.19.
  99. Returns
  100. -------
  101. out : (n_blocks_row, n_blocks_col, n_cells_row, n_cells_col, n_orient) ndarray
  102. HOG descriptor for the image. If `feature_vector` is True, a 1D
  103. (flattened) array is returned.
  104. hog_image : (M, N) ndarray, optional
  105. A visualisation of the HOG image. Only provided if `visualize` is True.
  106. Raises
  107. ------
  108. ValueError
  109. If the image is too small given the values of pixels_per_cell and
  110. cells_per_block.
  111. References
  112. ----------
  113. .. [1] https://en.wikipedia.org/wiki/Histogram_of_oriented_gradients
  114. .. [2] Dalal, N and Triggs, B, Histograms of Oriented Gradients for
  115. Human Detection, IEEE Computer Society Conference on Computer
  116. Vision and Pattern Recognition 2005 San Diego, CA, USA,
  117. https://lear.inrialpes.fr/people/triggs/pubs/Dalal-cvpr05.pdf,
  118. :DOI:`10.1109/CVPR.2005.177`
  119. .. [3] Lowe, D.G., Distinctive image features from scale-invatiant
  120. keypoints, International Journal of Computer Vision (2004) 60: 91,
  121. http://www.cs.ubc.ca/~lowe/papers/ijcv04.pdf,
  122. :DOI:`10.1023/B:VISI.0000029664.99615.94`
  123. .. [4] Dalal, N, Finding People in Images and Videos,
  124. Human-Computer Interaction [cs.HC], Institut National Polytechnique
  125. de Grenoble - INPG, 2006,
  126. https://tel.archives-ouvertes.fr/tel-00390303/file/NavneetDalalThesis.pdf
  127. Notes
  128. -----
  129. The presented code implements the HOG extraction method from [2]_ with
  130. the following changes: (I) blocks of (3, 3) cells are used ((2, 2) in the
  131. paper); (II) no smoothing within cells (Gaussian spatial window with sigma=8pix
  132. in the paper); (III) L1 block normalization is used (L2-Hys in the paper).
  133. Power law compression, also known as Gamma correction, is used to reduce
  134. the effects of shadowing and illumination variations. The compression makes
  135. the dark regions lighter. When the kwarg `transform_sqrt` is set to
  136. ``True``, the function computes the square root of each color channel
  137. and then applies the hog algorithm to the image.
  138. """
  139. image = np.atleast_2d(image)
  140. float_dtype = utils._supported_float_type(image.dtype)
  141. image = image.astype(float_dtype, copy=False)
  142. multichannel = channel_axis is not None
  143. ndim_spatial = image.ndim - 1 if multichannel else image.ndim
  144. if ndim_spatial != 2:
  145. raise ValueError(
  146. 'Only images with two spatial dimensions are '
  147. 'supported. If using with color/multichannel '
  148. 'images, specify `channel_axis`.'
  149. )
  150. """
  151. The first stage applies an optional global image normalization
  152. equalisation that is designed to reduce the influence of illumination
  153. effects. In practice we use gamma (power law) compression, either
  154. computing the square root or the log of each color channel.
  155. Image texture strength is typically proportional to the local surface
  156. illumination so this compression helps to reduce the effects of local
  157. shadowing and illumination variations.
  158. """
  159. if transform_sqrt:
  160. image = np.sqrt(image)
  161. """
  162. The second stage computes first order image gradients. These capture
  163. contour, silhouette and some texture information, while providing
  164. further resistance to illumination variations. The locally dominant
  165. color channel is used, which provides color invariance to a large
  166. extent. Variant methods may also include second order image derivatives,
  167. which act as primitive bar detectors - a useful feature for capturing,
  168. e.g. bar like structures in bicycles and limbs in humans.
  169. """
  170. if multichannel:
  171. g_row_by_ch = np.empty_like(image, dtype=float_dtype)
  172. g_col_by_ch = np.empty_like(image, dtype=float_dtype)
  173. g_magn = np.empty_like(image, dtype=float_dtype)
  174. for idx_ch in range(image.shape[2]):
  175. (
  176. g_row_by_ch[:, :, idx_ch],
  177. g_col_by_ch[:, :, idx_ch],
  178. ) = _hog_channel_gradient(image[:, :, idx_ch])
  179. g_magn[:, :, idx_ch] = np.hypot(
  180. g_row_by_ch[:, :, idx_ch], g_col_by_ch[:, :, idx_ch]
  181. )
  182. # For each pixel select the channel with the highest gradient magnitude
  183. idcs_max = g_magn.argmax(axis=2)
  184. rr, cc = np.meshgrid(
  185. np.arange(image.shape[0]),
  186. np.arange(image.shape[1]),
  187. indexing='ij',
  188. sparse=True,
  189. )
  190. g_row = g_row_by_ch[rr, cc, idcs_max]
  191. g_col = g_col_by_ch[rr, cc, idcs_max]
  192. else:
  193. g_row, g_col = _hog_channel_gradient(image)
  194. """
  195. The third stage aims to produce an encoding that is sensitive to
  196. local image content while remaining resistant to small changes in
  197. pose or appearance. The adopted method pools gradient orientation
  198. information locally in the same way as the SIFT [Lowe 2004]
  199. feature. The image window is divided into small spatial regions,
  200. called "cells". For each cell we accumulate a local 1-D histogram
  201. of gradient or edge orientations over all the pixels in the
  202. cell. This combined cell-level 1-D histogram forms the basic
  203. "orientation histogram" representation. Each orientation histogram
  204. divides the gradient angle range into a fixed number of
  205. predetermined bins. The gradient magnitudes of the pixels in the
  206. cell are used to vote into the orientation histogram.
  207. """
  208. s_row, s_col = image.shape[:2]
  209. c_row, c_col = pixels_per_cell
  210. b_row, b_col = cells_per_block
  211. n_cells_row = int(s_row // c_row) # number of cells along row-axis
  212. n_cells_col = int(s_col // c_col) # number of cells along col-axis
  213. # compute orientations integral images
  214. orientation_histogram = np.zeros(
  215. (n_cells_row, n_cells_col, orientations), dtype=float
  216. )
  217. g_row = g_row.astype(float, copy=False)
  218. g_col = g_col.astype(float, copy=False)
  219. _hoghistogram.hog_histograms(
  220. g_col,
  221. g_row,
  222. c_col,
  223. c_row,
  224. s_col,
  225. s_row,
  226. n_cells_col,
  227. n_cells_row,
  228. orientations,
  229. orientation_histogram,
  230. )
  231. # now compute the histogram for each cell
  232. hog_image = None
  233. if visualize:
  234. from .. import draw
  235. radius = min(c_row, c_col) // 2 - 1
  236. orientations_arr = np.arange(orientations)
  237. # set dr_arr, dc_arr to correspond to midpoints of orientation bins
  238. orientation_bin_midpoints = np.pi * (orientations_arr + 0.5) / orientations
  239. dr_arr = radius * np.sin(orientation_bin_midpoints)
  240. dc_arr = radius * np.cos(orientation_bin_midpoints)
  241. hog_image = np.zeros((s_row, s_col), dtype=float_dtype)
  242. for r in range(n_cells_row):
  243. for c in range(n_cells_col):
  244. for o, dr, dc in zip(orientations_arr, dr_arr, dc_arr):
  245. centre = tuple([r * c_row + c_row // 2, c * c_col + c_col // 2])
  246. rr, cc = draw.line(
  247. int(centre[0] - dc),
  248. int(centre[1] + dr),
  249. int(centre[0] + dc),
  250. int(centre[1] - dr),
  251. )
  252. hog_image[rr, cc] += orientation_histogram[r, c, o]
  253. """
  254. The fourth stage computes normalization, which takes local groups of
  255. cells and contrast normalizes their overall responses before passing
  256. to next stage. Normalization introduces better invariance to illumination,
  257. shadowing, and edge contrast. It is performed by accumulating a measure
  258. of local histogram "energy" over local groups of cells that we call
  259. "blocks". The result is used to normalize each cell in the block.
  260. Typically each individual cell is shared between several blocks, but
  261. its normalizations are block dependent and thus different. The cell
  262. thus appears several times in the final output vector with different
  263. normalizations. This may seem redundant but it improves the performance.
  264. We refer to the normalized block descriptors as Histogram of Oriented
  265. Gradient (HOG) descriptors.
  266. """
  267. n_blocks_row = (n_cells_row - b_row) + 1
  268. n_blocks_col = (n_cells_col - b_col) + 1
  269. if n_blocks_col <= 0 or n_blocks_row <= 0:
  270. min_row = b_row * c_row
  271. min_col = b_col * c_col
  272. raise ValueError(
  273. 'The input image is too small given the values of '
  274. 'pixels_per_cell and cells_per_block. '
  275. 'It should have at least: '
  276. f'{min_row} rows and {min_col} cols.'
  277. )
  278. normalized_blocks = np.zeros(
  279. (n_blocks_row, n_blocks_col, b_row, b_col, orientations), dtype=float_dtype
  280. )
  281. for r in range(n_blocks_row):
  282. for c in range(n_blocks_col):
  283. block = orientation_histogram[r : r + b_row, c : c + b_col, :]
  284. normalized_blocks[r, c, :] = _hog_normalize_block(block, method=block_norm)
  285. """
  286. The final step collects the HOG descriptors from all blocks of a dense
  287. overlapping grid of blocks covering the detection window into a combined
  288. feature vector for use in the window classifier.
  289. """
  290. if feature_vector:
  291. normalized_blocks = normalized_blocks.ravel()
  292. if visualize:
  293. return normalized_blocks, hog_image
  294. else:
  295. return normalized_blocks