OpenCV—反向投影直方图检测特定图像内容

Opencv计算机视觉编程攻略（第二版）在第4.5节介绍了反向投影直方图检测图像中的特定内容。本文为学习笔记、示例程序的实现方法以及自己的一些体会。本例的工程共包含3个头文件和一个源文件。3个头文件分别定义了1维直方图Histogram1D和3维直方图ColorHistogram操作方法，以及直方图检测方法ContentFinder的类。源文件为主程序。

#include "histogram.h"               //定义1维灰度直方图的获取与绘制方法类：class Histogram1D
#include "colorhistogram.h"          //定义3维彩色直方图的获取与绘制方法类：class ColorHistogram 
#include "contentFinder.h"           //定义反向投影直方图检测方法的类：class ContentFinder

本例中使用的图像如下：

1.直方图操作方法

（1）1维灰度直方图操作方法

头文件histogram.h中定义了操作1维灰度直方图的方法类class Histogram1D。histogram.h代码如下：

#if !defined HISTOGRAM
#define HISTOGRAM

#include <opencv2\\core\\core.hpp>
#include <opencv2\\imgproc\\imgproc.hpp>

// To create histograms of gray-level images
class Histogram1D { //将算法封装进类

  private:

    int histSize[1];         // number of bins in histogram直方图中箱子(bin)个数,[1]表示只有1维
	float hranges[2];    // range of values值范围,min和max共2个值，因此定义2维浮点数组
    const float* ranges[1];  // pointer to the different value ranges值范围的指针
    int channels[1];         // channel number to be examined要检查的通道数量

  public:

	Histogram1D() {

		// Prepare default arguments for 1D histogram
		histSize[0]= 256;   // 256 bins,只有1维，因此通过[0]来设置该维的箱子数
		hranges[0]= 0.0;    // from 0 (inclusive)直方图取值范围的min
		hranges[1]= 256.0;  // to 256 (exclusive)直方图取值范围的max
		ranges[0]= hranges; 
		channels[0]= 0;     // we look at channel 0，1维直方图暂时只看0通道
	}

	// Sets the channel on which histogram will be calculated.
	// By default it is channel 0.设置通道的方法
	void setChannel(int c) {

		channels[0]= c;
	}

	// Gets the channel used.获取通道的方法
	int getChannel() {     

		return channels[0];
	}

	// Sets the range for the pixel values.设置直方图值的范围
	// By default it is [0,256]
	void setRange(float minValue, float maxValue) {

		hranges[0]= minValue;
		hranges[1]= maxValue;
	}

	// Gets the min pixel value.
	float getMinValue() {

		return hranges[0];
	}

	// Gets the max pixel value.
	float getMaxValue() {

		return hranges[1];
	}

	// Sets the number of bins in histogram.设置直方图箱子数（统计多少个灰度级）
	// By default it is 256.构造函数中默认设置为256
	void setNBins(int nbins) {

		histSize[0]= nbins;
	}

	// Gets the number of bins in histogram.
	int getNBins() {

		return histSize[0];
	}

	// Computes the 1D histogram.自编函数计算1维直方图
	cv::Mat getHistogram(const cv::Mat &image) {//输入图像image

		cv::Mat hist;

		// Compute histogram
		cv::calcHist(&image, 
			1,		// histogram of 1 image only
			channels,	// the channel used
			cv::Mat(),	// no mask is used，不使用掩码
			hist,		// the resulting histogram
			1,		// it is a 1D histogram
			histSize,	// number of bins
			ranges		// pixel value range
		);

		return hist;
	}

    // Computes the 1D histogram and returns an image of it.直方图图像
	cv::Mat getHistogramImage(const cv::Mat &image, int zoom = 1){

		// Compute histogram first
		cv::Mat hist = getHistogram(image);

		// Creates image
        return Histogram1D::getImageOfHistogram(hist, zoom);
	}

    // Stretches the source image using min number of count in bins.
    cv::Mat stretch(const cv::Mat &image, int minValue = 0) {

        // Compute histogram first
        cv::Mat hist = getHistogram(image);

        // find left extremity of the histogram
        int imin = 0;
        for (; imin < histSize[0]; imin++) {
            // ignore bins with less than minValue entries
            if (hist.at<float>(imin) > minValue)
                break;
        }

        // find right extremity of the histogram
        int imax = histSize[0] - 1;
        for (; imax >= 0; imax--) {

            // ignore bins with less than minValue entries
            if (hist.at<float>(imax) > minValue)
                break;
        }

        // Create lookup table
        int dims[1] = { 256 };
        cv::Mat lookup(1, dims, CV_8U);

        for (int i = 0; i<256; i++) {

            if (i < imin) lookup.at<uchar>(i) = 0;
            else if (i > imax) lookup.at<uchar>(i) = 255;
            else lookup.at<uchar>(i) = cvRound(255.0*(i - imin) / (imax - imin));
        }

        // Apply lookup table
        cv::Mat result;
        result = applyLookUp(image, lookup);

        return result;
    }

