From 4a9534548956c80692f745b807041f138be08e97 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Johannes=20Sch=C3=B6nberger?= Date: Wed, 20 May 2015 17:51:41 -0700 Subject: [PATCH] Do not acquire GIL for texture functions --- skimage/feature/_texture.pyx | 279 +++++++++++++------------- skimage/feature/tests/test_texture.py | 3 + 2 files changed, 146 insertions(+), 136 deletions(-) diff --git a/skimage/feature/_texture.pyx b/skimage/feature/_texture.pyx index 5eb91486..d454a812 100644 --- a/skimage/feature/_texture.pyx +++ b/skimage/feature/_texture.pyx @@ -4,7 +4,7 @@ #cython: wraparound=False import numpy as np cimport numpy as cnp -from libc.math cimport sin, cos, abs +from libc.math cimport sin, cos, abs, NAN from .._shared.interpolation cimport bilinear_interpolation, round @@ -36,32 +36,33 @@ def _glcm_loop(cnp.uint8_t[:, ::1] image, double[:] distances, cnp.uint8_t i, j cnp.float64_t angle, distance - rows = image.shape[0] - cols = image.shape[1] + with nogil: + rows = image.shape[0] + cols = image.shape[1] - for a_idx in range(len(angles)): - angle = angles[a_idx] - for d_idx in range(len(distances)): - distance = distances[d_idx] - for r in range(rows): - for c in range(cols): - i = image[r, c] + for a_idx in range(angles.shape[0]): + angle = angles[a_idx] + for d_idx in range(distances.shape[0]): + distance = distances[d_idx] + for r in range(rows): + for c in range(cols): + i = image[r, c] - # compute the location of the offset pixel - row = r + round(sin(angle) * distance) - col = c + round(cos(angle) * distance) + # compute the location of the offset pixel + row = r + round(sin(angle) * distance) + col = c + round(cos(angle) * distance) - # make sure the offset is within bounds - if row >= 0 and row < rows and \ - col >= 0 and col < cols: - j = image[row, col] + # make sure the offset is within bounds + if row >= 0 and row < rows and \ + col >= 0 and col < cols: + j = image[row, col] - if i >= 0 and i < levels and \ - j >= 0 and j < levels: - out[i, j, d_idx, a_idx] += 1 + if i >= 0 and i < levels and \ + j >= 0 and j < levels: + out[i, j, d_idx, a_idx] += 1 -cdef inline int _bit_rotate_right(int value, int length): +cdef inline int _bit_rotate_right(int value, int length) nogil: """Cyclic bit shift to the right. Parameters @@ -132,124 +133,130 @@ def _local_binary_pattern(double[:, ::1] image, # To compute the variance features cdef double sum_, var_, texture_i - for r in range(image.shape[0]): - for c in range(image.shape[1]): - for i in range(P): - texture[i] = bilinear_interpolation(&image[0, 0], rows, cols, - r + rp[i], c + cp[i], - 'C', 0) - # signed / thresholded texture - for i in range(P): - if texture[i] - image[r, c] >= 0: - signed_texture[i] = 1 - else: - signed_texture[i] = 0 - - lbp = 0 - - # if method == 'var': - if method == 'V': - # Compute the variance without passing from numpy. - # Following the LBP paper, we're taking a biased estimate - # of the variance (ddof=0) - sum_ = 0.0 - var_ = 0.0 + with nogil: + for r in range(image.shape[0]): + for c in range(image.shape[1]): for i in range(P): - texture_i = texture[i] - sum_ += texture_i - var_ += texture_i * texture_i - var_ = (var_ - (sum_ * sum_) / P) / P - if var_ != 0: - lbp = var_ - else: - lbp = np.nan - # if method == 'uniform': - elif method == 'U' or method == 'N': - # determine number of 0 - 1 changes - changes = 0 - for i in range(P - 1): - changes += abs(signed_texture[i] - signed_texture[i + 1]) - if method == 'N': - # Uniform local binary patterns are defined as patterns - # with at most 2 value changes (from 0 to 1 or from 1 to - # 0). Uniform patterns can be caraterized by their number - # `n_ones` of 1. The possible values for `n_ones` range - # from 0 to P. - # Here is an example for P = 4: - # n_ones=0: 0000 - # n_ones=1: 0001, 1000, 0100, 0010 - # n_ones=2: 0011, 1001, 1100, 0110 - # n_ones=3: 0111, 1011, 1101, 1110 - # n_ones=4: 1111 - # - # For a pattern of size P there are 2 constant patterns - # corresponding to n_ones=0 and n_ones=P. For each other - # value of `n_ones` , i.e n_ones=[1..P-1], there are P - # possible patterns which are related to each other through - # circular permutations. The total number of uniform - # patterns is thus (2 + P * (P - 1)). - # Given any pattern (uniform or not) we must be able to - # associate a unique code: - # 1. Constant patterns patterns (with n_ones=0 and - # n_ones=P) and non uniform patterns are given fixed - # code values. - # 2. Other uniform patterns are indexed considering the - # value of n_ones, and an index called 'rot_index' - # reprenting the number of circular right shifts - # required to obtain the pattern starting from a - # reference position (corresponding to all zeros stacked - # on the right). This number of rotations (or circular - # right shifts) 'rot_index' is efficiently computed by - # considering the positions of the first 1 and the first - # 0 found in the pattern. - - if changes <= 2: - # We have a uniform pattern - n_ones = 0 # determines the number of ones - first_one = -1 # position was the first one - first_zero = -1 # position of the first zero - for i in range(P): - if signed_texture[i]: - n_ones += 1 - if first_one == -1: - first_one = i - else: - if first_zero == -1: - first_zero = i - if n_ones == 0: - lbp = 0 - elif n_ones == P: - lbp = P * (P - 1) + 1 - else: - if first_one == 0: - rot_index = n_ones - first_zero - else: - rot_index = P - first_one - lbp = 1 + (n_ones - 1) * P + rot_index - else: # changes > 2 - lbp = P * (P - 1) + 2 - else: # method != 'N' - if changes <= 2: - for i in range(P): - lbp += signed_texture[i] + texture[i] = bilinear_interpolation(&image[0, 0], rows, cols, + r + rp[i], c + cp[i], + 'C', 0) + # signed / thresholded texture + for i in range(P): + if texture[i] - image[r, c] >= 0: + signed_texture[i] = 1 else: - lbp = P + 1 - else: - # method == 'default' - for i in range(P): - lbp += signed_texture[i] * weights[i] + signed_texture[i] = 0 - # method == 'ror' - if method == 'R': - # shift LBP P times to the right and get minimum value - rotation_chain[0] = lbp - for i in range(1, P): - rotation_chain[i] = \ - _bit_rotate_right(rotation_chain[i - 1], P) - lbp = rotation_chain[0] - for i in range(1, P): - lbp = min(lbp, rotation_chain[i]) + lbp = 0 - output[r, c] = lbp + # if method == 'var': + if method == 'V': + # Compute the variance without passing from numpy. + # Following the LBP paper, we're taking a biased estimate + # of the variance (ddof=0) + sum_ = 0.0 + var_ = 0.0 + for i in range(P): + texture_i = texture[i] + sum_ += texture_i + var_ += texture_i * texture_i + var_ = (var_ - (sum_ * sum_) / P) / P + if var_ != 0: + lbp = var_ + else: + lbp = NAN + # if method == 'uniform': + elif method == 'U' or method == 'N': + # determine number of 0 - 1 changes + changes = 0 + for i in range(P - 1): + changes += (signed_texture[i] + - signed_texture[i + 1]) != 0 + if method == 'N': + # Uniform local binary patterns are defined as patterns + # with at most 2 value changes (from 0 to 1 or from 1 to + # 0). Uniform patterns can be characterized by their + # number `n_ones` of 1. The possible values for + # `n_ones` range from 0 to P. + # + # Here is an example for P = 4: + # n_ones=0: 0000 + # n_ones=1: 0001, 1000, 0100, 0010 + # n_ones=2: 0011, 1001, 1100, 0110 + # n_ones=3: 0111, 1011, 1101, 1110 + # n_ones=4: 1111 + # + # For a pattern of size P there are 2 constant patterns + # corresponding to n_ones=0 and n_ones=P. For each other + # value of `n_ones` , i.e n_ones=[1..P-1], there are P + # possible patterns which are related to each other + # through circular permutations. The total number of + # uniform patterns is thus (2 + P * (P - 1)). + + # Given any pattern (uniform or not) we must be able to + # associate a unique code: + # + # 1. Constant patterns patterns (with n_ones=0 and + # n_ones=P) and non uniform patterns are given fixed + # code values. + # + # 2. Other uniform patterns are indexed considering the + # value of n_ones, and an index called 'rot_index' + # reprenting the number of circular right shifts + # required to obtain the pattern starting from a + # reference position (corresponding to all zeros stacked + # on the right). This number of rotations (or circular + # right shifts) 'rot_index' is efficiently computed by + # considering the positions of the first 1 and the first + # 0 found in the pattern. + + if changes <= 2: + # We have a uniform pattern + n_ones = 0 # determines the number of ones + first_one = -1 # position was the first one + first_zero = -1 # position of the first zero + for i in range(P): + if signed_texture[i]: + n_ones += 1 + if first_one == -1: + first_one = i + else: + if first_zero == -1: + first_zero = i + if n_ones == 0: + lbp = 0 + elif n_ones == P: + lbp = P * (P - 1) + 1 + else: + if first_one == 0: + rot_index = n_ones - first_zero + else: + rot_index = P - first_one + lbp = 1 + (n_ones - 1) * P + rot_index + else: # changes > 2 + lbp = P * (P - 1) + 2 + else: # method != 'N' + if changes <= 2: + for i in range(P): + lbp += signed_texture[i] + else: + lbp = P + 1 + else: + # method == 'default' + for i in range(P): + lbp += signed_texture[i] * weights[i] + + # method == 'ror' + if method == 'R': + # shift LBP P times to the right and get minimum value + rotation_chain[0] = lbp + for i in range(1, P): + rotation_chain[i] = \ + _bit_rotate_right(rotation_chain[i - 1], P) + lbp = rotation_chain[0] + for i in range(1, P): + lbp = min(lbp, rotation_chain[i]) + + output[r, c] = lbp return np.asarray(output) diff --git a/skimage/feature/tests/test_texture.py b/skimage/feature/tests/test_texture.py index 788ccbc6..10444aa8 100644 --- a/skimage/feature/tests/test_texture.py +++ b/skimage/feature/tests/test_texture.py @@ -1,5 +1,6 @@ import numpy as np from skimage.feature import greycomatrix, greycoprops, local_binary_pattern +from skimage._shared.testing import test_parallel class TestGLCM(): @@ -10,6 +11,7 @@ class TestGLCM(): [0, 2, 2, 2], [2, 2, 3, 3]], dtype=np.uint8) + @test_parallel() def test_output_angles(self): result = greycomatrix(self.image, [1], [0, np.pi / 4, np.pi / 2, 3 * np.pi / 4], 4) assert result.shape == (4, 4, 1, 4) @@ -162,6 +164,7 @@ class TestLBP(): [ 0, 255, 30, 34, 255, 24], [146, 241, 255, 0, 189, 126]], dtype='double') + @test_parallel() def test_default(self): lbp = local_binary_pattern(self.image, 8, 1, 'default') ref = np.array([[ 0, 251, 0, 255, 96, 255],