Implementing Otsu binarization for faded document images
Taken from the Otsu method on Wikipedia
I = imread('cameraman.tif');
Step 1. Calculate the histogram and probabilities for each intensity level.
nbins = 256; % Number of bins
counts = imhist(I,nbins); % Each intensity increments the histogram from 0 to 255
p = counts / sum(counts); % Probabilities
Step 2. Set up initial omega_i (0) and mu_i (0)
omega1 = 0;
omega2 = 1;
mu1 = 0;
mu2 = mean(I(:));
Step 3. Step through all possible thresholds from 0 to maximum intensity (255)
Step 3.1 Update omega_i and mu_i
Step 3.2 Calculate sigma_b_squared
for t = 1:nbins
omega1(t) = sum(p(1:t));
omega2(t) = sum(p(t+1:end));
mu1(t) = sum(p(1:t).*(1:t)');
mu2(t) = sum(p(t+1:end).*(t+1:nbins)');
end
sigma_b_squared_wiki = omega1 .* omega2 .* (mu2-mu1).^2; % Eq. (14)
sigma_b_squared_otsu = (mu1(end) .* omega1-mu1) .^2 ./(omega1 .* (1-omega1)); % Eq. (18)
Step 4 The desired threshold matches the location of the sigma_b_squared maximum
[~,thres_level_wiki] = max(sigma_b_squared_wiki);
[~,thres_level_otsu] = max(sigma_b_squared_otsu);
There are some differences between the wiki version of eq. (14) in Otsu and eq. (18) and I don't know why. But thres_level_otsu
match MATLAB implementationgraythresh(I)
source to share
Since the function graythresh
in Matlab implements the Otsu method, you need to convert the image to grayscale and then use the function im2bw
to binarize the image using the threhsold level returned graythresh
.
To convert your image I
to grayscale, you can use the following code:
I = im2uint8(I);
if size(I,3) ~= 1
I = rgb2gray(I);
end;
To get a binary image Ib
using Otsu method use the following code:
Ib = im2bw(I, graythresh(I));
You should get the following output:
source to share
Starting from what your original question was doing OTSU, rearranging it true that the MATLAB function is graythresh
based on this method The OTSU method treats the threshold as a valley between two peaks, which is one of the foreground pixels and the other of the background pixels
Regarding your image, which looks like a historical manuscript, found this document comparing all the methods that can be used for document threshold images
You can also download and read the sauvola threshold from here
Good luck with its implementation =)
source to share
Fixed MATLAB implementation (for 2d matrix)
function [T] = myotsu(I,N); % create histogram nbins = N; [x,h] = hist(I(:),nbins); % calculate probabilities p = x./sum(x); % initialisation om1 = 0; om2 = 1; mu1 = 0; mu2 = mode(I(:)); for t = 1:nbins, om1(t) = sum(p(1:t)); om2(t) = sum(p(t+1:nbins)); mu1(t) = sum(p(1:t).*[1:t]); mu2(t) = sum(p(t+1:nbins).*[t+1:nbins]); end sigma = (mu1(nbins).*om1-mu1).^2./(om1.*(1-om1)); idx = find(sigma == max(sigma)); T = h(idx(1));
source to share