    // Stretches the source image using percentile.
    cv::Mat stretch(const cv::Mat &image, float percentile) {

        // number of pixels in percentile
        float number= image.total()*percentile;

        // Compute histogram first
        cv::Mat hist = getHistogram(image);

        // find left extremity of the histogram
        int imin = 0;
        for (float count=0.0; imin < histSize[0]; imin++) {
            // number of pixel at imin and below must be > number
            if ((count+=hist.at<float>(imin)) >= number)
                break;
        }

        // find right extremity of the histogram
        int imax = histSize[0] - 1;
        for (float count=0.0; imax >= 0; imax--) {
            // number of pixel at imax and below must be > number
            if ((count += hist.at<float>(imax)) >= number)
                break;
        }

        // Create lookup table
        int dims[1] = { 256 };
        cv::Mat lookup(1, dims, CV_8U);

        for (int i = 0; i<256; i++) {

            if (i < imin) lookup.at<uchar>(i) = 0;
            else if (i > imax) lookup.at<uchar>(i) = 255;
            else lookup.at<uchar>(i) = cvRound(255.0*(i - imin) / (imax - imin));
        }

        // Apply lookup table
        cv::Mat result;
        result = applyLookUp(image, lookup);

        return result;
    }

    // static methods

    // Create an image representing a histogram
    static cv::Mat getImageOfHistogram(const cv::Mat &hist, int zoom) {

        // Get min and max bin values
        double maxVal = 0;
        double minVal = 0;
        cv::minMaxLoc(hist, &minVal, &maxVal, 0, 0);

        // get histogram size
        int histSize = hist.rows;

        // Square image on which to display histogram
        cv::Mat histImg(histSize*zoom, histSize*zoom, CV_8U, cv::Scalar(255));

        // set highest point at 90% of nbins (i.e. image height)
        int hpt = static_cast<int>(0.9*histSize);

        // Draw vertical line for each bin
        for (int h = 0; h < histSize; h++) {

            float binVal = hist.at<float>(h);
            if (binVal>0) {
                int intensity = static_cast<int>(binVal*hpt / maxVal);
                cv::line(histImg, cv::Point(h*zoom, histSize*zoom),
                    cv::Point(h*zoom, (histSize - intensity)*zoom), cv::Scalar(0), zoom);
            }
        }

        return histImg;
    }

    // Equalizes the source image.
    static cv::Mat equalize(const cv::Mat &image) {

		cv::Mat result;
		cv::equalizeHist(image,result);

		return result;
	}

	// Applies a lookup table transforming an input image into a 1-channel image
    static cv::Mat applyLookUp(const cv::Mat& image, // input image
      const cv::Mat& lookup) { // 1x256 uchar matrix

      // the output image
      cv::Mat result;

      // apply lookup table
      cv::LUT(image,lookup,result);

      return result;
    }
};

#endif

（2）3维彩色直方图操作方法

头文件ColorHistogram.h中定义了操作彩色直方图的方法类class ColorHistogram。colorhistogram.h代码如下：

#if !defined COLHISTOGRAM
#define COLHISTOGRAM

#include <opencv2\\core\\core.hpp>
#include <opencv2\\imgproc\\imgproc.hpp>

class ColorHistogram {

  private:

    int histSize[3];         // size of each dimension
	float hranges[2];    // range of values
    const float* ranges[3];  // array of ranges for each dimension
    int channels[3];         // channel to be considered

  public:

	ColorHistogram() {

		// Prepare default arguments for a color histogram
		// each dimension has equal size and range
		histSize[0]= histSize[1]= histSize[2]= 256;
		hranges[0]= 0.0;    // BRG range from 0 to 256
		hranges[1]= 256.0;
		ranges[0]= hranges; // in this class,  
		ranges[1]= hranges; // all channels have the same range
		ranges[2]= hranges; 
		channels[0]= 0;	    // the three channels 
		channels[1]= 1; 
		channels[2]= 2; 
	}

	// set histogram size for each dimension
	void setSize(int size) {

		// each dimension has equal size 
		histSize[0]= histSize[1]= histSize[2]= size;
	}

