Merge pull request #1 from amueller/hog

Hog
This commit is contained in:
Brian Holt
2011-10-05 04:55:10 -07:00
9 changed files with 236 additions and 89 deletions
+3
View File
@@ -63,3 +63,6 @@
- Brian Holt
Histograms of Oriented Gradients
- David-Warde Farley, Sturla Molden
Bresenheim line drawing, from snippets on numpy-discussion.
+4
View File
@@ -122,3 +122,7 @@ Implement Algorithms
- Graph cut segmentation
- Probabilistic Hough transform
Drawing
```````
- Wu's algorithm for lines and circles
+1
View File
@@ -0,0 +1 @@
from draw import *
+66
View File
@@ -0,0 +1,66 @@
import numpy as np
cimport numpy as np
cimport cython
cdef extern from "math.h":
int abs(int i)
@cython.boundscheck(False)
@cython.wraparound(False)
def bresenham(int y, int x, int y2, int x2):
"""
Generate line pixel coordinates.
Parameters
----------
y, x : int
Starting position (row, column).
y2, x2 : int
End position (row, column).
Returns
-------
rr, cc : (N,) ndarray of int
Indices of pixels that belong to the line.
May be used to directly index into an array, e.g.
``img[rr, cc] = 1``.
"""
cdef np.ndarray[np.int32_t, ndim=1, mode="c"] rr, cc
cdef int steep = 0
cdef int dx = abs(x2 - x)
cdef int dy = abs(y2 - y)
cdef int sx, sy, d, i
if (x2 - x) > 0: sx = 1
else: sx = -1
if (y2 - y) > 0: sy = 1
else: sy = -1
if dy > dx:
steep = 1
x,y = y,x
dx,dy = dy,dx
sx,sy = sy,sx
d = (2 * dy) - dx
rr = np.zeros(int(dx) + 1, dtype=np.int32)
cc = np.zeros(int(dx) + 1, dtype=np.int32)
for i in range(dx):
if steep:
rr[i] = x
cc[i] = y
else:
rr[i] = y
cc[i] = x
while d >= 0:
y = y + sy
d = d - (2 * dx)
x = x + sx
d = d + (2 * dy)
rr[dx] = y2
cc[dx] = x2
return rr, cc
+6
View File
@@ -0,0 +1,6 @@
"""
Methods to draw on arrays.
"""
from _draw import bresenham
+30
View File
@@ -0,0 +1,30 @@
#!/usr/bin/env python
import os
from scikits.image._build import cython
base_path = os.path.abspath(os.path.dirname(__file__))
def configuration(parent_package='', top_path=None):
from numpy.distutils.misc_util import Configuration, get_numpy_include_dirs
config = Configuration('draw', parent_package, top_path)
config.add_data_dir('tests')
cython(['_draw.pyx'], working_path=base_path)
config.add_extension('_draw', sources=['_draw.c'],
include_dirs=[get_numpy_include_dirs()])
return config
if __name__ == '__main__':
from numpy.distutils.core import setup
setup(maintainer = 'Scikits-image developers',
author = 'Scikits-image developers',
maintainer_email = 'scikits-image@googlegroups.com',
description = 'Drawing',
url = 'https://github.com/scikits-image/scikits.image',
license = 'SciPy License (BSD Style)',
**(configuration(top_path='').todict())
)
+52
View File
@@ -0,0 +1,52 @@
from numpy.testing import assert_array_equal
import numpy as np
from scikits.image.draw import bresenham
def test_bresenham_horizontal():
img = np.zeros((10, 10))
rr, cc = bresenham(0, 0, 0, 9)
img[rr, cc] = 1
img_ = np.zeros((10, 10))
img_[0, :] = 1
assert_array_equal(img, img_)
def test_bresenham_vertical():
img = np.zeros((10, 10))
rr, cc = bresenham(0, 0, 9, 0)
img[rr, cc] = 1
img_ = np.zeros((10, 10))
img_[:, 0] = 1
assert_array_equal(img, img_)
def test_reverse():
img = np.zeros((10, 10))
rr, cc = bresenham(0, 9, 0, 0)
img[rr, cc] = 1
img_ = np.zeros((10, 10))
img_[0, :] = 1
assert_array_equal(img, img_)
def test_diag():
img = np.zeros((5, 5))
rr, cc = bresenham(0, 0, 4, 4)
img[rr, cc] = 1
img_ = np.eye(5)
assert_array_equal(img, img_)
if __name__ == "__main__":
from numpy.testing import run_module_suite
+70 -87
View File
@@ -1,51 +1,43 @@
"""
:author: Brian Holt, 2011
:license: modified BSD
"""
import numpy as np
from scipy import sqrt, pi, arctan2, cos, sin
from scipy.ndimage import uniform_filter
# XXX Replace with integral after merge
from ..transform import sat_sum
def hog(image, n_orientations=9, pixels_per_cell=(8, 8),
def hog(image, orientations=9, pixels_per_cell=(8, 8),
cells_per_block=(3, 3), visualise=False, normalise=False):
""" Extract Histogram of Oriented Gradients (HOG) for a given image.
"""Extract Histogram of Oriented Gradients (HOG) for a given image.
Compute a Histogram of Oriented Gradients (HOG) by
1) (optional) global image normalisation
2) computing the gradient image in x and y
3) computing gradient histograms
3) normalise across blocks
4) flatten into a feature vector
3) normalising across blocks
4) flattening into a feature vector
Parameters
----------
image: ndarray, 2D
2D image (greyscale)
n_orientations : int
number of orientation bins
pixels_per_cell : 2 tuple (int,int)
pixels per cell, size in pixels of a cell
image : (M, N) ndarray
Input image (greyscale).
orientations : int
Number of orientation bins.
pixels_per_cell : 2 tuple (int, int)
Size (in pixels) of a cell.
cells_per_block : 2 tuple (int,int)
cells per block, number of cells in each block
Number of cells in each block.
visualise : bool, optional
return an image of the HOG
Also return an image of the HOG.
normalise : bool, optional
apply power law compression to normalise the image before
processing
Apply power law compression to normalise the image before
processing.
Returns
-------
newarr : ndarray
HOG for the image as a 1D (flattened) array.
hog_image : PIL Image (if visualise=True)
A visualisation of the HOG image
hog_image : ndarray (if visualise=True)
A visualisation of the HOG image.
References
----------
@@ -54,8 +46,8 @@ def hog(image, n_orientations=9, pixels_per_cell=(8, 8),
* Dalal, N and Triggs, B, Histograms of Oriented Gradients for
Human Detection, IEEE Computer Society Conference on Computer
Vision and Pattern Recognition 2005 San Diego, CA, USA
"""
"""
image = np.atleast_2d(image)
"""
@@ -68,9 +60,9 @@ def hog(image, n_orientations=9, pixels_per_cell=(8, 8),
shadowing and illumination variations.
"""
if image.ndim == 3:
# replace RGB with locally dominant colour channel
pass # TODO
if image.ndim > 3:
raise ValueError("Currently only supports grey-level images")
if normalise:
image = sqrt(image)
@@ -91,40 +83,22 @@ def hog(image, n_orientations=9, pixels_per_cell=(8, 8),
"""
The third stage aims to produce an encoding that is sensitive to
local image content while remaining resistant to small changes in pose
or appearance. The adopted method pools gradient orientation information
locally in the same way as the SIFT [Lowe 2004] feature. The image window
is divided into small spatial regions, called "cells". For each cell we
accumulate a local 1-D histogram of gradient or edge orientations over
all the pixels in the cell. This combined cell-level 1-D histogram
forms the basic "orientation histogram" representation. Each orientation
histogram divides the gradient angle range into a fixed number of
predetermined bins. The gradient magnitudes of the pixels in the cell
are used to vote into the orientation histogram.
local image content while remaining resistant to small changes in
pose or appearance. The adopted method pools gradient orientation
information locally in the same way as the SIFT [Lowe 2004]
feature. The image window is divided into small spatial regions,
called "cells". For each cell we accumulate a local 1-D histogram
of gradient or edge orientations over all the pixels in the
cell. This combined cell-level 1-D histogram forms the basic
"orientation histogram" representation. Each orientation histogram
divides the gradient angle range into a fixed number of
predetermined bins. The gradient magnitudes of the pixels in the
cell are used to vote into the orientation histogram.
"""
magnitude = sqrt(gx ** 2 + gy ** 2)
orientation = arctan2(gy, (gx + 1e-15)) * (180 / pi) + 90
# compute n_orientations integral images
integral_images = []
for i in range(0, n_orientations):
#create new integral image for this orientation
# isolate orientations in this range
temp_ori = np.where(orientation < 180 / n_orientations * (i + 1),
orientation, 0)
temp_ori = np.where(orientation >= 180 / n_orientations * i,
temp_ori, 0)
# select magnitudes for those orientations
cond2 = temp_ori > 0
temp_mag = np.where(cond2, magnitude, 0)
#compute integral image
integral = np.cumsum(np.cumsum(temp_mag, axis=0, dtype=float),
axis=1, dtype=float)
integral_images.append(integral)
sx, sy = image.shape
cx, cy = pixels_per_cell
bx, by = cells_per_block
@@ -132,34 +106,43 @@ def hog(image, n_orientations=9, pixels_per_cell=(8, 8),
n_cellsx = int(np.floor(sx // cx)) # number of cells in x
n_cellsy = int(np.floor(sy // cy)) # number of cells in y
# compute orientations integral images
orientation_histogram = np.zeros((n_cellsx, n_cellsy, orientations))
for i in range(orientations):
#create new integral image for this orientation
# isolate orientations in this range
temp_ori = np.where(orientation < 180 / orientations * (i + 1),
orientation, 0)
temp_ori = np.where(orientation >= 180 / orientations * i,
temp_ori, 0)
# select magnitudes for those orientations
cond2 = temp_ori > 0
temp_mag = np.where(cond2, magnitude, 0)
orientation_histogram[:,:,i] = uniform_filter(temp_mag, size=(cx, cy))[cx/2::cx, cy/2::cy].T
# now for each cell, compute the histogram
orientation_histogram = np.zeros((n_cellsx, n_cellsy, n_orientations))
#orientation_histogram = np.zeros((n_cellsx, n_cellsy, orientations))
radius = min(cx, cy) // 2 - 1
hog_image = None
if visualise:
import Image
import ImageDraw
hog_image = Image.new("F", (sy, sx))
draw = ImageDraw.Draw(hog_image)
hog_image = np.zeros((sy, sx), dtype=float)
for x in range(0, n_cellsx):
for y in range(0, n_cellsy):
for o in range(0, n_orientations):
# compute the histogram from integral image
orientation_histogram[x, y, o] = sat_sum(integral_images[o],
y * cy,
x * cx,
(y + 1) * cy - 1,
(x + 1) * cx - 1)
if visualise:
if visualise:
from scikits.image import draw
for x in range(n_cellsx):
for y in range(n_cellsy):
for o in range(orientations):
centre = tuple([y * cy + cy // 2, x * cx + cx // 2])
dx = radius * cos(float(o) / n_orientations * np.pi)
dy = radius * sin(float(o) / n_orientations * np.pi)
draw.line([(centre[0] - dx, centre[1] - dy),
(centre[0] + dx, centre[1] + dy)],
fill=orientation_histogram[x, y, o])
dx = radius * cos(float(o) / orientations * np.pi)
dy = radius * sin(float(o) / orientations * np.pi)
rr, cc = draw.bresenham(centre[0] - dx, centre[1] - dy,
centre[0] + dx, centre[1] + dy)
hog_image[rr, cc] += orientation_histogram[x, y, o]
"""
The fourth stage computes normalisation, which takes local groups of
@@ -179,18 +162,18 @@ def hog(image, n_orientations=9, pixels_per_cell=(8, 8),
n_blocksx = (n_cellsx - bx) + 1
n_blocksy = (n_cellsy - by) + 1
normalised_blocks = np.zeros((n_blocksx, n_blocksy,
bx, by, n_orientations))
bx, by, orientations))
for x in range(0, n_blocksx):
for y in range(0, n_blocksy):
for x in range(n_blocksx):
for y in range(n_blocksy):
block = orientation_histogram[x:x + bx, y:y + by, :]
eps = 1e-5
normalised_blocks[x, y, :] = block / sqrt(block.sum() ** 2 + eps)
"""
The final step collects the HOG descriptors from all blocks of a dense
overlapping grid of blocks covering the detection window into a combined
feature vector for use in the window classifier
The final step collects the HOG descriptors from all blocks of a dense
overlapping grid of blocks covering the detection window into a combined
feature vector for use in the window classifier.
"""
if visualise:
+4 -2
View File
@@ -4,10 +4,12 @@ import scipy
from scikits.image.feature import hog
def test_histogram_of_oriented_gradients():
img = scipy.lena().astype(np.int8)
# Replace with scikits.image.data.lena() after merge
img = scipy.misc.lena().astype(np.int8)
fd = hog(img, n_orientations=9, pixels_per_cell=(8, 8),
fd = hog(img, orientations=9, pixels_per_cell=(8, 8),
cells_per_block=(1, 1))
assert len(fd) == 9 * (512//8) ** 2
if __name__ == '__main__':