如何使用直方图的峰值自动选择阈值？

Question

通过 OpenCV 库，我想对这样的图像进行阈值处理：

threshold(image, thresh, 220, 255, THRESH_BINARY_INV)

但我想自动找到阈值（220）。

我使用 Otsu 来估计阈值。但它在我的情况下不起作用。

因此，我应该使用直方图峰值技术。我想在直方图中找到与图像的背景和对象相对应的两个峰值。它会在两个峰值之间自动设置阈值。 我使用这本书（第 117 和 496-505 页）：Dwayne Phillips (http://homepages.inf.ed.ac.uk/rbf/BOOKS/PHILLIPS/) 的“C 中的图像处理”。我使用源代码在直方图中找到与图像的背景和对象相对应的两个峰值。这是我的图片：

这是我的 C++ 代码：

#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <opencv2/opencv.hpp>

#include <iostream>
#include <stdio.h>
#include <fstream>

using namespace std;
using namespace cv;


int main()
{
    Mat image0 = imread("C:/Users/Alireza/Desktop/contrast950318/2.bmp");
    imshow("image0", image0);

    Mat image, thresh, Tafrigh;
    cvtColor(image0, image, CV_RGB2GRAY);


    int N = image.rows*image.cols;

    int histogram[256];

    for (int i = 0; i < 256; i++) {
        histogram[i] = 0;
    }

    //create histo
    for (int i = 0; i < image.rows; i++){
        for (int j = 0; j < image.cols; j++){
            histogram[((int)image.at<uchar>(i, j))]++;
        }
    }


    int peak1, peak2;
    #define PEAKS 30
    int distance[PEAKS], peaks[PEAKS][2];

    int i, j = 0, max = 0, max_place = 0;

    for (int i = 0; i<PEAKS; i++){
        distance[i] = 0;
        peaks[i][0] = -1;
        peaks[i][1] = -1;
    }

    for (i = 0; i <= 255; i++){
        max = histogram[i];
        max_place = i;

        //insert_into_peaks(peaks, max, max_place);
        //int max, max_place, peaks[PEAKS][2];
        //int i, j;
        /* first case */
        if (max > peaks[0][0]){
            for (i = PEAKS - 1; i > 0; i--){
                peaks[i][0] = peaks[i - 1][0];
                peaks[i][1] = peaks[i - 1][1];
            }
            peaks[0][0] = max;
            peaks[0][1] = max_place;
        } /* ends if */

        /* middle cases */
        for (j = 0; j < PEAKS - 3; j++){
            if (max < peaks[j][0] && max > peaks[j + 1][0]){
                for (i = PEAKS - 1; i > j + 1; i--){
                    peaks[i][0] = peaks[i - 1][0];
                    peaks[i][1] = peaks[i - 1][1];
                }

                peaks[j + 1][0] = max;
                peaks[j + 1][1] = max_place;
            } /* ends if */
        } /* ends loop over j */
        /* last case */
        if (max < peaks[PEAKS - 2][0] &&
            max > peaks[PEAKS - 1][0]){
            peaks[PEAKS - 1][0] = max;
            peaks[PEAKS - 1][1] = max_place;
        } /* ends if */

}/* ends loop over i */



    for (int i = 1; i<PEAKS; i++){
        distance[i] = peaks[0][1] - peaks[i][1];
        if (distance[i] < 0)
            distance[i] = distance[i] * (-1);
    }

    peak1 = peaks[0][1];
    cout << "  peak1= " << peak1;

    for (int i = PEAKS - 1; i > 0; i--){
        if (distance[i] > 1)
            peak2 = peaks[i][1];

    }
    cout << "  peak2= " << peak2;


    int mid_point;
    //int peak1, peak2;
    short hi, low;

    unsigned long sum1 = 0, sum2 = 0;
    if (peak1 > peak2)
        mid_point = ((peak1 - peak2) / 2) + peak2;
    if (peak1 < peak2)
        mid_point = ((peak2 - peak1) / 2) + peak1;

    for (int i = 0; i<mid_point; i++)
        sum1 = sum1 + histogram[i];
    for (int i = mid_point; i <= 255; i++)
        sum2 = sum2 + histogram[i];
    if (sum1 >= sum2){
        low = mid_point;
        hi = 255;
    }
    else{
        low = 0;
        hi = mid_point;
    }

    cout << "  low= " << low << "  hi= " << hi;



double  threshnum = 0.5* (low  + hi);
threshold(image, thresh, threshnum, hi, THRESH_BINARY_INV);


    waitKey(0);
    return 0;
}

但我不知道这段代码是否正确。如果正确，为什么阈值是 202？

您对如何解决此任务有什么建议？或者我可以在 Internet 上的哪些资源上找到帮助？

Answer 1

您也可以使用最大熵。在某些情况下，只使用熵的高频可能会更好

int maxentropie(const cv::Mat1b& src)
{
    // Histogram
    cv::Mat1d hist(1, 256, 0.0);
    for (int r=0; r<src.rows; ++r)
        for (int c=0; c<src.cols; ++c)
            hist(src(r,c))++;

    // Normalize
    hist /= double(src.rows * src.cols);

    // Cumulative histogram
    cv::Mat1d cumhist(1, 256, 0.0);
    float sum = 0;
    for (int i = 0; i < 256; ++i)
    {
        sum += hist(i);
        cumhist(i) = sum;
    }

    cv::Mat1d hl(1, 256, 0.0);
    cv::Mat1d hh(1, 256, 0.0);

    for (int t = 0; t < 256; ++t)
    {
        // low range entropy
        double cl = cumhist(t);
        if (cl > 0)
        {
            for (int i = 0; i <= t; ++i)
            {
                if (hist(i) > 0)
                {
                    hl(t) = hl(t) - (hist(i) / cl) * log(hist(i) / cl);
                }
            }
        }

        // high range entropy
        double ch = 1.0 - cl;  // constraint cl + ch = 1
        if (ch > 0)
        {
            for (int i = t+1; i < 256; ++i)
            {
                if (hist(i) > 0)
                {
                    hh(t) = hh(t) - (hist(i) / ch) * log(hist(i) / ch);
                }
            }
        }
    }

    // choose best threshold

    cv::Mat1d entropie(1, 256, 0.0);
    double h_max = hl(0) + hh(0);
    int threshold = 0;
    entropie(0) = h_max;

    for (int t = 1; t < 256; ++t)
    {
        entropie(t) = hl(t) + hh(t);
        if (entropie(t) > h_max)
        {
            h_max = entropie(t);
            threshold = uchar(t);
        }
    }
    if(threshold==0) threshold=255;
    return threshold;
}