大津法---OTSU算法

簡介：

大津法（OTSU）是一種確定圖像二值化分割閾值的算法，由日本學者大津於1979年提出。從大津法的原理上來講，該方法又稱作最大類間方差法，因為按照大津法求得的閾值進行圖像二值化分割后，前景與背景圖像的類間方差最大（何為類間方差？原理中有介紹）。

OTSU算法

OTSU算法也稱最大類間差法，有時也稱之為大津算法，由大津於1979年提出，被認為是圖像分割中閾值選取的最佳算法，計算簡單，不受圖像亮度和對比度的影響，因此在數字圖像處理上得到了廣泛的應用。它是按圖像的灰度特性,將圖像分成背景和前景兩部分。因方差是灰度分布均勻性的一種度量,背景和前景之間的類間方差越大,說明構成圖像的兩部分的差別越大,當部分前景錯分為背景或部分背景錯分為前景都會導致兩部分差別變小。因此,使類間方差最大的分割意味着錯分概率最小。

原理：

對於圖像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)
將式(5)代入式(6),得到等價公式:
　　　　　　g=ω0ω1(μ0-μ1)^2 　　 (7)　這就是類間方差
采用遍歷的方法得到使類間方差g最大的閾值T,即為所求。

matlab函數：

matlab中函數graythresh既是使用大津法求得分割閾值Ｔ。用法如下：

　　　　　　T = graythresh(img);

　　　　　　BW = im2bw(img,T);

大津法的形象理解：

對於直方圖有兩個峰值的圖像，大津法求得的Ｔ近似等於兩個峰值之間的低谷。

　　　　　　imhist(img);

T = graythresh(img);

如下圖為圖像的直方圖，使大津法求得的Ｔ＝0.5294，轉換在［0,255］之間為134.9970，只好是兩個峰值之間低谷的位置。

OpenCV的二值化操作中，有一種“大津閾值處理”的方法，使用函數cvThreshold(image,image2,0,255,CV_THRESH_OTSU) 實現，該函數就會使用大律法OTSU得到的全局自適應閾值來進行二值化圖片，而參數中的threshold不再起作用。

以下是一段在OpenCV中實現的C語言程序，即一個使用OTSU算法提取圖像閾值的函數，輸入參數為一個圖像指針，返回分割該圖像的最佳閾值。

其中的變量說明：當分割的閾值為t時

w0為背景像素點占整幅圖像的比例

u0為w0平均灰度

w1為前景像素點占整幅圖像的比例

u1為w1平均灰度

u為整幅圖像的平均灰度

公式:g = w0*pow((u-u0),2) + w1*pow((u-u1),2)

int MyAutoFocusDll::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]  
    float w0, w1, u0tmp, u1tmp, u0, u1, u, deltaTmp, deltaMax = 0;
    for (i = 0; i < GrayScale; i++)     // 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;
}