Max энтропии с использованием OpenCV пороговой

Question

Я пытаюсь преобразовать код для использования максимальной энтропии порогового из this matlab code:

%************************************************************************** 
%************************************************************************** 
% 
% 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("here\n"); 
       hl.push_back(hltmp.at<float>(0)-(histogram.at<float>  (i)/cl)*log(histogram.at<float>(i)/cl)); 
        printf("here\n"); 
       cout<<hl.at<float>(i); 
        printf("here\n"); 
       hltmp.push_back(hl.row(t)); 
        printf("here\n"); 
      } 
     } 
    } 
    printf("here\n"); 
    //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

Ваша проблема в том, что вы читаете float элементы из CV_8U (aka uchar) Mat.

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

Вы должны знать, как использовать отладчик, и вы узнаете, эти проблемы очень скоро.

---

Поскольку вы чрезмерно сложные вещи в вашей реализации, сделал несколько ошибок, и код обрастает от отладки гравюр, я предлагаю один ниже вместо точечно корректировать ваши (не много, но в основном концептуальный) ошибки. Вы можете видеть, что при правильном написании почти 1: 1 конверсия из Matlab в OpenCV.

#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; 
}