在图像处理中Otsu方法，是以 Nobuyuki otsu 的名字命名的（日本人，大津展之），常用于基于图像分割的聚类。该算法的理论依据是：假定图像包含两类像素（前景像素和背景像素），直方图为双峰直方图，然后计算使得两类像素能分开的最佳阈值（类内方差），或等价的间类间方差最大。
Otsu算法原理：
对于图像 I(x,y)，前景(即目标)和背景的分割阈值记作 T，属于前景的像素点数占整幅图像的比例记为 ω0，平均灰度为 μ0；背景像素点数占整幅图像的比例为 ω1，平均灰度为 μ1；整幅图像的平均灰度记为μ，类间方差记为g。
假设图像大小为M×N，图像中像素的灰度值小于阈值 T 的像素个数为 N0，像素灰度大于阈值T的像素个数为 N1，那么：
ω0=N0/ M×N (1)
ω1=N1/ M×N (2)
N0+N1=M×N (3)
ω0+ω1=1 (4)
μ=ω0*μ0+ω1*μ1 (5)
g=ω0(μ0-μ)^2+ω1(μ1-μ)^2 (6)
g=ω0ω1(μ0-μ1)^2 (7)
采用遍历的方法使得类间方差g最大的阈值T，即为所求。Ostu方法可以形象地理解为：求取直方图有两个峰值的图像中那两个峰值之间的低谷值 T 。
Otsu算法实现：
matlab函数：
matlab中有现成的函数实现，函数名为： graythresh， 该函数便是用Ostu方法求分割阈值Ｔ。用法如下：
imgScr=imread('..');
T = graythresh(imgScr);
BW = im2bw(imgScr,T);我这里po下源码，方便理解：
function [level em] = graythresh(I)
%GRAYTHRESH Global image threshold using Otsu's method.
% LEVEL = GRAYTHRESH(I) computes a global threshold (LEVEL) that can be
% used to convert an intensity image to a binary image with IM2BW. LEVEL
% is a normalized intensity value that lies in the range [0, 1].
% GRAYTHRESH uses Otsu's method, which chooses the threshold to minimize
% the intraclass variance of the thresholded black and white pixels.
%
% [LEVEL EM] = GRAYTHRESH(I) returns effectiveness metric, EM, as the
% second output argument. It indicates the effectiveness of thresholding
% of the input image and it is in the range [0, 1]. The lower bound is
% attainable only by images having a single gray level, and the upper
% bound is attainable only by two-valued images.
%
% Class Support
% -------------
% The input image I can be uint8, uint16, int16, single, or double, and it
% must be nonsparse. LEVEL and EM are double scalars.
%
% Example
% -------
% I = imread('coins.png');
% level = graythresh(I);
% BW = im2bw(I,level);
% figure, imshow(BW)
%
narginchk(1,1);
validateattributes(I,{'uint8','uint16','double','single','int16'},{'nonsparse'}, ...mfilename,'I',1);if ~isempty(I)% Convert all N-D arrays into a single column. Convert to uint8 for% fastest histogram computation.I = im2uint8(I(:));num_bins = 256;counts = imhist(I,num_bins);% Variables names are chosen to be similar to the formulas in% the Otsu paper.p = counts / sum(counts);omega = cumsum(p);mu = cumsum(p .* (1:num_bins)');mu_t = mu(end);sigma_b_squared = (mu_t * omega - mu).^2 ./ (omega .* (1 - omega));% Find the location of the maximum value of sigma_b_squared.% The maximum may extend over several bins, so average together the% locations. If maxval is NaN, meaning that sigma_b_squared is all NaN,% then return 0.maxval = max(sigma_b_squared);isfinite_maxval = isfinite(maxval);if isfinite_maxvalidx = mean(find(sigma_b_squared == maxval));% Normalize the threshold to the range [0, 1].level = (idx - 1) / (num_bins - 1);elselevel = 0.0;end
elselevel = 0.0;isfinite_maxval = false;
end% compute the effectiveness metric
if nargout > 1if isfinite_maxvalem = maxval/(sum(p.*((1:num_bins).^2)') - mu_t^2);elseem = 0;end
end
opencv函数实现：
po下opencv的源码：
getThreshVal_Otsu_8u( const Mat& _src ) { Size size = _src.size(); if( _src.isContinuous() ) { size.width *= size.height; size.height = 1; } const int N = 256; int i, j, h[N] = {0}; for( i = 0; i < size.height; i++ ) { const uchar* src = _src.data + _src.step*i; j = 0; #if CV_ENABLE_UNROLLED for( ; j <= size.width - 4; j += 4 ) { int v0 = src[j], v1 = src[j+1]; h[v0]++; h[v1]++; v0 = src[j+2]; v1 = src[j+3]; h[v0]++; h[v1]++; } #endif for( ; j < size.width; j++ ) h[src[j]]++; } double mu = 0, scale = 1./(size.width*size.height); for( i = 0; i < N; i++ ) mu += i*(double)h[i]; mu *= scale; double mu1 = 0, q1 = 0; double max_sigma = 0, max_val = 0; for( i = 0; i < N; i++ ) { double p_i, q2, mu2, sigma; p_i = h[i]*scale; mu1 *= q1; q1 += p_i; q2 = 1. - q1; if( std::min(q1,q2) < FLT_EPSILON || std::max(q1,q2) > 1. - FLT_EPSILON ) continue; mu1 = (mu1 + i*p_i)/q1; mu2 = (mu - q1*mu1)/q2; sigma = q1*q2*(mu1 - mu2)*(mu1 - mu2); if( sigma > max_sigma ) { max_sigma = sigma; max_val = i; } } return max_val; }
int otsuThreshold(IplImage *frame)
{const int GrayScale = 256;int width = frame->width;int height = frame->height;int pixelCount[GrayScale];float pixelPro[GrayScale];int i, j, pixelSum = width * height, threshold = 0;uchar* data = (uchar*)frame->imageData;for (i = 0; i < GrayScale; i++){pixelCount[i] = 0;pixelPro[i] = 0;}//统计灰度级中每个像素在整幅图像中的个数 for (i = 0; i < height; i++){for (j = 0; j < width; j++){pixelCount[(int)data[i * width + j]]++; }}//计算每个像素在整幅图像中的比例 float maxPro = 0.0;int kk = 0;for (i = 0; i < GrayScale; i++){pixelPro[i] = (float)pixelCount[i] / pixelSum;if (pixelPro[i] > maxPro){maxPro = pixelPro[i];kk = i;}}//遍历灰度级[0,255]， i作为阈值 float w0, w1, u0tmp, u1tmp, u0, u1, u, deltaTmp, deltaMax = 0;for (i = 0; i < GrayScale; i++) {w0 = w1 = u0tmp = u1tmp = u0 = u1 = u = deltaTmp = 0;for (j = 0; j < GrayScale; j++){if (j <= i) //背景部分 {w0 += pixelPro[j];u0tmp += j * pixelPro[j];}else //前景部分 {w1 += pixelPro[j];u1tmp += j * pixelPro[j];}}u0 = u0tmp / w0;u1 = u1tmp / w1;u = u0tmp + u1tmp;deltaTmp = w0 * pow((u0 - u), 2) + w1 * pow((u1 - u), 2);if (deltaTmp > deltaMax){deltaMax = deltaTmp;threshold = i;}}return threshold;
}参考：
https://en.wikipedia.org/wiki/Otsu%27s_method
http://www.cnblogs.com/ranjiewen/p/6385564.html
http://blog.csdn.net/glouds/article/details/38976573
