使用OpenCV的最大熵阈值

Question

我正在尝试从这个matlab代码转换使用最大熵阈值的代码：

%**************************************************************************
%**************************************************************************
%   
% maxentropie is a function for thresholding using Maximum Entropy
% 
% 
% input = I ==> Image in gray level 
% output =
%           I1 ==> binary image
%           threshold ==> the threshold choosen by maxentropie
%  
% F.Gargouri
%
%
%**************************************************************************
%**************************************************************************

function [threshold I1]=maxentropie(I)
    [n,m]=size(I);
    h=imhist(I);
    %normalize the histogram ==>  hn(k)=h(k)/(n*m) ==> k  in [1 256]
    hn=h/(n*m);
    %Cumulative distribution function
    c(1) = hn(1);
    for l=2:256
        c(l)=c(l-1)+hn(l);
    end

    hl = zeros(1,256);
    hh = zeros(1,256);
    for t=1:256
        %low range entropy
        cl=double(c(t));
        if cl>0
            for i=1:t
                if hn(i)>0
                    hl(t) = hl(t)- (hn(i)/cl)*log(hn(i)/cl);                      
                end
            end
        end
        %high range entropy
        ch=double(1.0-cl);  %constraint cl+ch=1
        if ch>0
            for i=t+1:256
                if hn(i)>0
                    hh(t) = hh(t)- (hn(i)/ch)*log(hn(i)/ch);
                end
            end
        end
    end
    % choose best threshold
    h_max =hl(1)+hh(1)
    threshold = 0;
    entropie(1)=h_max;
    for t=2:256
        entropie(t)=hl(t)+hh(t);
        if entropie(t)>h_max
            h_max=entropie(t);
            threshold=t-1;
        end
    end
    % Display    
    I1 = zeros(size(I));
    I1(I<threshold) = 0;
    I1(I>threshold) = 255;
    %imshow(I1)      
end 

问题是我在代码中出现了浮点excpetion错误，我无法理解为什么

这是我的实现：

#include <iostream>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <math.h>
using namespace cv;
using namespace std;
int main(){
cout.setf(std::ios_base::fixed, std::ios_base::floatfield);
cout.precision(9);
Mat old_image=imread("2.png",CV_LOAD_IMAGE_GRAYSCALE);
double minval, maxval;
minMaxLoc(old_image,&minval, &maxval);
cout<<minval<<" "<<maxval<<endl;
Mat image;
old_image.convertTo(image, CV_8UC1, 255.0/(maxval-minval), -minval*255.0/(maxval-minval));
minMaxLoc(image,&minval, &maxval);
cout<<minval<<" "<<maxval;
int k=0;
imshow("im",image);
waitKey(0);
for(int y=0; y<image.rows;y++){
    for(int x=0; x<image.cols;x++){
        if((int) image.at<uchar>(y,x)==0){
            k++;
        }
    }
}
cout<<k<<endl<<endl;

int i, l, j, t;
int histSize = 256;
float range[] = { 0, 255 };
const float *ranges[] = { range };
Mat hist, histogram, c, ctmp, hl,  hh, hhtmp, entropy;
calcHist( &image, 1, 0, Mat(), hist, 1, &histSize, ranges, true, false );
for( int h = 1; h < histSize; h++){
    histogram.push_back(hist.at<float>(h,0));
    cout<<histogram.rows<<endl;
    cout<<histogram.row(h-1)<<endl;
    cout<<hist.row(h)<<endl;
}
histogram=histogram/(image.rows*image.cols-hist.at<float>(0,0));
//cumulative distribution function
float cl,ch;
ctmp.push_back(histogram.row(0));
c.push_back(histogram.row(0));
cout<<c.row(0)<<endl;
for(l=1;l<255;l++){
    c.push_back(ctmp.at<float>(0)+histogram.at<float>(l));
    ctmp.push_back(c.row(l));
    cout<<c.at<float>(l)<<endl;
    //c.row(l)=c.row(l-1)+histogram.row(l);
}
 Mat hltmp= Mat::zeros(1,256,CV_8U);
// THE PROBLEM IS IN THIS TWO FOR CYCLES 
for(t=0;t<255;t++){
    //low range entropy
    cl=c.at<float>(t);
    if(cl>0){
        for(i=0;i<=t;i++){
            if(histogram.at<float>(t)>0){
                 printf("heren");
                hl.push_back(hltmp.at<float>(0)-(histogram.at<float>        (i)/cl)*log(histogram.at<float>(i)/cl));
                   printf("heren");
                 cout<<hl.at<float>(i);
                  printf("heren");
                 hltmp.push_back(hl.row(t));
                     printf("heren");
            }
        }
    }
    printf("heren");
    //high range entropy
    ch=1.0-cl;
    if(ch>0){
        for(i=t+1;i<255;i++){
            if(histogram.at<float>(i)>0){
                hh.push_back(hh.at<float>(t)-(histogram.at<float>   (i)/ch)*log(histogram.at<float>(i)/ch));
            }
        }
    }
}
 //choose the best threshold
float h_max=hl.at<float>(0,0)+hh.at<float>(0,0);
float threshold=0;
entropy.at<float>(0,0)=h_max;
for(t=1;t<255;t++){
    entropy.at<float>(t)=hl.at<float>(t)+hh.at<float>(t);
    if(entropy.at<float>(t)>h_max){
        h_max=entropy.at<float>(t);
        threshold=t-1;
    }
    cout<<threshold<<endl;
}
//display
Mat I1= Mat::zeros(image.rows,image.cols,CV_8UC1);
for(int y=0; y<image.rows;y++){
    for(int x=0; x<image.cols;x++){
        if((int) image.at<uchar>(y,x)<threshold){
            I1.at<uchar>(y,x)=0;
        }
        else{
            I1.at<uchar>(y,x)=255;
        }
    }
}
imshow("image",I1);
waitKey(0);*/
return 0;
}

Answer 1

您的问题是从CV_8U（也称为uchar）Mat读取float元素。

Mat hltmp = Mat::zeros(1, 256, CV_8U);  
...
hltmp.at<float>(0)

您应该学习如何使用调试器，很快就会发现这些问题。

---

由于您在实现中处理了复杂的事情，犯了一些错误，并且调试打印的代码很混乱，所以我建议使用下面的方法，而不是及时纠正您的（不多，但主要是概念性的）错误。您可以看到，如果编写得当，从Matlab到OpenCV的转换几乎是1:1。

#include <opencv2/opencv.hpp>
#include <iostream>
using namespace std;
using namespace cv;
uchar maxentropie(const Mat1b& src, Mat1b& dst)
{
    // Histogram
    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
    Mat1d cumhist(1, 256, 0.0);
    float sum = 0;
    for (int i = 0; i < 256; ++i)
    {
        sum += hist(i);
        cumhist(i) = sum;
    }
    Mat1d hl(1, 256, 0.0);
    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
    Mat1d entropie(1, 256, 0.0);
    double h_max = hl(0) + hh(0);
    uchar 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);
        }
    }
    // Create output image
    dst = src > threshold;
    return threshold;
}
int main()
{
    Mat1b img = imread("path_to_image", IMREAD_GRAYSCALE);
    Mat1b res;
    uchar th = maxentropie(img, res);
    imshow("Original", img);
    imshow("Result", res);
    waitKey();
    return 0;
}