Revision history [back]

Mat GrayHistMatching(Mat I, Mat R) {

/* Histogram Matching of a gray image with a reference*/
// accept two images I (input image)  and R (reference image)

Mat Result;     // The Result image

int L = 256;    // Establish the number of bins
if(I.channels()!=1)
{
    cout<<"Please use Gray image"<<endl;
    return Mat::zeros(I.size(),CV_8UC1);
}
Mat G,S,F; //G is the reference CDF, S the CDF of the equlized given image, F is the map from S->G
if(R.cols>1)
{
    if(R.channels()!=1)
    {
        cout<<"Please use Gray reference"<<endl;
        return Mat::zeros(I.size(),CV_8UC1);
    }


    Mat R_hist, Rp_hist; //R_hist the counts of pixels for each level, Rp_hist is the PDF of each level
    /// Set the ranges ( for B,G,R) )
    float range[] = { 0, 256 } ;
    const float* histRange = { range };
    bool uniform = true; bool accumulate = false;
    //calcHist(const Mat* images, int nimages, const int* channels, InputArray mask, OutputArray hist, int dims, const int* histSize, const float**         ranges, bool uniform=true, bool accumulate=false )
    calcHist( &R, 1, 0, Mat(), R_hist, 1, &L, &histRange, uniform, accumulate );
    //Calc PDF of the image
    Rp_hist=R_hist/(I.rows*I.cols);
    //calc G=(L-1)*CDF(p)
    Mat CDF=cumsum(Rp_hist);
    G=(L-1)*CDF;
    for(int i=0;i<G.rows;i++)
        G.at<Point2f>(i,0).x=(float)cvRound(G.at<Point2f>(i,0).x);//round G
}
else
{
    //else, the given R is the reference PDF
    Mat CDF=cumsum(R);
    G=(L-1)*CDF;
    for(int i=0;i<G.rows;i++)
        G.at<Point2f>(i,0).x=(float)cvRound(G.at<Point2f>(i,0).x);//round G
}
/// Establish the number of bins
Mat S_hist, Sp_hist; //S_hist the counts of pixels for each level, Sp_hist is the PDF of each level
/// Set the ranges ( for B,G,R) )
float range[] = { 0, 256 } ;
const float* histRange = { range };

bool uniform = true; bool accumulate = false;

//calcHist(const Mat* images, int nimages, const int* channels, InputArray mask, OutputArray hist, int dims, const int* histSize, const float**         ranges, bool uniform=true, bool accumulate=false )
calcHist( &I, 1, 0, Mat(), S_hist, 1, &L, &histRange, uniform, accumulate );

//Calc PDF of the image
Sp_hist=S_hist/(I.rows*I.cols);
//calc s=(L-1)*CDF(p)
Mat CDF=cumsum(Sp_hist);
S=(L-1)*CDF;
for(int i=0;i<S.rows;i++)
    S.at<Point2f>(i,0).x=(float)cvRound(S.at<Point2f>(i,0).x);//round S

F=Mat::zeros(S.size(),CV_32F);
int minIndex=-1;
double T,min=100000;
for(int i=0;i<S.rows;i++)
{
    for(int j=0;j<G.rows;j++)
    {
        T=abs(double(S.at<Point2f>(i,0).x)-double(G.at<Point2f>(j,0).x));
        if (T==0)
        {
            minIndex=j;
            break;
        }
        else
            if(T<min)
            {
                minIndex=j;
                min=T;
            }
    }
    F.at<Point2f>(i,0).x=(float)minIndex;
    minIndex=-1;
    min=1000000;
}
uchar table[256]; 
for(int i=0;i<256;i++)
{
    table[i]=(int)F.at<Point2f>(i,0).x;
}
Result= ScanImageAndReduceC(I, table); 

return Result;

}