-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathphansalkar.py
78 lines (69 loc) · 3.35 KB
/
phansalkar.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
import numpy as np
import skimage.transform
import skimage.filters
import itertools
from scipy import ndimage as ndi
from collections.abc import Iterable
from skimage import util
from math import e
# Fixed Parameters
q = 10 # q is a constant that determines how the exponential term affects the threshold value
R = .5 # R is the dynamic range of standard deviation (0.5 for normalized images)
# Slight modification of preexisting algorithm used for Niblack and Sauvola Thresholding, based of the scipy
# implementation. The real trick is using their version of _mean_std which efficiently computes the mean, std for
# each pixel group (sliding window). These values are then plugged into phansalkers equation and compared to the original
# to create a binary image. The values k and p are tunable parameters that change how the algorithm operate.
def pfilter(img, k=.25, p=3):
og_img = np.array(img)/255 # Normalize the images
m, s = _mean_std(og_img, 3)
return og_img > m * (1 + (p * (e ** (-q * m))) + (k * ((s / R) - 1)))
# Exact code as found in the filters package of scipy with the exception of removing the line related to valid window
# sizes, as I am not going to let the user specify window size.
def _mean_std(image, w):
"""Return local mean and standard deviation of each pixel using a
neighborhood defined by a rectangular window size ``w``.
The algorithm uses integral images to speedup computation. This is
used by :func:`threshold_niblack` and :func:`threshold_sauvola`.
Parameters
----------
image : ndarray
Input image.
w : int, or iterable of int
Window size specified as a single odd integer (3, 5, 7, …),
or an iterable of length ``image.ndim`` containing only odd
integers (e.g. ``(1, 5, 5)``).
Returns
-------
m : ndarray of float, same shape as ``image``
Local mean of the image.
s : ndarray of float, same shape as ``image``
Local standard deviation of the image.
References
----------
.. [1] F. Shafait, D. Keysers, and T. M. Breuel, "Efficient
implementation of local adaptive thresholding techniques
using integral images." in Document Recognition and
Retrieval XV, (San Jose, USA), Jan. 2008.
:DOI:`10.1117/12.767755`
"""
if not isinstance(w, Iterable):
w = (w,) * image.ndim
# _validate_window_size(w) My only edit!
pad_width = tuple((k // 2 + 1, k // 2) for k in w)
padded = np.pad(image.astype('float'), pad_width,
mode='reflect')
padded_sq = padded * padded
integral = skimage.transform.integral.integral_image(padded)
integral_sq = skimage.transform.integral.integral_image(padded_sq)
kern = np.zeros(tuple(k + 1 for k in w))
for indices in itertools.product(*([[0, -1]] * image.ndim)):
kern[indices] = (-1) ** (image.ndim % 2 != np.sum(indices) % 2)
total_window_size = np.prod(w)
sum_full = ndi.correlate(integral, kern, mode='constant')
m = skimage.util.crop(sum_full, pad_width) / total_window_size
sum_sq_full = ndi.correlate(integral_sq, kern, mode='constant')
g2 = skimage.util.crop(sum_sq_full, pad_width) / total_window_size
# Note: we use np.clip because g2 is not guaranteed to be greater than
# m*m when floating point error is considered
s = np.sqrt(np.clip(g2 - m * m, 0, None))
return m, s