According to the HISTEQ function of Matlab rewrite the histogram under the OpenCV to set the C source code!

It is said that the histogram of the image of the regulation of the use of more than the histogram equalization, but it is very strange that OpenCV actually no image histogram rules of the source code! Therefore, it is necessary for me to write an image histogram OpenCV function to facilitate future use.

I found on the internet a few histogram are dilute the source code, and based on OpenCV to rewrite these source code, the effect is not as good as the Histeq function of MATLAB, which rewrite the hardships and cumbersome will not elaborate. Finally, no way, had to learn matalb Histeq function source code, and on it based on OpenCV rewrite.

Although I finally rewrite the success, but the algorithm still do not understand, only according to the help of MATLAB in the document to explain the Histeq, roughly understand its principle, the explanation is as follows:

When you supply a desired histogram hgram, Histeq chooses the grayscale transformation T to Minimize


where C0 is the cumulative histogram of A, C1 are the cumulative sum of Hgram for all intensities K. This minimization was subject to the constraints that T must be monotonic (monotonic) and C1 (T (a)) cannot overshoot (over) C0 (a) by M Ore than half the distance between the histogram counts at A. Histeq uses the transformation B = T (a) to map the gray Leve LS in X (or the ColorMap) to their new values.

On the basis of rewriting success, I will try to translate the above sentence:

For the HISTEQ function if you provide a second argument, and the second parameter is a histogram vector, then the HISTEQ function will calculate a mapping T, using this mapping to map the pixel values of the original image to the values after the histogram match.

This mapping t satisfies the following requirements:

The ① guarantee is minimal, where c0 is the histogram accumulation function of the histogram-defined image, and C1 is the histogram accumulation function of the standard histogram,

②t are monotonically increasing.

③ after T mapping, the value of C1 at point A cannot be straight that C0 is half the value of a point

The principle is that so much above, the following to the source:

The download links for the required images in the source code are:

Coins.png Http://pan.baidu.com/s/1slilbPF

Rice.png http://pan.baidu.com/s/1cDNVx8

First give the MATLAB Histeq function as a histogram of the source code extracted from the source bar:

Clear all;
Close all;
CLC.
I1 = imread (' COINS. PNG ');
I2 = imread (' rice. PNG ');

A = I1.

Hgram = imhist (I2);

N = 256;
Hgram = hgram * (numel (a)/sum (hgram)); % Set sum = numel(a)
M = length (hgram);

% [nn, cum] = computeCumulativeHistogram (a, n);
Nn = imhist (a, n) ';
Cum = cumsum (nn);

A, % T = createTransformationToIntensityImage (hgram, m, n, nn, cum).
% Create the transformation to the an intensity image by minimizing the error
% between desired and actual cumulative histogram.

Cumd = cumsum (hgram * numel (a)/sum (hgram));

Tol = ones (m, 1) * min ([nn (1: n - 1), 0, 0, nn (2: n]]) / 2;
Err = ones (cumd (:) * (1, n) - 'ones (m, 1) * cum (:)') + tol;
D = find(err < -numel(a)* SQRT (eps));
If ~isempty(d) % if d isempty then the value of isempty is 1 and the value of ~isempty is 0, the if statement is not executed
Err (d) = numel (a) * ones (size (d));
The end
[dum, T] = min (err); %
T = (T - 1)/(m - 1); This changes the grayscale from 1 to 256 to 0 to 255
T = T * 255; The gray-level mapping is stored in %T. For example, if the value of 50 is 44, it means that the new gray-level of the pixel in the original image is 44

can be out of the final histogram after the processing of the image of the MATLAB source code as follows:

Clear all;
Close all;
CLC.
I1 = imread (' COINS. PNG ');
I2 = imread (' rice. PNG ');
Histv = imhist (I2);
Histv = histv ';
Fid = fopen (' C: \ Users \ \ Administrator \ Documents \ MATLAB histv3. TXT ', 'wt);
Fprintf (fid, '% g, histv);
The fclose (fid);
Histv = histv '; And then I'm going to transpose it back
J = histeq (I1, histv);
Ranghistv = histv. / 65536;

HistI1 = imhist (I1);

Subplot (1, 2, 1);
Imshow (I1);
Subplot (1,2,2);
Imshow (J);

Then give the OPENCV environment C code Source:

# include < opencv2 / opencv. HPP >
# include < opencv2 / legacy/compat HPP >
# include < fstream >
Using namespace STD.
# pragma comment (would, "/ subsystem: \ Windows \" mainCRTStartup "/ entry: \" \ "")


Void mycvCalcHist(IplImage *img,double out_hist[256]
{
Int I = 0, j = 0;
Double temp1 = 0;
Int temp2 = 0;
Const int hist_sz = 256; //0 to 255, a total of 256 grayscale values
Double hist [hist_sz];
Memset (hist, 0, sizeof (hist));
For (I = 0; I < img - > height; I++)
{
For (j = 0; J < img - > width; J++)
{temp1 = cvGet2D (img, I, j). Val [0].
Temp2 = int (temp1); // type conversion
Hist [temp2] + +; // the frequency of each gray value stored in hist is realized here
}
}

Memcpy (out_hist, hist, sizeof (hist)); // someone must ask why you don't pass the array name as an argument to change the value of the real argument set
// this method is generally ok, I have also tested, but here is not, I estimate and
/ / memset (hist, 0, sizeof (hist)); This sentence is relevant

}

Void my_hist_specification(IplImage *a,IplImage *I2)//a is a pointer to the source image, and I2 is a pointer to the image corresponding to the standard histogram required by a for histogram specification
{

Int I = 0;
Int j = 0;
Int n = 256;

// calculate the gray histogram of the image
Double hgram [256].
MycvCalcHist (I2, hgram);
Double nn [256].
MycvCalcHist (a, nn);

// is equivalent to hgram in M = hgram*(numel(a)/sum(hgram));
Double sum_hgram = 0;
For (I = 0; I < 256; I++)
Sum_hgram = sum_hgram + hgram [I];
Double M_a N_a;
M_a = a - > height;
N_a = a - > width;
Double numel_a;
Numel_a = M_a * N_a;
Double numel_aDivsum_hgram;
Numel_aDivsum_hgram = numel_a/sum_hgram;
For (I = 0; I < 256; I++)
Hgram [I] = hgram [I] * numel_aDivsum_hgram;

// is equivalent to m = length(hgram);
Int m;
M = 256;

Cum = cumsum(nn);
Double cum [256].
Double val_1 = 0;
Int index;
For (index = 0; The index < 256; Index++)
{
Val_1 = val_1 + nn [index];
Cum [index] = val_1;
}

Cumd = cumsum(hgram*numel(a)/sum(hgram));
Double cumd [256].
Double cumd_temp1 [256].
Double cumd_temp2;
Sum_hgram = 0;
For (I = 0; I < 256; I++)
Sum_hgram = sum_hgram + hgram [I];
Cumd_temp2 = numel_a/sum_hgram;
For (I = 0; I < 256; I++)
Cumd_temp1 [I] = hgram [I] * cumd_temp2;
Val_1 = 0;
For (index = 0; The index < 256; Index++)
{
Val_1 = val_1 + cumd_temp1 [index];
Cumd [index] = val_1;
}

/ / equivalent to tol = ones (m, 1) * min ([nn (1: n - 1), 0; 0, nn (2: n]]) / 2;
Double tol_temp1 [256].
For (index = 0; The index < 256; Index++)
Tol_temp1 [index] = nn [index] / 2;
Tol_temp1 [0] = 0;
Tol_temp1 [255] = 0;
The static double tol [256] [256].
For (index = 0; The index < 256; Index++)
Memcpy (tol [index], tol_temp1, sizeof (tol_temp1)); // here the value is correct, equivalent to a row copy of tol_temp1;

/ / equivalent to err = ones (cumd (:) * (1, n) - 'ones (m, 1) * cum (:)') + tol;
The static double err_1 [256] [256].
For (I = 0; I < 256; I++)
For (j = 0; J < 256; J++)
[j] = cumd err_1 [I] [I]; // here the value is correct, and cumd is copied as a column; The second check is also correct

The static double err_2 [256] [256].
For (index = 0; The index < 256; Index++)
Memcpy (err_2 [index], cum, sizeof (cum)); // here the value is correct and the row copy is cum
The static double err [256] [256].
For (I = 0; I < 256; I++)
For (j = 0; J < 256; J++)
Err [I] [j] = err_1 [I] [j] - err_2 [I] [j] + tol [I] [j];

// equivalent to d = find(err < -numel(a)* SQRT (eps));
Double matlab_eps e-008 = 1.4901;
Double find_err_temp1 = - numel_a * matlab_eps;
Int k = 0;
Int breakvari = 0;
Double watch_err;
For (I = 0; I < 256; I++) // this is column traversal, not row traversal
{
For (j = 0; J < 256; J++)
If (err [j] [I] < find_err_temp1)
{watch_err = err [j] [I];
K++;
}
}
Double *d = (double *)malloc(k*sizeof(double));
K = 0;
For (I = 0; I < 256; I++) // this is column traversal, not row traversal
{
For (j = 0; J < 256; J++)
If (err [j] [I] < find_err_temp1)
{watch_err = err [j] [I];
D [k] = I * 256 + j;
K++;
}
}

/ *
The equivalent of
If ~isempty(d) % if d isempty then the value of isempty is 1 and the value of ~isempty is 0, the if statement is not executed
Err (d) = numel (a) * ones (size (d));
The end
* /
Int d_size;
D_size = k;
Char isempty_flag = 1;
For (I = 0; I < d_size; I++)
{
If (d [I]! = 0)
{
Isempty_flag = 0;
Break;
}
}
Double err_temp1 = 0;
Int d_m d_n;
If (! Isempty_flag)
{
For (I = 0; I < d_size; I++)
{
[I] d_m = int (d) % 256;
[I] d_n = int (d) / 256;
Err [d_m] [d_n] = numel_a;

}

}
// is equivalent to [dum,T] = min(err);
Int T [256];
Double err_min = 0;
For (j = 0; J < 256; J++)
{err_min = err [0] [j];
T [j] = 0;
For (I = 0; I < 256; I++)
{if (err [I] [j] < err_min)
{
T [j] = I;
Err_min = err [I] [j]; //T is the pixel value mapping after histogram matching!
}
}
}

// the following program maps each pixel value in the original image to the mapping relation in T
T [0] = 0; // this is the grayscale value mapping I calculated through the previous program, no matter how you map,0 must be 0, right
Int s1_temp1;
CvScalar s1.
CvScalar s2.
For (I = 0; I < M_a; I++)
{
For (j = 0; J < N_a; J++)
{
S1 = cvGet2D (a, I, j);
S1_temp1 = int (s1. Val [0]).
S2. Val [s1_temp1] [0] = T;
CvSet2D (a, I, j, s2);
}

}

}

Int main ()
{
// load the original image from the file
IplImage *pSrcImage_coins = cvLoadImage(" minutes.png ", CV_LOAD_IMAGE_UNCHANGED);

IplImage *a = cvLoadImage(" minutes.png ", CV_LOAD_IMAGE_UNCHANGED);
IplImage *I2 = cvLoadImage("rice.png", CV_LOAD_IMAGE_UNCHANGED);


My_hist_specification (a, I2);

Const char *pstrWindowsATitle = "original ";
Const char *pstrWindowsBTitle = "histogram specified source map ";
// create a window
CvNamedWindow (pstrWindowsATitle, CV_WINDOW_AUTOSIZE);
CvNamedWindow (pstrWindowsBTitle, CV_WINDOW_AUTOSIZE);
// displays the image in the specified window
CvShowImage (pstrWindowsATitle pSrcImage_coins);
CvShowImage (pstrWindowsBTitle, a);
// wait for key events
CvWaitKey ();
CvDestroyWindow (pstrWindowsATitle);
CvDestroyWindow (pstrWindowsBTitle);
cvReleaseImage(&a);    
cvReleaseImage(&I2); 

return 0;
}

The final run results are as follows:

I tested it with another picture, and there was no problem, as follows:

Visible, the results are consistent, the program test success! In addition, in my backup source histogram_specification_05.cpp, download link: http://pan.baidu.com/s/1i5g0Kzr full of my entire rewrite process, I set a lot of intermediate variables and arrays, You can see if the program is correct!

Welcome everyone to join the image recognition technology Exchange Group: 271891601, in addition, the special welcome Chengdu engaged in image recognition work of friends, my QQ number 248787278

According to the HISTEQ function of Matlab rewrite the histogram under the OpenCV to set the C source code!