	// Computes the histogram.
	cv::Mat getHistogram(const cv::Mat &image) {

		cv::Mat hist;

		// BGR color histogram
		hranges[0]= 0.0;    // BRG range
		hranges[1]= 256.0;
		channels[0]= 0;		// the three channels 
		channels[1]= 1; 
		channels[2]= 2; 

		// Compute histogram
		cv::calcHist(&image, 
			1,		// histogram of 1 image only
			channels,	// the channel used
			cv::Mat(),	// no mask is used
			hist,		// the resulting histogram
			3,		// it is a 3D histogram
			histSize,	// number of bins
			ranges		// pixel value range
		);
		return hist;
	}

	// Computes the histogram.
	cv::SparseMat getSparseHistogram(const cv::Mat &image) {

		cv::SparseMat hist(3,        // number of dimensions
			           histSize, // size of each dimension
				       CV_32F);

		// BGR color histogram
		hranges[0]= 0.0;    // BRG range
		hranges[1]= 256.0;
		channels[0]= 0;	    // the three channels 
		channels[1]= 1; 
		channels[2]= 2; 

		// Compute histogram
		cv::calcHist(&image, 
			1,		// histogram of 1 image only
			channels,	// the channel used
			cv::Mat(),	// no mask is used
			hist,		// the resulting histogram
			3,		// it is a 3D histogram
			histSize,	// number of bins
			ranges		// pixel value range
		);
		return hist;
	}

	// Computes the 1D Hue histogram with a mask.
	// BGR source image is converted to HSV
	// Pixels with low saturation are ignored
	cv::Mat getHueHistogram(const cv::Mat &image, 
                             int minSaturation=0) {
		cv::Mat hist;

		// Convert to HSV colour space
		cv::Mat hsv;
		cv::cvtColor(image, hsv, CV_BGR2HSV);

		// Mask to be used (or not)
		cv::Mat mask;

		if (minSaturation>0) {
			// Spliting the 3 channels into 3 images
			std::vector<cv::Mat> v;
			cv::split(hsv,v);

			// Mask out the low saturated pixels
			cv::threshold(v[1],mask,minSaturation,255,
                                 cv::THRESH_BINARY);
		}

		// Prepare arguments for a 1D hue histogram
		hranges[0]= 0.0;    // range is from 0 to 180
		hranges[1]= 180.0;
		channels[0]= 0;    // the hue channel 

		// Compute histogram
		cv::calcHist(&hsv, 
			1,		// histogram of 1 image only
			channels,	// the channel used
			mask,		// binary mask
			hist,		// the resulting histogram
			1,		// it is a 1D histogram
			histSize,	// number of bins
			ranges		// pixel value range
		);
		return hist;
	}

	// Computes the 2D ab histogram.
	// BGR source image is converted to Lab
	cv::Mat getabHistogram(const cv::Mat &image) {

		cv::Mat hist;

		// Convert to Lab color space
		cv::Mat lab;
		cv::cvtColor(image, lab, CV_BGR2Lab);

		// Prepare arguments for a 2D color histogram
		hranges[0]= 0;
		hranges[1]= 256.0;
		channels[0]= 1; // the two channels used are ab 
		channels[1]= 2; 

		// Compute histogram
		cv::calcHist(&lab, 
			1,		    // histogram of 1 image only
			channels,	    // the channel used
			cv::Mat(),	    // no mask is used
			hist,		    // the resulting histogram
			2,		    // it is a 2D histogram
			histSize,	    // number of bins
			ranges		    // pixel value range
		);
		return hist;
	}
};
#endif

2.计算反向投影直方图的方法

计算反向投影直方图的过程：

• 从ROI区域的归一化的直方图中读取概率值

• 把输入图像（代检测图）中的每一个像素替换成与之对应的归一化中方图中的概率值

• 把替换成概率值（0~1）的像素值再从0~1映射到0~255

• 灰度值越大的像素越有可能是ROI的成分。

本例中，将反向投影直方图封装成一个类ContentFinder，contentfinder.h代码如下：

#if !defined OFINDER
#define OFINDER

#include <opencv2\\core\\core.hpp>
#include <opencv2\\imgproc\\imgproc.hpp>

class ContentFinder {

  private:

	// histogram parameters
	float hranges[2];
    const float* ranges[3];
    int channels[3];

	float threshold;           // decision threshold
	cv::Mat histogram;         // histogram can be sparse 输入直方图
	cv::SparseMat shistogram;  // or not
	bool isSparse;

  public:

