PageRenderTime 42ms CodeModel.GetById 17ms RepoModel.GetById 0ms app.codeStats 0ms

/skimage/feature/_texture.pyx

https://github.com/nitishb/scikit-image
Cython | 182 lines | 136 code | 31 blank | 15 comment | 47 complexity | 42fc20030018d8b99ec6e8c1b1c609f7 MD5 | raw file
Possible License(s): BSD-3-Clause
  1. #cython: cdivision=True
  2. #cython: boundscheck=False
  3. #cython: nonecheck=False
  4. #cython: wraparound=False
  5. import numpy as np
  6. cimport numpy as np
  7. from libc.math cimport sin, cos, abs
  8. from skimage._shared.interpolation cimport bilinear_interpolation
  9. def _glcm_loop(np.ndarray[dtype=np.uint8_t, ndim=2,
  10. negative_indices=False, mode='c'] image,
  11. np.ndarray[dtype=np.float64_t, ndim=1,
  12. negative_indices=False, mode='c'] distances,
  13. np.ndarray[dtype=np.float64_t, ndim=1,
  14. negative_indices=False, mode='c'] angles,
  15. int levels,
  16. np.ndarray[dtype=np.uint32_t, ndim=4,
  17. negative_indices=False, mode='c'] out
  18. ):
  19. """Perform co-occurrence matrix accumulation.
  20. Parameters
  21. ----------
  22. image : ndarray
  23. Input image, which is converted to the uint8 data type.
  24. distances : ndarray
  25. List of pixel pair distance offsets.
  26. angles : ndarray
  27. List of pixel pair angles in radians.
  28. levels : int
  29. The input image should contain integers in [0, levels-1],
  30. where levels indicate the number of grey-levels counted
  31. (typically 256 for an 8-bit image)
  32. out : ndarray
  33. On input a 4D array of zeros, and on output it contains
  34. the results of the GLCM computation.
  35. """
  36. cdef:
  37. np.int32_t a_inx, d_idx
  38. np.int32_t r, c, rows, cols, row, col
  39. np.int32_t i, j
  40. rows = image.shape[0]
  41. cols = image.shape[1]
  42. for a_idx, angle in enumerate(angles):
  43. for d_idx, distance in enumerate(distances):
  44. for r in range(rows):
  45. for c in range(cols):
  46. i = image[r, c]
  47. # compute the location of the offset pixel
  48. row = r + <int>(sin(angle) * distance + 0.5)
  49. col = c + <int>(cos(angle) * distance + 0.5);
  50. # make sure the offset is within bounds
  51. if row >= 0 and row < rows and \
  52. col >= 0 and col < cols:
  53. j = image[row, col]
  54. if i >= 0 and i < levels and \
  55. j >= 0 and j < levels:
  56. out[i, j, d_idx, a_idx] += 1
  57. cdef inline int _bit_rotate_right(int value, int length):
  58. """Cyclic bit shift to the right.
  59. Parameters
  60. ----------
  61. value : int
  62. integer value to shift
  63. length : int
  64. number of bits of integer
  65. """
  66. return (value >> 1) | ((value & 1) << (length - 1))
  67. def _local_binary_pattern(np.ndarray[double, ndim=2] image,
  68. int P, float R, char method='D'):
  69. """Gray scale and rotation invariant LBP (Local Binary Patterns).
  70. LBP is an invariant descriptor that can be used for texture classification.
  71. Parameters
  72. ----------
  73. image : (N, M) double array
  74. Graylevel image.
  75. P : int
  76. Number of circularly symmetric neighbour set points (quantization of the
  77. angular space).
  78. R : float
  79. Radius of circle (spatial resolution of the operator).
  80. method : {'D', 'R', 'U', 'V'}
  81. Method to determine the pattern.
  82. * 'D': 'default'
  83. * 'R': 'ror'
  84. * 'U': 'uniform'
  85. * 'V': 'var'
  86. Returns
  87. -------
  88. output : (N, M) array
  89. LBP image.
  90. """
  91. # texture weights
  92. cdef np.ndarray[int, ndim=1] weights = 2 ** np.arange(P, dtype=np.int32)
  93. # local position of texture elements
  94. rp = - R * np.sin(2 * np.pi * np.arange(P, dtype=np.double) / P)
  95. cp = R * np.cos(2 * np.pi * np.arange(P, dtype=np.double) / P)
  96. cdef np.ndarray[double, ndim=2] coords = np.round(np.vstack([rp, cp]).T, 5)
  97. # pre allocate arrays for computation
  98. cdef np.ndarray[double, ndim=1] texture = np.zeros(P, np.double)
  99. cdef np.ndarray[char, ndim=1] signed_texture = np.zeros(P, np.int8)
  100. cdef np.ndarray[int, ndim=1] rotation_chain = np.zeros(P, np.int32)
  101. output_shape = (image.shape[0], image.shape[1])
  102. cdef np.ndarray[double, ndim=2] output = np.zeros(output_shape, np.double)
  103. cdef int rows = image.shape[0]
  104. cdef int cols = image.shape[1]
  105. cdef double lbp
  106. cdef int r, c, changes, i
  107. for r in range(image.shape[0]):
  108. for c in range(image.shape[1]):
  109. for i in range(P):
  110. texture[i] = bilinear_interpolation(<double*>image.data,
  111. rows, cols, r + coords[i, 0], c + coords[i, 1], 'C', 0)
  112. # signed / thresholded texture
  113. for i in range(P):
  114. if texture[i] - image[r, c] >= 0:
  115. signed_texture[i] = 1
  116. else:
  117. signed_texture[i] = 0
  118. lbp = 0
  119. # if method == 'uniform' or method == 'var':
  120. if method == 'U' or method == 'V':
  121. # determine number of 0 - 1 changes
  122. changes = 0
  123. for i in range(P - 1):
  124. changes += abs(signed_texture[i] - signed_texture[i + 1])
  125. if changes <= 2:
  126. for i in range(P):
  127. lbp += signed_texture[i]
  128. else:
  129. lbp = P + 1
  130. if method == 'V':
  131. var = np.var(texture)
  132. if var != 0:
  133. lbp /= var
  134. else:
  135. lbp = np.nan
  136. else:
  137. # method == 'default'
  138. for i in range(P):
  139. lbp += signed_texture[i] * weights[i]
  140. # method == 'ror'
  141. if method == 'R':
  142. # shift LBP P times to the right and get minimum value
  143. rotation_chain[0] = <int>lbp
  144. for i in range(1, P):
  145. rotation_chain[i] = \
  146. _bit_rotate_right(rotation_chain[i - 1], P)
  147. lbp = rotation_chain[0]
  148. for i in range(1, P):
  149. lbp = min(lbp, rotation_chain[i])
  150. output[r, c] = lbp
  151. return output