mirror of
https://github.com/wassname/scikit-image.git
synced 2026-07-19 11:27:45 +08:00
Rebased to master, incorporated view_as_blocks, cleaned up interpolate method and example
This commit is contained in:
Steven Silvester
@@ -123,10 +123,6 @@
|
||||
|
||||
- Luis Pedro Coelho
|
||||
imread plugin
|
||||
=======
|
||||
Polygon, circle and ellipse drawing functions
|
||||
Adaptive thresholding
|
||||
Implementation of Matlab's `regionprops`
|
||||
|
||||
- Steven Silvester, Karel Zuiderveld
|
||||
Adaptive Histogram Equalization
|
||||
|
||||
@@ -18,7 +18,7 @@ that fall within the 2nd and 98th percentiles [2]_.
|
||||
|
||||
"""
|
||||
|
||||
from skimage import data
|
||||
from skimage import data, img_as_ubyte
|
||||
from skimage.util.dtype import dtype_range
|
||||
from skimage import exposure
|
||||
|
||||
@@ -39,7 +39,7 @@ def plot_img_and_hist(img, axes, bins=256):
|
||||
ax_img.set_axis_off()
|
||||
|
||||
# Display histogram
|
||||
ax_hist.hist(img.ravel(), bins=bins)
|
||||
ax_hist.hist(img.ravel(), bins=bins, histtype='stepfilled')
|
||||
ax_hist.ticklabel_format(axis='y', style='scientific', scilimits=(0, 0))
|
||||
ax_hist.set_xlabel('Pixel intensity')
|
||||
|
||||
@@ -63,10 +63,11 @@ img_rescale = exposure.rescale_intensity(img, in_range=(p2, p98))
|
||||
|
||||
# Equalization
|
||||
img_eq = exposure.equalize(img)
|
||||
img_eq = img_as_ubyte(img_eq)
|
||||
|
||||
# Adaptive Equalization
|
||||
img_adapteq = exposure.adapthist(img, clip_limit=0.03)
|
||||
|
||||
img_adapteq = img_as_ubyte(img_adapteq)
|
||||
|
||||
# Display results
|
||||
f, axes = plt.subplots(2, 4, figsize=(8, 4))
|
||||
|
||||
@@ -25,6 +25,7 @@ import numpy as np
|
||||
import skimage
|
||||
from skimage import color
|
||||
from skimage.exposure import rescale_intensity
|
||||
from skimage.util import view_as_blocks
|
||||
|
||||
|
||||
MAX_REG_X = 16 # max. # contextual regions in x-direction */
|
||||
@@ -58,8 +59,6 @@ def adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, nbins=256):
|
||||
* The algorithm relies on an image whose rows and columns are even
|
||||
multiples of the number of tiles, so the extra rows and columns are left
|
||||
at their original values, thus preserving the input image shape.
|
||||
* For grayscale images, CLAHE is performed on one channel,
|
||||
and a grayscale is returned
|
||||
* For color images, the following steps are performed:
|
||||
- The image is converted to LAB color space
|
||||
- The CLAHE algorithm is run on the L channel
|
||||
@@ -127,7 +126,7 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128):
|
||||
if clip_limit == 1.0:
|
||||
return image # is OK, immediately returns original image.
|
||||
|
||||
map_array = np.zeros((ntiles_x * ntiles_y, nbins), dtype=int)
|
||||
map_array = np.zeros((ntiles_y, ntiles_x, nbins), dtype=int)
|
||||
|
||||
y_res = image.shape[0] - image.shape[0] % ntiles_y
|
||||
x_res = image.shape[1] - image.shape[1] % ntiles_x
|
||||
@@ -146,21 +145,17 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128):
|
||||
bin_size = 1 + NR_OF_GREY / nbins
|
||||
aLUT = np.arange(NR_OF_GREY)
|
||||
aLUT /= bin_size
|
||||
img_blocks = view_as_blocks(image, (y_size, x_size))
|
||||
# Calculate greylevel mappings for each contextual region
|
||||
ystart = 0
|
||||
for y in range(ntiles_y):
|
||||
xstart = 0
|
||||
for x in range(ntiles_x):
|
||||
sub_img = image[ystart: ystart + y_size,
|
||||
xstart: xstart + x_size]
|
||||
sub_img = img_blocks[y, x]
|
||||
hist = aLUT[sub_img.ravel()]
|
||||
hist = np.bincount(hist)
|
||||
hist = np.append(hist, np.zeros(nbins - hist.size, dtype=int))
|
||||
hist = clip_histogram(hist, clip_limit)
|
||||
hist = map_histogram(hist, 0, NR_OF_GREY, n_pixels)
|
||||
map_array[y * ntiles_x + x] = hist
|
||||
xstart += x_size
|
||||
ystart += y_size
|
||||
map_array[y, x] = hist
|
||||
# Interpolate greylevel mappings to get CLAHE image
|
||||
ystart = 0
|
||||
for y in range(ntiles_y + 1):
|
||||
@@ -190,11 +185,13 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128):
|
||||
xstep = x_size
|
||||
xL = x - 1
|
||||
xR = xL + 1
|
||||
mapLU = map_array[yU * ntiles_x + xL]
|
||||
mapRU = map_array[yU * ntiles_x + xR]
|
||||
mapLB = map_array[yB * ntiles_x + xL]
|
||||
mapRB = map_array[yB * ntiles_x + xR]
|
||||
interpolate(image, xstart, xstep, ystart, ystep,
|
||||
mapLU = map_array[yU, xL]
|
||||
mapRU = map_array[yU, xR]
|
||||
mapLB = map_array[yB, xL]
|
||||
mapRB = map_array[yB, xR]
|
||||
xslice = np.arange(xstart, xstart + xstep)
|
||||
yslice = np.arange(ystart, ystart + ystep)
|
||||
interpolate(image, xslice, yslice,
|
||||
mapLU, mapRU, mapLB, mapRB, aLUT)
|
||||
xstart += xstep # set pointer on next matrix */
|
||||
ystart += ystep
|
||||
@@ -280,7 +277,7 @@ def map_histogram(hist, min_val, max_val, n_pixels):
|
||||
return out.astype(int)
|
||||
|
||||
|
||||
def interpolate(image, xstart, xstep, ystart, ystep,
|
||||
def interpolate(image, xslice, yslice,
|
||||
mapLU, mapRU, mapLB, mapRB, aLUT):
|
||||
'''Find the new grayscale level for a region using bilinear interpolation
|
||||
|
||||
@@ -288,10 +285,8 @@ def interpolate(image, xstart, xstep, ystart, ystep,
|
||||
----------
|
||||
image : np.ndarray
|
||||
full image
|
||||
xstart, xstop : int
|
||||
indices of xslice
|
||||
ystart, ystop : int
|
||||
indices of yslice
|
||||
xslice, yslice : array-like
|
||||
indices of the region
|
||||
map* : np.ndarray
|
||||
mappings of greylevels from histograms
|
||||
aLUT : np.ndarray
|
||||
@@ -309,18 +304,18 @@ def interpolate(image, xstart, xstep, ystart, ystep,
|
||||
This is done by a bilinear interpolation between four different
|
||||
mappings in order to eliminate boundary artifacts.
|
||||
'''
|
||||
norm = xstep * ystep # Normalization factor
|
||||
|
||||
norm = xslice.size * yslice.size # Normalization factor
|
||||
# interpolation weight matrices
|
||||
x_coef, y_coef = np.meshgrid(np.arange(xstep),
|
||||
np.arange(ystep))
|
||||
x_coef, y_coef = np.meshgrid(np.arange(xslice.size),
|
||||
np.arange(yslice.size))
|
||||
x_inv_coef, y_inv_coef = x_coef[:, ::-1] + 1, y_coef[::-1] + 1
|
||||
|
||||
im_slice = image[ystart: ystart + ystep, xstart: xstart + xstep]
|
||||
im_slice = aLUT[im_slice]
|
||||
view = image[yslice[0]: yslice[-1] + 1, xslice[0]: xslice[-1] + 1]
|
||||
im_slice = aLUT[view]
|
||||
new = ((y_inv_coef * (x_inv_coef * mapLU[im_slice]
|
||||
+ x_coef * mapRU[im_slice])
|
||||
+ y_coef * (x_inv_coef * mapLB[im_slice]
|
||||
+ x_coef * mapRB[im_slice]))
|
||||
/ norm)
|
||||
image[ystart: ystart + ystep, xstart: xstart + xstep] = new
|
||||
view[:, :] = new
|
||||
return image
|
||||
|
||||
@@ -1,11 +1,13 @@
|
||||
import numpy as np
|
||||
from numpy.testing import assert_array_almost_equal as assert_close
|
||||
import skimage
|
||||
from skimage import io
|
||||
from skimage import data
|
||||
from skimage import exposure
|
||||
from skimage.color import rgb2gray
|
||||
from skimage.util.dtype import dtype_range
|
||||
|
||||
io.use_plugin('qt')
|
||||
|
||||
# Test histogram equalization
|
||||
# ===========================
|
||||
@@ -81,7 +83,7 @@ def test_adapthist_scalar():
|
||||
img = skimage.img_as_ubyte(data.moon())
|
||||
adapted = exposure.adapthist(img, clip_limit=0.02)
|
||||
assert adapted.min() == 0
|
||||
assert adapted.max() == (1<<16) - 1
|
||||
assert adapted.max() == (1 << 16) - 1
|
||||
assert img.shape == adapted.shape
|
||||
full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img))
|
||||
assert_almost_equal = np.testing.assert_almost_equal
|
||||
@@ -90,6 +92,7 @@ def test_adapthist_scalar():
|
||||
0.0410010)
|
||||
return img, adapted
|
||||
|
||||
|
||||
def test_adapthist_grayscale():
|
||||
'''Test a grayscale float image
|
||||
'''
|
||||
|
||||
Reference in New Issue
Block a user