	ContentFinder() : threshold(0.1f), isSparse(false) {

		// in this class,
		// all channels have the same range
		ranges[0]= hranges;  
		ranges[1]= hranges; 
		ranges[2]= hranges; 
	}
   
	// Sets the threshold on histogram values [0,1]
	void setThreshold(float t) {

		threshold= t;
	}

	// Gets the threshold
	float getThreshold() {

		return threshold;
	}

	// Sets the reference histogram
	void setHistogram(const cv::Mat& h) {

		isSparse= false;
		cv::normalize(h,histogram,1.0);
	}

	// Sets the reference histogram
	void setHistogram(const cv::SparseMat& h) {

		isSparse= true;
		cv::normalize(h,shistogram,1.0,cv::NORM_L2);
	}

	// All channels used, with range [0,256]
	cv::Mat find(const cv::Mat& image) {

		cv::Mat result;

		hranges[0]= 0.0;	// default range [0,256]
		hranges[1]= 256.0;
		channels[0]= 0;		// the three channels 
		channels[1]= 1; 
		channels[2]= 2; 

		return find(image, hranges[0], hranges[1], channels);
	}

	// Finds the pixels belonging to the histogram
	cv::Mat find(const cv::Mat& image, float minValue, float maxValue, int *channels) {

		cv::Mat result;

		hranges[0]= minValue;
		hranges[1]= maxValue;

		if (isSparse) { // call the right function based on histogram type

		   for (int i=0; i<shistogram.dims(); i++)
			  this->channels[i]= channels[i];

		   cv::calcBackProject(&image,
                      1,            // we only use one image at a time
                      channels,     // vector specifying what histogram dimensions belong to what image channels
                      shistogram,   // the histogram we are using
                      result,       // the resulting back projection image
                      ranges,       // the range of values, for each dimension
                      255.0         // the scaling factor is chosen such that a histogram value of 1 maps to 255
		   );

		} else {

		   for (int i=0; i<histogram.dims; i++)
			  this->channels[i]= channels[i];
//某对象的this指针，指向被调用函数所在的对象，此处对象为ContentFinder类
		   //this->channels[i]即ContentFinder类的私有成员channels[3]
		   //对ContentFinder类各成员的访问均通过this进行
		   cv::calcBackProject(&image,
                      1,            // we only use one image at a time
                      channels,     // 向量表示哪个直方图维度属于哪个图像通道
                      histogram,    // 用到的直方图
                      result,       // 反向投影的图像
                      ranges,       // 每个维度值的范围
                      255.0         // 选用的换算系数
		   );
		}
        // Threshold back projection to obtain a binary image阈值分割反向投影图像得到二值图
		if (threshold>0.0)// 设置的阈值>0时，才进行阈值分割
			cv::threshold(result, result, 255.0*threshold, 255.0, cv::THRESH_BINARY);
		return result;
	}
};
#endif

3.提取感兴趣区域ROI的直方图

以上定义了直方图操作方法和反向投影直方图检测方法。下面将在main函数中调用上述方法来做反向投影直方图检测的实验。

要查找图像中特定的内容（例如在下图中检测出天空中的云彩），首先要选择一个包括所需样本的兴趣区域，即画一个矩形框选出ROI。

（1）读取图像

补充：imread(const string& filename, int flags=1)函数读取图像色彩空间参数flags:

enum
{
    // 8bit, color or not
    IMREAD_UNCHANGED  =-1,
    // 8bit, gray
    IMREAD_GRAYSCALE  =0,
    // color
    IMREAD_COLOR      =1,
    // any depth,
    IMREAD_ANYDEPTH   =2,
    // any color
    IMREAD_ANYCOLOR   =4
};

在main函数中运行程序主干。以灰度图的形式读取图像。

// Read input image
    cv::Mat image= cv::imread("f:\\\\images\\\\waves.jpg",0);// /灰度图方式读取图像
	if (!image.data)
		return 0; 

（2）设置ROI

设置读取图像中的云层区域为ROI。

	// define image ROI
	cv::Mat imageROI;
	imageROI= image(cv::Rect(406,146,30,24)); // Cloud  region

// Find histogram of reference
	Histogram1D h;
	cv::Mat hist= h.getHistogram(imageROI);
	cv::namedWindow("Reference Hist");
	cv::imshow("Reference Hist",h.getHistogramImage(imageROI));
	waitKey(0);

4.调用反向投影直方图检测灰度图像

前文定义了1维直方图的操作方法Histogram1D类，定义了计算反向投影直方图的方法ContentFinder类，在main函数中读取了图像、设置好ROI，获取了ROI直方图后，就可以调用ContentFinder类来计算反向投影了。

注意为了能将处理过程的中间结果显示出来，设置阈值为-0.1<0，用ContentFinder类在find()方法中当设置的阈值<0时不进行阈值分割，因此输出的结果为灰度图像，便于显示观察。下面的代码先使用find()获取反向投影结果result1，然后使用result1.convertTo(tmp,CV_8U,-1.0,255.0);将result1转换为反向投影图像tmp。

// Create the content finder
	ContentFinder finder;

	// set histogram to be back-projected
	finder.setHistogram(hist);//将ROI直方图hist传入反向投影计算类ContentFinder finder
	finder.setThreshold(-1.0f);// /设置阈值为-0.1，float型。阈值<0，不使用find方法中的阈值分割，输出result为灰度图 

	// Get back-projection
	cv::Mat result1;
	result1= finder.find(image);

	// Create negative image and display result
	cv::Mat tmp;
	result1.convertTo(tmp,CV_8U,-1.0,255.0);
	cv::namedWindow("Backprojection result");
	cv::imshow("Backprojection result",tmp);

为了搞清楚convertTo()的操作，我们用如下代码扫描图像像素的方法读取result1中一个15×15矩形区域的值在屏幕上输出：

if(1)//for debug
	{
		cv::namedWindow("Backprojection result1");
		cv::imshow("Backprojection result1",result1);
		cout<<"pixel vaule of result1(421,141)~(436,156),15x15 region"<<endl<<endl;
		for ( int j = 140;j<155;j++)//从第j=101行开始
		{
			uchar* data= result1.ptr<uchar>(j); 
			for (int i=420;i<435;i++)
			{
				cout<<(int)data[i]<<" ";//打印输出第j行的所有数据
			}
			cout<<endl;//输出1行后换行
		}
	}

result1图像该区域的像素值输出结果如下图所示

图像result1从点(421,141)到点(436,156)的矩形区域的像素值

使用使用上述方法输出tmp图像同样位置的15×15矩形区域的像素值如下图所示。比较转换前后的图像，发现对源图中的像素逐点I(x,y)做了反向操作I(x,y)=255-I(x,y).

图像tmp从点(421,141)到点(436,156)的矩形区域的像素值

显示出反向前后的全图如下2图所示。

阈值分割后的反向投影检测结果

5.调用反向投影直方图检测彩色图像

从上图检测结果中可以看到，本来是检测云彩区域，但得到的检测结果还包括了底部许多沙滩区域和部分海水区域。一个改进的方法就是使用彩色图像的直方图来进行反向投影检测。使用彩色信息检测需要用到前面提到的3维彩色直方图操作ColorHistogram类。在main函数中通过如下代码来进行彩色图像的反向投影直方图检测。

// Load color image
	ColorHistogram hc; //声明一个彩色直方图类hc用于下面获取ROI直方图的操作
    cv::Mat color= cv::imread("f:\\\\images\\\\waves.jpg");//这里要使用本机的图像路径

	// extract region of interest
	imageROI= color(cv::Rect(406,154,30,24)); //设置检测云彩区域

	// Get 3D colour histogram (8 bins per channel)
	hc.setSize(8); // 8x8x8 降低bin的数量使用8,16,64分别计算，结果如下图所示。
	cv::Mat shist= hc.getHistogram(imageROI);

	// set histogram to be back-projected
	finder.setHistogram(shist);
	finder.setThreshold(0.05f);

	// Get back-projection of color histogram
	result1= finder.find(color);
	cv::namedWindow("Color Detection Result");
	cv::imshow("Color Detection Result",result1);

下图是输入的彩色图像与ROI检测区域

彩色图像与选择的ROI区域

在程序中，可以通过hc.setSize()来设置直方图中每个通道bin个数，本实验中，选择了8/16/64三个参数，检测结果如下图所示。从检测结果可以看出来，使用彩色直方图后，检测效果大大优于灰度直方图，并且连水中云的倒影都能检测出来。而设置的bin个数越多，检测精度越高。并且本例中使用16个bin得到的检测结果优于8和64个bin，因为精度过高将排除掉部分与ROI区域类似的云彩区域，所以bin不是越多越好，合适的bin的个数需要根据实际情况而定。

彩色图检测结果（16个bin，基本检测出了全部的云彩）

彩色图检测结果（64个bin，有部分云彩漏检）

// Convert to Lab space
	cv::Mat lab;
	cv::cvtColor(color, lab, CV_BGR2Lab);

	// Convert to HSV space
	cv::Mat hsv;
	cv::cvtColor(color, hsv, CV_BGR2HSV);