精华内容
下载资源
问答
  • 2019-07-05 23:22:01

    初衷

    果然halcon用顺手了,人就变懒了,正好有项目需要自己写个形状匹配的程序,就拿来练练手,程序不是很复杂,速度上感觉和halcon里面find_scaled_shape_model还是有差距,目前也不知道如何进一步改进,暂时就先这样了。大佬们如果有什么改进的想法可以评论一下,这样才能不断进步。

    思路

    主要的想法还是基于散点的重合度,虽然openCV自带matchTemplate和matchShape函数,但其内部是进行遍历像素矩计算,由于检测图中会有干扰边缘的存在,所以这些函数对于形状匹配(局部像素点集匹配)来说是不适合的。当然这些函数运行速度还是非常快的,目前我也不知道用什么办法可以达到同数量级的运算速度。

    算法实现过程中主要遇到了几个问题:

    1.生成旋转模版的问题

    openCV没有自带的旋转函数,需要通过getRotationMatrix2D生成旋转矩阵,再用warpAffine进行仿射变换。但是warpAffine生成的图形会自动进行裁剪,会带来一些问题,因此需要稍微修改(详情见imrotate函数)。

    2.遍历边缘点还是进行卷积的问题

    由于openCV和halcon对于取像素点的操作都特别慢,所以对比之下只能采取用filter2D进行遍历,不得不说openCV自带的硬件加速库还是高效。

    3.遍历运算速度过慢

    openCV一些函数的源代码在内部还是进行了hal库的调用,用户没有办法修改,因此没有办法在函数内部进行算法修改。尝试用ppl库进行并行加速,效果还可以。

    程序使用说明

    (1)优先修改金字塔层数PyrLevel,层数越大,搜索速度越快,主函数中有缩放后的边缘图片,根据边缘是否清晰来修改PyrLevel和Canny函数的阈值,保证待检测图和模版图片的边缘足够清晰完整。
    (2)如果找不到,适当减小minScore,越接近0,搜索相似容忍度越大,但运行速度同时也越慢。

    函数片段

    具体代码工程见
    https://download.csdn.net/download/sillykog/10727775

    /*************************************************
    Function:       //	CreateScaledShapeModel
    Description:    //	创建模版集
    Input:          //	Template:模版图片的边缘图
    					PyrLevel:金字塔缩小层数
    					AngleStart,AngleExtent,AngleStep:旋转角度候选参数
    					ScaleMin,ScaleMax,ScaleStep:缩放比例候选参数
    Output:         //	pModelImageSet,pModelPointSet,pScaleSet,pAngleSet:模版集指针
    Return:         //  无
    Others:         //	
    *************************************************/
    void CreateScaledShapeModel(Mat Template, int PyrLevel, int AngleStart, int AngleExtent, int AngleStep, float ScaleMin, float ScaleMax, float ScaleStep, \
    	vector<Mat>* pModelImageSet, vector<int>* pModelPointSet, vector<float>* pScaleSet, vector<float>* pAngleSet)
    {
    	vector<Mat> ModelImageSet;
    	vector<int> ModelPointSet;
    	vector<float> AngleSet;
    	vector<float> ScaleSet;
    	while (ScaleMin <= ScaleMax)
    	{
    		cout << ScaleMax << endl;
    		ScaleSet.push_back(ScaleMax);
    		ScaleMax -= ScaleStep;
    	}
    	while (AngleStart <= AngleExtent)
    	{
    		cout << AngleExtent << endl;
    		AngleSet.push_back(AngleExtent);
    		AngleExtent -= AngleStep;
    	}
    	//模版生成	
    	for (int level = 0; level <= PyrLevel; level++)
    	{
    		Mat pyrmodelImage = Template;
    		for (int i = 0; i < level; i++)
    		{
    			pyrDown(pyrmodelImage, pyrmodelImage);
    		}
    		//缩放
    		for (int i = 0; i < ScaleSet.size(); i++)
    		{
    			Mat scaleImage;
    			resize(pyrmodelImage, scaleImage, Size(round(pyrmodelImage.cols*ScaleSet[i]), round(pyrmodelImage.cols*ScaleSet[i])), 0, 0, INTER_LINEAR);
    			//旋转
    			for (int j = 0; j < AngleSet.size(); j++)
    			{
    				Mat rotateImage;
    				imrotate(scaleImage, rotateImage, AngleSet[j]);
    				//threshold(rotateImage, rotateImage, 1, 255, 0);
    				Canny(rotateImage, rotateImage, 50, 100, 3, false);
    				rotateImage /= 255;
    				//imshow("旋转", rotateImage);
    				//imwrite("旋转.jpg", rotateImage);
    				//waitKey(0);
    				ModelImageSet.push_back(rotateImage);
    				int pointNum = 0;
    				for (int i = 0; i < rotateImage.cols; i++)
    				{
    					for (int j = 0; j < rotateImage.rows; j++)
    					{
    						if (rotateImage.at<uchar>(Point(i, j)) != 0)
    							pointNum++;
    					}
    				}
    				ModelPointSet.push_back(pointNum);
    				rotateImage.release();
    			}
    			scaleImage.release();
    		}
    	}
    	*pModelImageSet = ModelImageSet;
    	*pModelPointSet = ModelPointSet;
    	*pAngleSet = AngleSet;
    	*pScaleSet = ScaleSet;
    }
    
    
    
    /*************************************************
    Function:       //	FindScaledShapeModel
    Description:    //	在一张图片中搜索与模版相似的图形
    Input:          //	Image:待检测图片
    					ModelImageSet,ModelPointSet,ScaleSet,AngleSet:模版集
    					PyrLevel:金字塔缩小层数
    					MinScore:筛选相似度阈值		
    Output:         //	pRow,pCol,pScale,pAngle,pScore:输出匹配到的元素参数集合的指针
    Return:         //  无
    Others:         //	使用该函数前需要先调用CreateScaledShapeModel
    *************************************************/
    void FindScaledShapeModel(Mat Image, vector<Mat> ModelImageSet, vector<int> ModelPointSet, vector<float> ScaleSet, vector<float> AngleSet, int PyrLevel, float MinScore,\
    	vector<int>* pRow, vector<int> * pCol, vector<float>* pScale, vector<float>* pAngle, vector<float>* pScore)
    {
    	mutex mt;
    	Mat modelImage = ModelImageSet[0];
    	vector<int> Row;
    	vector<int> Col;
    	vector<float> Scale;
    	vector<float> Angle;
    	vector<float> Score;
    	bool findFlag = false;
    	//金字塔分层匹配
    	for (int level = PyrLevel; !findFlag && level >= PyrLevel; level--)
    	{		
    		Mat pyrsrcImage = Image;
    		for (int i = 0; i < level; i++)
    		{
    			pyrDown(pyrsrcImage, pyrsrcImage);
    		}
    
    		int kernSize = floor(sqrt(min(pyrsrcImage.rows / 100, pyrsrcImage.cols / 100)));		
    		Mat kern = Mat::ones(2 * kernSize + 1, 2 * kernSize + 1, CV_8U);
    			
    		Mat blurImage;
    		filter2D(pyrsrcImage, blurImage, pyrsrcImage.depth(), kern);		
    		//imshow("糊化原图", blurImage);
    		//moveWindow("糊化原图", 0, 0);
    		//waitKey(10);
    
    		Mat tempblurImage;
    		blurImage.convertTo(tempblurImage, CV_8U);
    		tempblurImage /= 255;
    		int parallelnum = ScaleSet.size() *AngleSet.size();
    
    		parallel_for(0, parallelnum, [&](int k)	
    		{
    			Mat scoreImage(tempblurImage.size(), CV_16U);
    			Mat tempmodelImage = ModelImageSet.at(ModelImageSet.size()- 1 - k);
    			int temppointNum = ModelPointSet.at(ModelPointSet.size() - 1 - k);
    			float max_score = 0;		
    			/*imshow("模版", tempmodelImage);
    			resizeWindow("模版", tempmodelImage.rows, tempmodelImage.cols);		
    			moveWindow("模版", blurImage.cols,0);
    			waitKey(10);*/
    			//double start = static_cast<double>(getTickCount());
    
    			filter2D(tempblurImage, scoreImage, scoreImage.depth(), tempmodelImage);
    			//double time = ((double)getTickCount() - start) / getTickFrequency();
    			//cout << "所用时间为:" << time << "秒" << endl;
    			mt.lock();
    			while (1)
    			{
    				double v_min, v_max;
    				int idx_min[2] = { 255,255 }, idx_max[2] = { 255, 255 };
    				minMaxIdx(scoreImage, &v_min, &v_max, idx_min, idx_max);
    				scoreImage.at<ushort>(idx_max[0], idx_max[1]) = 0;
    				
    				max_score = v_max / temppointNum;
    
    				//cout << "第" << level << "层" << "第" << k + 1 << "个成绩:" << max_score << endl;
    				if (max_score > MinScore)
    				{
    					float scale = ScaleSet[ScaleSet.size() - 1 - (k) / AngleSet.size()];
    					float angle = AngleSet[AngleSet.size() - 1 - (k) % AngleSet.size()];
    					//int selectx = (idx_max[1] - (tempmodelImage.cols - 1) / 2)*pow(2, level);
    					//int selecty = (idx_max[0] - (tempmodelImage.rows - 1) / 2)*pow(2, level);
    					//int pyrselectx = idx_max[1] - tempmodelImage.cols / 2;
    					//int pyrselecty = idx_max[0] - tempmodelImage.rows / 2;
    					Row.push_back(idx_max[1] * pow(2, level));
    					Col.push_back(idx_max[0] * pow(2, level));
    					Scale.push_back(scale);
    					Angle.push_back(angle);
    					Score.push_back(max_score);
    					//cout << Point(selectx, selecty) << " " << Point(pyrselectx, pyrselecty) << endl;
    					//rectangle(blurImage, Rect(pyrselectx, pyrselecty, tempmodelImage.cols, tempmodelImage.rows), 255, 2, 8);
    					//rectangle(Image, Rect(selectx, selecty, modelImage.cols / scale, modelImage.rows / scale), 255, 2, 8);
    					//imshow("缩放图位置", blurImage);
    					//imshow("原图位置", Image);
    					//waitKey(0);
    					//findFlag = true;				
    				}		
    				else break;
    			}
    			tempmodelImage.release();
    			scoreImage.release();
    			mt.unlock();
    		}
    		);
    		for (int m = 0; m < Row.size() ; m++)
    		{
    			for (int n = m+1; n < Row.size() ; n++)
    			{
    				if (abs(Row[n] - Row[m])<modelImage.rows*0.3  && abs(Col[n] - Col[m])< modelImage.cols*03)
    				{
    					if (Score[n] < Score[m])
    					{
    						Row.erase(Row.begin() + n);
    						Col.erase(Col.begin() + n);
    						Scale.erase(Scale.begin() + n);
    						Angle.erase(Angle.begin() + n);
    						Score.erase(Score.begin() + n);
    						n--;
    					}
    					else if (Score[n] >= Score[m])
    					{
    						swap(Row[m], Row[n]);
    						swap(Col[m], Col[n]);
    						swap(Scale[m], Scale[n]); 
    						swap(Angle[m], Angle[n]); 
    						swap(Score[m], Score[n]);
    						Row.erase(Row.begin() + n);
    						Col.erase(Col.begin() + n);
    						Scale.erase(Scale.begin() + n);
    						Angle.erase(Angle.begin() + n);
    						Score.erase(Score.begin() + n);
    						n = m ;
    					}
    				
    				}
    			}
    		}
    		*pRow = Row;
    		*pCol = Col;
    		*pScale = Scale;
    		*pAngle = Angle;
    		*pScore = Score;
    	}
    }
    void imrotate(Mat& img, Mat& newIm, double angle) 
    {
    	int heightNew = int(img.cols*fabs(sin(angle*3.14 / 180)) + img.rows * fabs(cos(angle*3.14 / 180)));
    	int widthNew = int(img.rows*fabs(sin(angle*3.14 / 180)) + img.cols * fabs(cos(angle*3.14 / 180)));
    	int len = max(img.cols, img.rows);
    	Point2f pt(len / 2., len / 2.);
    	Mat r = getRotationMatrix2D(pt, angle, 1.0);
    	r.at<double>(0, 2) += (widthNew - img.cols) / 2;
    	r.at<double>(1, 2) += (heightNew - img.rows) / 2;	
    	warpAffine(img, newIm, r, Size(widthNew, heightNew), INTER_LINEAR, BORDER_REPLICATE);
    }
    
    
    更多相关内容
  • 用于形状匹配和点对应的形状上下文。 忠实地实现了Belongie,Malik和Puzicha的“形状上下文”。 基于玫瑰直方图之间的卡方距离。 入门 这个小项目是Belongie,Malik和Puzicha的“形状上下文”的忠实实现。 它们可...
  • QT+Opencv相机标定获取去畸变矫正源码 QT+OpenCV 9点标定源码 QT+OpenCV 圆拟合源码 QT+Opencv模板匹配源码 QT+Opencv做匹配+旋转+缩放 QT+Opencv做形状匹配源码 QT+Opencv特征匹配源码 QT+Opencv线点匹配源码
  • 针对现有形状匹配算法匹配速度慢、可靠性差的问题,提出了一种基于轮廓矢量化的形状匹配算法。将轮廓曲线点集离散化,使用多组向量对轮廓线性进行逼近。匹配时,结合线段在几何中的匹配方法,通过加权求和,计算源...
  • opencv 模板匹配&&形状匹配

    千次阅读 2022-04-28 20:17:23
    找圆垫子1.1 得到模板1.2 形状匹配2. 找瓜子 这是第四次作业要求 所以今天就趁机会讲讲模板匹配,正好之前的项目有一部分重要工作就是和模板匹配紧密相关,对于今天作业来说,之前的项目难度更大,因为涉及到许多...


    这是第四次作业要求

    image-20220428142137548

    所以今天就趁机会讲讲模板匹配,正好之前的项目有一部分重要工作就是和模板匹配紧密相关,对于今天作业来说,之前的项目难度更大,因为涉及到许多要考虑的因素,还要考虑效率实时性等问题。太详细的我也不方便展开,下面先看看之前的效果

    17_50_13_9

    17_50_13_9

    当然也有其他的车型,视频就不放了直接上结果图

    18_00_13_1

    17_50_13_3

    17_50_13_8

    现在就正式开始解决今天这个作业

    1. 找圆垫子

    9_1

    在上面视频拼车里模板都是在每一帧中取得,并且为了模板匹配能匹配得上还要预处理每一帧图片,调整模板区域、大小以及匹配方法。但是在拼车中我们只需要匹配一次,而作业中明显是需要多目标匹配的,并且待匹配目标的角度差异很大,在这种情况下使用opencv中的matchTemplate()效果肯定很差,至于为什么这可以从原理上解释。

    image-20220428193221641

    有兴趣的可以去看看相关资料,我就根据上图简单来说说模板匹配流程

    1. 得到原图和模板
    2. 滑动模板,与原图作对比
    3. 比较依据也就是评价指标如下图

    image-20220428193436767

    平方差方法匹配值越小越好,相关性方法匹配值越大越好,这样就会得到一个result矩阵,其中有最亮/最暗的地方就是我们匹配到的地方。利用minMaxLoc()找到这些地方,那也就代表找到了待匹配目标。

    但是,我在上面也说过了,这些待匹配目标和模板的角度差异特别大,你用模板匹配最多也只能匹配出和模板一样大小、方向的待匹配目标。或者利用一堆复杂的操作,把模板图片进行角度变换搞一个任意角度都有的模板组然后再一一匹配,这个倒是可行,但是效率上不敢说。所以干脆换一种思路,不用模板匹配而是用形状去匹配。

    匹配嘛就是根据相似性进行比较,上面已经说了因为角度导致模板匹配失效,那就得想办法排除掉角度干扰,最容易想到的就是利用某些角度不变特性的算子,所以很自然就想到HuMoments,关于Hu不变矩可以看看opencv的解释

    image-20220428195311262

    利用归一化中心矩构造了不变特征矩,这个对于平移,旋转都可以保持不变性,这不就是我们需要的吗?然后利用这个开始实现作业。

    image-20220428195703459

    image-20220428195730804

    这是opencv中形状匹配的函数以及参数,我之前也没用过所以放在这一起学习一下。

    1.1 得到模板

    我们只有原图,并没有模板图。自己想去截图但是待匹配图之间挨得比较近,截下来的模板图不太容易找轮廓,因此直接在待匹配图操作,找到其中一个的轮廓作为模板。

    vector<Point> TemplateContours(Mat img_template)
    {
    	//灰度化
    	Mat gray_img_template;
    	cvtColor(img_template, gray_img_template, COLOR_BGR2GRAY);
    
    	//二值化
    	Mat thresh_img_template;
    	threshold(gray_img_template, thresh_img_template, 0, 255, THRESH_OTSU);
    
    	imshow("b", thresh_img_template);
    
    	//膨胀处理
    	Mat ellipse = getStructuringElement(MORPH_ELLIPSE, Size(15, 15));
    	morphologyEx(thresh_img_template, thresh_img_template, MORPH_OPEN, ellipse, Point(-1, -1), 1);
    
    	//寻找边界
    	vector<vector<Point>> contours_template;
    	vector<Vec4i> hierarchy;
    	findContours(thresh_img_template, contours_template, hierarchy,RETR_LIST, CHAIN_APPROX_SIMPLE, Point());
    
    	//绘制边界,查看是否是我们的目标轮廓
    	drawContours(img_template, contours_template, 0, Scalar(0, 0, 255), 1, 8, hierarchy);
    	imshow("img_template", img_template);
    	waitKey(0);
    	return contours_template[0];
    }
    

    image-20220428200330063

    1.2 形状匹配

    那个红色的圈就是每个圆形垫子的共有特征,然后根据这个圈作为模板轮廓对图像进行形状匹配,并根据设定的阈值排除掉不是相同形状物体的干扰。

    void ShapeTemplateMatch(Mat image, vector<Point> imgTemplatecontours, double thresh)
    {
    	//灰度化
    	Mat gray_img;
    	cvtColor(image, gray_img, COLOR_BGR2GRAY);
    
    	//二值化
    	Mat thresh_img;
    	threshold(gray_img, thresh_img, 0, 255, THRESH_OTSU);
    
    	//寻找边界
    	vector<vector<Point>> contours_img;
    	vector<Vec4i> hierarchy;
    	findContours(thresh_img, contours_img, hierarchy, RETR_LIST, CHAIN_APPROX_NONE, Point());
    
    	//根据形状模板进行匹配
    	int min_pos = -1;
    
    	for (int i = 0; i < contours_img.size(); i++)
    	{
    		//计算轮廓面积,筛选掉一些没必要的小轮廓
    		// 根据计算面积以及可视化结果确定下面的阈值
    		//cout << contourArea(contours_img[i]) << endl;
    		if (contourArea(contours_img[i]) > 120)
    		{
    			//进行形状匹配,形状越相似value越小
    			double value = matchShapes(contours_img[i], imgTemplatecontours, CONTOURS_MATCH_I3, 0.0);
    
    			//查看value以及边界图得到thresh阈值
    			//cout << value << endl;
    
    			//将匹配分值与设定阈值比较 
    			if (value < thresh)
    			{
    				min_pos = i;
    				//绘制目标边界
    				drawContours(image, contours_img, min_pos, Scalar(0, 255,0), 5, 8, hierarchy, 0);
    				imshow("1", image);
    			}
    		}
    	}
    	waitKey(0);
    }
    

    image-20220428200455886

    2. 找瓜子

    9_2

    和上面的函数一样,只是改了一下阈值

    int main() {
    	Mat src = imread("D:/9_2.jpg");
    	resize(src, src, Size(src.cols / 4, src.rows / 4));
    	vector<Point> tempContours =TemplateContours(src.clone());
    
    	//圆垫阈值0.03,瓜子阈值0.2
    	ShapeTemplateMatch(src.clone(),tempContours, 0.2);
    	return 0;
    }
    

    image-20220428201530908

    image-20220428201554768

    展开全文
  • 在本文中,我们提出了一种基于视觉感知的形状分解和匹配的新方法,该方法使用多个尺度下的凸角。 我们的方法由多尺度形状分解过程和轮廓线段匹配过程组成。 首先,形态学运算用于更稳定和视觉感知的拐角,以在多个...
  • 为了对轮毂型号进行识别, 提出一种基于形状匹配及纹理筛选的轮型识别算法。首先, 确定一个轮辐形状为标准模板并得出其边缘图, 把模板作为移动窗口在待识别轮毂图片中移动, 逐一计算模板到轮毂图片各感兴趣区域( ...
  • DFT的matlab源代码一小段代码说明形状上下文匹配。 实现的代码通常非常简单,可以在中遵循。 我们使用OpenCV库执行大多数I / O。 我对这些点进行装箱,然后使用辅助库进行加权二分匹配。 为了方便起见,我从文件中...
  • DFT的matlab源代码一小段代码说明形状上下文匹配。 实现的代码通常非常简单,可以在中遵循。 我们使用OpenCV库执行大多数I / O。 我对这些点进行装箱,然后使用辅助库进行加权二分匹配。 为了方便起见,我从文件中...
  • 该方法通过分析图像特征描述符特性,考虑静脉图像的全局信息,确定提取具有平移、旋转、比例和仿射不变性的边界和区域的形状特征进行粗匹配;为了考虑局部特征使形状描述符具有良好的易分辨能力,依据小波矩具有小波多...
  • HALCON基于形状匹配详解

    千次阅读 2019-05-06 13:43:51
    去年有过一段时间的集中学习,做了许多的练习和实验,并对基于HDevelop的形状匹配算法的参数优化进行了研究,写了一篇《基于HDevelop的形状匹配算法参数的优化研究》文章,总结了在形状匹配过程中哪些参数影响到模板...

    HALCON基于形状的模板匹配详细说明

    很早就想总结一下前段时间学习HALCON的心得,但由于其他的事情总是抽不出时间。去年有过一段时间的集中学习,做了许多的练习和实验,并对基于HDevelop的形状匹配算法的参数优化进行了研究,写了一篇《基于HDevelop的形状匹配算法参数的优化研究》文章,总结了在形状匹配过程中哪些参数影响到模板的搜索和匹配,又如何来协调这些参数来加快匹配过程,提高匹配的精度,这篇paper放到了中国论文在线了,需要可以去下载。

    德国MVTec公司开发的HALCON机器视觉开发软件,提供了许多的功能,在这里我主要学习和研究了其中的形状匹配的算法和流程。HDevelop开发环境中提供的匹配的方法主要有三种,即Component-Based、Gray-Value-Based、Shape-Based,分别是基于组件(或成分、元素)的匹配,基于灰度值的匹配和基于形状的匹配。这三种匹配的方法各具特点,分别适用于不同的图像特征,但都有创建模板和寻找模板的相同过程。这三种方法里面,我主要就第三种-基于形状的匹配,做了许多的实验,因此也做了基于形状匹配的物体识别,基于形状匹配的视频对象分割和基于形状匹配的视频对象跟踪这些研究,从中取得较好的效果,简化了用其他工具,比如VC++来开发的过程。在VC下往往针对不同的图像格式,就会弄的很头疼,更不用说编写图像特征提取、模板建立和搜寻模板的代码呢,我想其中间过程会很复杂,效果也不一定会显著。下面我就具体地谈谈基于HALCON的形状匹配算法的研究和心得总结。

    1.   Shape-Based matching的基本流程
      

    HALCON提供的基于形状匹配的算法主要是针对感兴趣的小区域来建立模板,对整个图像建立模板也可以,但这样除非是对象在整个图像中所占比例很大,比如像视频会议中人体上半身这样的图像,我在后面的视频对象跟踪实验中就是针对整个图像的,这往往也是要牺牲匹配速度的,这个后面再讲。基本流程是这样的,如下所示:
    ⑴ 首先确定出ROI的矩形区域,这里只需要确定矩形的左上点和右下点的坐标即可,gen_rectangle1()这个函数就会帮助你生成一个矩形,利用area_center()找到这个矩形的中心;

    ⑵ 然后需要从图像中获取这个矩形区域的图像,reduce_domain()会得到这个ROI;这之后就可以对这个矩形建立模板,而在建立模板之前,可以先对这个区域进行一些处理,方便以后的建模,比如阈值分割,数学形态学的一些处理等等;

    ⑶ 接下来就可以利用create_shape_model()来创建模板了,这个函数有许多参数,其中金字塔的级数由Numlevels指定,值越大则找到物体的时间越少,AngleStart和AngleExtent决定可能的旋转范围,AngleStep指定角度范围搜索的步长;这里需要提醒的是,在任何情况下,模板应适合主内存,搜索时间会缩短。对特别大的模板,用Optimization来减少模板点的数量是很有用的;MinConstrast将模板从图像的噪声中分离出来,如果灰度值的波动范围是10,则MinConstrast应当设为10;Metric参数决定模板识别的条件,如果设为’use_polarity’,则图像中的物体和模板必须有相同的对比度;创建好模板后,这时还需要监视模板,用inspect_shape_model()来完成,它检查参数的适用性,还能帮助找到合适的参数;另外,还需要获得这个模板的轮廓,用于后面的匹配,get_shape_model_contours()则会很容易的帮我们找到模板的轮廓;

    ⑷ 创建好模板后,就可以打开另一幅图像,来进行模板匹配了。这个过程也就是在新图像中寻找与模板匹配的图像部分,这部分的工作就由函数find_shape_model()来承担了,它也拥有许多的参数,这些参数都影响着寻找模板的速度和精度。这个的功能就是在一幅图中找出最佳匹配的模板,返回一个模板实例的长、宽和旋转角度。其中参数SubPixel决定是否精确到亚像素级,设为’interpolation’,则会精确到,这个模式不会占用太多时间,若需要更精确,则可设为’least_square’,’lease_square_high’,但这样会增加额外的时间,因此,这需要在时间和精度上作个折中,需要和实际联系起来。比较重要的两个参数是MinSocre和Greediness,前一个用来分析模板的旋转对称和它们之间的相似度,值越大,则越相似,后一个是搜索贪婪度,这个值在很大程度上影响着搜索速度,若为0,则为启发式搜索,很耗时,若为1,则为不安全搜索,但最快。在大多数情况下,在能够匹配的情况下,尽可能的增大其值。
    ⑸ 找到之后,还需要对其进行转化,使之能够显示,这两个函数vector_angle_to_rigid()和affine_trans_contour_xld()在这里就起这个作用。前一个是从一个点和角度计算一个刚体仿射变换,这个函数从匹配函数的结果中对构造一个刚体仿射变换很有用,把参考图像变为当前图像。

    其详细的流程图和中间参数,如下图所示:(无法上传)

    1.   基于形状匹配的参数关系与优化
      
      在HALCON的说明资料里讲到了这些参数的作用以及关系,在上面提到的文章中也作了介绍,这里主要是重复说明一下这些参数的作用,再强调一下它们影响匹配速度的程度;

    在为了提高速度而设置参数之前,有必要找出那些在所有测试图像中匹配成功的设置,这时需考虑以下情况:

    ① 必须保证物体在图像边缘处截断,也就是保证轮廓的清晰,这些可以通过形态学的一些方法来处理;

    ② 如果Greediness值设的太高,就找不到其中一些可见物体,这时最后将其设为0来执行完全搜索;
    ③ 物体是否有封闭区域,如果要求物体在任何状态下都能被识别,则应减小MinScore值;

    ④ 判断在金字塔最高级上的匹配是否失败,可以通过find_shape_model()减小NumLevels值来测试;

    ⑤ 物体是否具有较低的对比度,如果要求物体在任何状态下都能被识别,则应减小MinContrast值;

    ⑥ 判断是否全局地或者局部地转化对比度极性,如果需要在任何状态下都能被识别,则应给参数Metric设置一个合适的值;
    ⑦ 物体是否与物体的其他实例重叠,如果需要在任何状态下都能识别物体,则应增加MaxOverlap值;

    ⑧ 判断是否在相同物体上找到多个匹配值,如果物体几乎是对称的,则需要控制旋转范围;

    如何加快搜索匹配,需要在这些参数中进行合理的搭配,有以下方法可以参考:

    ① 只要匹配成功,则尽可能增加参数MinScore的值;
    ② 增加Greediness值直到匹配失败,同时在需要时减小MinScore值;

    ③ 如果有可能,在创建模板时使用一个大的NumLevels,即将图像多分几个金字塔级;

    ④ 限定允许的旋转范围和大小范围,在调用find_shape_model()时调整相应的参数;

    ⑤ 尽量限定搜索ROI的区域;
    除上面介绍的以外,在保证能够匹配的情况下,尽可能的增大Greediness的值,因为在后面的实验中,用模板匹配进行视频对象跟踪的过程中,这个值在很大程度上影响到匹配的速度。

    当然这些方法都需要跟实际联系起来,不同图像在匹配过程中也会有不同的匹配效果,在具体到某些应用,不同的硬件设施也会对这个匹配算法提出新的要求,所以需要不断地去尝试。在接下来我会结合自己做的具体的实验来如何利用HALCON来进行实验,主要是在视频对象分割和视频对象的跟踪方面。待续…………_

    例子:区分硬币 Example: solution_guide/basics/matching_coins.hdev 见附件trainmodel.hdev

    点击查看原图

    create_shape_model(

    const Hobject& Template , //reduce_domain后的模板图像

    Hlong NumLevels, //金字塔的层数,可设为“auto”或0—10的整数

    Double AngleStart, //模板旋转的起始角度

    Double AngleExtent, //模板旋转角度范围, >=0

    Double AngleStep, //旋转角度的步长, >=0 and <=pi/16

    const char* Optimization, //设置模板优化和模板创建方法

    const char* Metric, //匹配方法设置

    Hlong Contrast, //设置对比度

    Hlong MinContrast , //设置最小对比度

    Hlong* ModelID ) //输出模板句柄

    进一步分析:

    NumLevels越大,找到匹配使用的时间就越小。另外必须保证最高层的图像具有足够的信息(至少四个点)。可以通过inspect_shape_model函数查看设置的结果。如果最高层金字塔的消息太少,算法内部会自动减少金字塔层数,如果最底层金字塔的信息太少,函数就会报错。如果设为auto,算法会自动计算金字塔的层数,我们可以通过get_shape_model_params函数查看金字塔的层数。如果金字塔的层数太大,模板不容易识别出来,这是需要将find_shape_model函数中MinScore和Greediness参数设置的低一些。如果金字塔层数太少找到模板的时间会增加。可以先使用inspect_shape_model函数的输出结果来选择一个较好的金字塔层数。

    参数AngleStart、AngleExtent定义了模板可能发生旋转的范围。注意模板在find_shape_model函数中只能找到这个范围内的匹配。参数AngleStep定义了旋转角度范围内的步长。 如果在find_shape_model函数中没有指定亚像素精度,这个参数指定的精度是可以实现find_shape_mode函数中的角度的。参数AngleStep的选择是基于目标的大小的,如果模板图像太小不能产生许多不同离散角度的图像,因此对于较小的模板图像AngleStep应该设置的比较大。如果AngleExtent不是AngleStep的整数倍, 将会相应的修改AngleStep。

    如果选择 complete pregeneration ,不同角度的模板图像将会产生并保存在内存中。用来存储模板的内存与旋转角度的数目和模板图像的的点数是成正比的。 因此,如果AngleStep太小或是AngleExtent太大, 将会出现该模型不再适合(虚拟)内存的情况。在任何情况下,模型是完全适合主存储器的,因为这避免了操作系统的内存分页,使得寻找匹配模板的时间变短。由于find_shape_model函数中的角度可以使用亚像素精度,一个直径小于200像素的模板可以选择AngleStep>= 1. 如果选择AngleStep=‘auto’ (or 0 向后兼容),create_shape_model将会基于模板的大小自动定义一个合适的角度步长. 自动计算出来的AngleStep可以使用get_shape_model_params函数查看。

    如果没有选择complete pregeneration, 该模型会在每一层金字塔上建立在一个参考的位置。这样在find_shape_model函数运行时,该模型必须转化为不同的角度和尺度在运行时在。正因为如此,匹配该模型可能需要更多的时间。

    对于特别大的模板图像,将参数Optimization设置为不同于’none’的其他数值是非常有用的。如果Optimization= ‘none’, 所有的模型点将要存储。在其他情况下, 按照Optimization的数值会将模型的点数减少. 如果模型点数变少了,必须在find_shape_model函数中将参数Greediness设为一个比较小的值, 比如:0.7、0.8。对于比较小的模型, 减少模型点数并不能提高搜索速度,因为这种情况下通常显着更多的潜在情况的模型必须进行检查。如果Optimization设置为’auto’, create_shape_model自动确定模型的点数。

    Optimization的第二个值定义了模型是否进行预处理(pregenerated completely),是通过选择’pregeneration’或者’no_pregeneration’来设置的。如果不使用第二个值(例如:仅仅设置了第一个值), 默认的是系统中的设置,是通过set_system(‘pregenerate _shape_models’,…)来设置的,对于默认值是 (‘pregenerate_shape_models’ = ‘false’), 模型没有进行预处理. 模型的预处理设置通常会导致比较低的运行时间,因为模型不需要 在运行时间时转换。然而在这种情况下,内存的要求和创建模板所需要的时间是比较高的。还应该指出,不能指望这两个模式返回完全相同的结果,因为在运行时变换一定会导致变换模型和预处理变换模型之间不同的内部数据。比如,如果模型没有 completely pregenerated,在find_shape_model函数中通常返回一个较低的scores,这可能需要将MinScore设置成一个较低的值。此外,在两个模型中插值法获得的位置可能略有不同。如果希望是最高精确度,应该使用最小二乘调整得到模型位置。

    参数Contras决定着模型点的对比度。对比度是用来测量目标与背景之间和目标不同部分之间局部的灰度值差异。Contrast的选择应该确保模板中的主要特征用于模型中。Contrast也可以是两个数值,这时模板使用近似edges_image函数中滞后阈值的算法进行分割。这里第一个数值是比较低的阈值,第二个数值是比较高的阈值。Contrast也可以包含第三个,这个数值是在基于组件尺寸选择重要模型组件时所设置的阈值,比如,比指定的最小尺寸的点数还少的组件将被抑制。这个最小尺寸的阈值会在每相邻的金字塔层之间除以2。如果一个小的模型组件被抑制,但是不使用滞后阈值,然而在Contrast中必须指定三个数值,在这种情况下前两个数值设置成相同的数值。这个参数的设置可以在inspect_shape_model函数中查看效果。如果Contrast设置为’auto’,create_shape_model将会自动确定三个上面描述的数值。或者仅仅自动设置对比度(‘auto_contrast’),滞后阈值(‘auto_contrast_hyst’)或是最小尺寸(‘auto_min_size’)中一个。其他没有自动设置的数值可以按照上面的格式再进行设置。可以允许各种组合,例如:如果设置 [‘auto_contrast’,‘auto_min_size’],对比度和最小尺寸自动确定;如果设置 [‘auto_min_size’,20,30],最小尺寸会自动设定,而滞后阈值被设为20和30。有时候可能对比度阈值自动设置的结果是不满意的,例如,由于一些具体应用的原因当某一个模型组件是被包含或是被抑制时,或是目标包含几种不同的对比度时,手动设置这些参数效果会更好。因此对比度阈值可以使用determine_shape_model_params函数自动确定,也可以在调用create_shape_model之前使用inspect_shape_mode函数检查效果。

    MinContrast用来确定在执行find_shape_model函数进行识别时模型的哪一个对比度必须存在,也就是说,这个参数将模型从噪声图像中分离出来。因此一个好的选择应该是在图像中噪声所引起的灰度变化范围。例如,如果灰度浮动在10个灰度级内,MinContrast应该设置成10。如果模板和搜索图像是多通道图像,Metric参数设置成’ignore_color_polarity’,在一个通道中的噪声必须乘以通道个数的平方根再去设置MinContrast。例如,如果灰度值在一个通道的浮动范围是10个灰度级,图像是三通道的,那么MinContrast应该设置为17。很显然,MinContrast必须小于Contrast。如果要在对比度较低的图像中识别模板,MinContrast必须设置为一个相对较小的数值。如果要是模板即使严重遮挡(occluded)也能识别出来,MinContrast应该设置成一个比噪声引起的灰度浮动范围略大的数值,这样才能确保在find_shape_model函数中提取出模板准确的位置和旋转角度。如果MinContrast设置为’auto’,最小对比度会基于模板图像中的噪声自动定义。因此自动设定仅仅在搜索图像和模板图像噪声近似时才可以使用。此外,在某些情况下为了遮挡的鲁棒性,采用自动设定数值是比较好的。使用get_shape_model_params函数可以查询自动计算的最小对比度。

    参数Metric定义了在图像中匹配模板的条件。如果Metric= ‘use_polarity’,图像中的目标必须和模型具有一样的对比度。例如,如果模型是一个亮的目标在一个暗的背景上,那么仅仅那些比背景亮的目标可以找到。如果Metric= ‘ignore_global_polarity’,在两者对比度完全相反时也能找到目标。在上面的例子中,如果目标是比背景暗的也能将目标找到。find_shape_model函数的运行时间在这种情况下将会略微增加。如果Metric= ‘ignore_local_polarity’, 即使局部对比度改变也能找到模型。例如,当目标包含一部分中等灰度,并且其中部分比较亮部分比较暗时,这种模式是非常有用的。由于这种模式下find_shape_model函数的运行时间显著增加,最好的方法是使用create_shape_model创建几个反映目标可能的对比度变化的模型,同时使用find_shape_models去匹配他们。上面三个metrics仅仅适用于单通道图像。如果是多通道图像作为模板图像或搜索图像,仅仅第一个通道被使用。如果Metric='ignore_color_pol

    arity’, 即使颜色对比度局部变化也能找到模型。例如,当目标的部分区域颜色发生变化(e.g.从红到绿)的情况。如果不能提前知道目标在哪一个通道是可见的这种模式是非常有用的。在这种情况下find_shape_model函数的运行时间也会急剧增加。'ignore_color_polarity’可以使用于具有任意通道数目的图像中。如果使用于单通道图像,他的效果和’ignore_loc al_polarity’是完全相同的。当Metric= ‘ignore_color_polarity’ 时,

    create_shape_model创建的模板通道数目和find_shape_model中的图像通道数目可以是不同的。例如,可以使用综合生成的单通道图像创建模型。另外,这些通道不需要是经过光谱细分(像RGB图像)的。这些通道还可以包括具有在不同方向照亮同一个目标所获得的图像。

    模型图像Template的domain区域的重心是模板的初始位置,可以在set_shape_model_origin函数中设置不同的初始位置。

    LIntExport Herror find_shape_model(

    const Hobject& Image, //搜索图像

    Hlong ModelID, //模板句柄

    Double AngleStart, // 搜索时的起始角度

    Double AngleExtent, //搜索时的角度范围,必须与创建模板时的有交集。

    Double MinScore, // 输出的匹配的质量系数Score 都得大于该值

    Hlong NumMatches, // 定义要输出的匹配的最大个数

    Double MaxOverlap, // 当找到的目标存在重叠时,且重叠大于该值时选择一个好的输出

    const char* SubPixel, // 计算精度的设置,五种模式,多选2,3

    Hlong NumLevels, // 搜索时金字塔的层数

    Double Greediness , //贪婪度,搜索启发式,一般都设为0.9,越高速度快容易出现找不到的情况

    Halcon::HTuple* Row, //输出匹配位置的行坐标

    Halcon::HTuple* Column, //输出匹配位置的列坐标

    Halcon::HTuple* Angle, //输出匹配角度

    Halcon::HTuple* Score ) //输出匹配质量

    进一步分析:

    注意Row、Column的坐标并不是模板在搜索图像中的精确位置,因此不能直接使用他们。这些数值是为了创建变换矩阵被优化后的,你可以用这个矩阵的匹配结果完成各种任务,比如调整后续步骤的ROI。

    Score是一个0到1之间的数,是模板在搜索图像中可见比例的近似测量。如果模板的一半被遮挡,该值就不能超过0.5。

    Image的domain定义了模型参考点的搜索区域,模型参考点是在create_shape_model中用来创建模型的图像的domain区域的重心。不考虑使用函数set_shape_model_origin设置不同的初始位置。在图像domain区域的这些点内搜索模型,其中模型完全属于这幅图像。这意味着如果模型超出图像边界,即使获得的质量系数(score)大于MinScore也不能找到模型。这种性能可以通过set_system(‘border_shape_models’,‘true’)改变,这样那些超出图像边界,质量系数大于MinScore的模型也能找到。这时那些在图像外面的点看作是被遮挡了,可以降低质量系数。在这种模式下搜索的时间将要增加。

    参数AngleStart和AngleExtent确定了模型搜索的旋转角度,如果有必要,旋转的范围会被截取成为create_shape_model函数中给定的旋转范围。这意味着创建模型和搜索时的角度范围必须真正的重叠。在搜索时的角度范围不会改变为模2*pi的。为了简化介绍,在该段落剩下的部分所有角度都用度来表示,而在find_shape_model函数中使用弧度来设置的。因此,如果创建模板时,AngleStart=-20°、AngleExtent=40°,在搜索模板函数find_shape_model中设置AngleStart=350°、AngleExtent=20°,尽管角度模360后是重叠的,还是会找不到模板的。为了找到模板,在这个例子中必须将AngleStart=350°改为AngleStart=-10°。

    参数MinScore定义模板匹配时至少有个什么样的质量系数才算是在图像中找到模板。MinScore设置的越大,搜索的就越快。如果模板在图像中没有被遮挡,MinScore可以设置为0.8这么高甚至0.9。

    NumMatches定义了在图像上找到模板的最大的个数。如果匹配时的质量系数大于MinScore的目标个数多于NumMatches,仅仅返回质量系数最好的NumMatches个目标位置。如果找的匹配目标不足NumMatches,那么就只返回找到的这几个。参数MinScore优于NumMatches。

    如果模型具有对称性,会在搜索图像的同一位置和不同角度上找到多个与目标匹配的区域。参数MaxOverlap是0到1之间的,定义了找到的两个目标区域最多重叠的系数,以便于把他们作为两个不同的目标区域分别返回。如果找到的两个目标区域彼此重叠并且大于MaxOverlap,仅仅返回效果最好的一个。重叠的计算方法是基于找到的目标区域的任意方向的最小外接矩形(看smallest_rectangle2)。如果MaxOverlap=0, 找到的目标区域不能存在重叠, 如果MaxOverla p=1,所有找到的目标区域都要返回。

    SubPixel确定找到的目标是否使用亚像素精度提取。如果SubPixel设置为’none’(或者’false’ 背景兼容),模型的位置仅仅是一个像素精度和在create_shape_model中定义的角度分辨率。如果SubPixel设置为’interpo lation’(或’true’),位置和角度都是亚像素精度的。在这种模式下模型的位置是在质量系数函数中插入的,这种模式几乎不花费计算时间,并且能达到足够高的精度,被广泛使用。然而在一些精度要求极高的应用中,模板的位置应该通过最小二乘调整决定,比如通过最小化模板点到相关图像点的距离。与 ‘interpolation’相比,这种模式需要额外的计算时间。对于最小二乘调整的模式有:‘least_squares’, ‘least_squares_high’, 和’least_squares_very_high’。他们可用来定义被搜索的最小距离的精度,选择的精度越高,亚像素提取的时间越长。然而,通常SubPixel设置为’interpolation’。如果希望设置最小二乘就选择’least_squares’, 因为这样才能确保运行时间和精度的权衡。

    NumLevels是在搜索时使用的金字塔层数,如有必要,层数截成创建模型时的范围。如果NumLevels=0,使用创建模板时金字塔的层数。另外NumLevels还可以包含第二个参数,这个参数定义了找到匹配模板的最低金字塔层数。NumLevels=[4,2]表示匹配在第四层金字塔开始,在第二层金字塔找到匹配(最低的设为1)。可以使用这种方法降低匹配的运行时间,但是这种模式下位置精度是比正常模式下低的,所谓正常模式是在金字塔最底层匹配。因此如果需要较高的精度,应该设置SubPixel至少为’least_squares’。如果金字塔最底层设置的过大,可能不会达到期望的精度,或者找到一个不正确的匹配区域。这是因为在较高层的金字塔上模板是不够具体的,不足以找到可靠的模板最佳匹配。在这种情况下最低金字塔层数应设为最小值。

    参数Greediness确定在搜索时的“贪婪程度”。如果Greediness=0,使用一个安全的搜索启发式,只要模板在图像中存在就一定能找到模板,然而这种方式下搜索是相对浪费时间的。如果Greediness=1,使用不安全的搜索启发式,这样即使模板存在于图像中,也有可能找不到模板,但只是少数情况。如果设置Greediness=0.9,在几乎所有的情况下,总能找到模型的匹配。

    转载自halcon学习网Trevan帖子。

    展开全文
  • matlab人脸匹配代码ASM人脸特征点匹配 简介任务:通过结合使用主动形状模型(ASM)算法和PCA方法来训练可用于检测测试人脸图片关键点的回归器。 然后,进行Procrustes分析以匹配那些图像的重建关键点; 环境:Matlab...
  • 光谱匹配 matlab代码
  • 2019论文中描述的用于3D表面形状匹配的Local Point Signature的MATLAB实现: 王一群,郭建伟,严东明,王凯,张晓鹏。 如果此代码/程序(或其部分)对您的项目有利,请考虑引用以上文章。 用法 请使用MATLAB运行“ ...
  • 图像矩阵matlab代码该软件包包括用于与层次投影不变上下文进行形状匹配的所有必要源代码。 其中一些功能是通过直接使用或略微更改Ling和Jacobs在线提供的内部距离形状上下文(IDSC)的源代码以及Xiang Bai在线提供的...
  • 用过halcon形状匹配的都知道,这个算子贼好用,随便截一个ROI做模板就可以在搜索图像中匹配到相似的区域,并且能输出搜索图像的位置,匹配尺度,匹配角度。现在我们就要利用opencv在C++的环境下复现这个效果。 我们...

    前言

    用过halcon形状匹配的都知道,这个算子贼好用,随便截一个ROI做模板就可以在搜索图像中匹配到相似的区域,并且能输出搜索图像的位置,匹配尺度,匹配角度。现在我们就要利用opencv在C++的环境下复现这个效果。

    我们先看下复现的效果图,提升下学习的欲望(要在搜索图像中找到所有的K字母)。

    下图是模板图像,为一个"K"字母。

    下图是待搜索的图像,其中的K字符存在旋转,缩放,残缺遮挡,要利用上面的"K"字母模板在下图中找到所以的K字母。并输出它的位置,旋转角度,尺度,相似度。

    下图是经过形状匹配后的结果图像,可以看到匹配的结果,用蓝色画出来,输出相似度,匹配角度,尺度。匹配时间为0.201040s。

    原理简介

    不具备多角度,多尺度的形状匹配原理

    1.对模板图像进行特征提取,并存储特征信息

    2.对搜索图像进行特征提取

    3.将模板图像的特征信息在搜索图像上进行特征相似度比对,然后滑动窗口继续比对

    4.直到比对完所有的搜索图像区域,则生成相似度矩阵

    5.根据具体需求对相似度矩阵进行操作

    具备多角度,多尺度的形状匹配原理

    1.对模板图像进行特征提取,并存储特征信息,然后对模板图像进行多角度变换(旋转)后进行特征提前,然后存储,再对各个多角度变换后的模板图像进行多尺度变换(缩放),再存储。则一个旋转范围为0~360°,旋转步长为1°,缩放范围为0.9~1.1,缩放步长为0.1的模板个数为((360-0)/1) * ((1.1-0.9)/0.1) = 360 * 20 = 7200个模板,制作模板就是存储这7200个模板特征

    2.对搜索图像进行特征提取

    3.将所有的模板图像(比如这边的7200个)的特征信息依次在搜索图像上进行特征相似度比对,然后滑动窗口继续比对

    4.直到所有的模板比对完所有的搜索图像区域,则生成相似度矩阵(7200个模板就是7200个相似度矩阵)

    5.根据具体需求对所有的相似度矩阵进行操作

    流程简介

    主要流程就是两步

    1.制作模板

    2.开始匹配(使用模板在搜索图像上进行匹配)

    制作模板

    对所有的模板进行特征提取,这边使用类似于canny边缘提取的算法来提取特征。

    开始匹配

    对搜索图像提取特征,将所有的模板在搜索图像上进行相似度计算,然后滑动窗口直到匹配完成。

    代码实现与分析

    我们将实现三个主要的函数功能:制作模板(MakingTemplates),加载模板(LoadModel),匹配(Matching)。

    /*
    @model: 输入图像
    @angle_range: 角度范围
    @scale_range: 尺度范围
    @num_features: 特征数
    @weak_thresh:弱阈值
    @strong_thresh: 强阈值
    @mask: 掩码
    */
    void MakingTemplates(Mat model, AngleRange angle_range, ScaleRange scale_range,
    	int num_features, float weak_thresh = 30.0f, float strong_thresh = 60.0f,
    	Mat mask = Mat());
    /*
    加载模型
    */
    void LoadModel();
    /*
    @source: 输入图像
    @score_thresh: 匹配分数阈值
    @overlap: 重叠阈值
    @mag_thresh: 最小梯度阈值
    @greediness: 贪婪度,越小匹配越快,但是可能无法匹配到目标
    @pyrd_level: 金字塔层数,越大匹配越快,但是可能无法匹配到目标
    @T: T参数
    @top_k: 最多匹配多少个
    @strategy: 精确匹配(0), 普通匹配(1), 粗略匹配(2)
    @mask: 匹配掩码
    */
    vector<Match> Matching(Mat source, float score_thresh = 0.9f, float overlap = 0.4f,
    	float mag_thresh = 30.f, float greediness = 0.8f, PyramidLevel pyrd_level = PyramidLevel_3,
    	int T = 2, int top_k = 0, MatchingStrategy strategy = Strategy_Accurate, const Mat mask = Mat());

    制作模板

    先来看两个结构体AngleRange(角度范围(起始角度,终止角度,角度步长))和ScaleRange(尺度范围(起始尺度,终止尺度,角度尺度))。

    struct MatchRange
    {
    	float begin;
    	float end;
    	float step;
    
    	MatchRange() : begin(0.f), end(0.f), step(0.f) {}
    	MatchRange(float b, float e, float s);
    };
    inline MatchRange::MatchRange(float b, float e, float s) : begin(b), end(e), step(s) {}
    typedef struct MatchRange AngleRange; // 角度范围(起始角度,终止角度,角度步长)
    typedef struct MatchRange ScaleRange; // 尺度范围(起始尺度,终止尺度,角度尺度)

    制作模板的代码

    void KcgMatch::MakingTemplates(Mat model, AngleRange angle_range, ScaleRange scale_range,
    	int num_features,float weak_thresh, float strong_thresh, Mat mask) {
    
    	ClearModel();
    	// padding模板和模板掩码
    	// 为什么要padding呢?因为在制作旋转模板的时候可能丢失有效区域
    	PaddingModelAndMask(model, mask, scale_range.end);
    	// 初始化角度,尺度范围
    	angle_range_ = angle_range;
    	scale_range_ = scale_range;
    	// 生成所有的模板信息
    	vector<ShapeInfo> shape_infos = ProduceShapeInfos(angle_range, scale_range);
    	vector<Mat> l0_mdls; l0_mdls.clear();
    	vector<Mat> l0_msks; l0_msks.clear();
    	// 生成所有模板的底层金字塔图像
    	for (int s = 0; s < shape_infos.size(); s++) {
    
    		l0_mdls.push_back(MdlOf(model, shape_infos[s]));
    		l0_msks.push_back(MskOf(mask, shape_infos[s]));
    	}
    	// 对所有层的金字塔图像进行特征提取(这边最多8层,只制作8层,足够用)
    	for (int p = 0; p <= PyramidLevel_7; p++) {
    
    		// 某层金字塔的所有角度,尺度图像
    		for (int s = 0; s < shape_infos.size(); s++) {
    
    			Mat mdl_pyrd = l0_mdls[s];
    			Mat msk_pyrd = l0_msks[s];
    			if (p > 0) {
    
    				Size sz = Size(l0_mdls[s].cols >> 1, l0_mdls[s].rows >> 1);
    				pyrDown(l0_mdls[s], mdl_pyrd, sz);
    				pyrDown(l0_msks[s], msk_pyrd, sz);
    			}
    			// 为什么要erode?因为有效特征信息可能在边缘
    			erode(msk_pyrd, msk_pyrd, Mat(), Point(-1, -1), 1, BORDER_REPLICATE);
    			l0_mdls[s] = mdl_pyrd;
    			l0_msks[s] = msk_pyrd;
    
    			// 计算某层金字塔需要的特征数
    			int features_pyrd = (int)((num_features >> p) * shape_infos[s].scale);
    
    			Mat mag8, angle8, quantized_angle8;
    			// 量化边缘特征(8个方向)
    			QuantifyEdge(mdl_pyrd, angle8, quantized_angle8, mag8, weak_thresh, false);
    			// 提取模板信息(8个方向)
    			Template templ = ExtractTemplate(	angle8, quantized_angle8, mag8, 
    												shape_infos[s], PyramidLevel(p),
    												weak_thresh, strong_thresh,
    												features_pyrd, msk_pyrd);
    			templ_all_[p].push_back(templ);
    
    			Mat mag180, angle180, quantized_angle180;
    			// 量化边缘特征(180个方向)
    			QuantifyEdge(mdl_pyrd, angle180, quantized_angle180, mag180, weak_thresh, true);
    			// 提取模板信息(180个方向)
    			templ = ExtractTemplate(	angle180, quantized_angle180, mag180,
    										shape_infos[s], PyramidLevel(p),
    										weak_thresh, strong_thresh,
    										features_pyrd, msk_pyrd);
    			templ_all_[p + 8].push_back(templ);
    
    			// 画出提取过程
    			Mat draw_mask;
    			msk_pyrd.copyTo(draw_mask);
    			DrawTemplate(draw_mask, templ, Scalar(0));
    			imshow("draw_mask", draw_mask);
    			waitKey(1);
    		}
    		cout << "train pyramid level " << p << " complete." << endl;
    	}
    	// 保存模板
    	SaveModel();
    }

    加载模板

    void KcgMatch::LoadModel() {
    
    	// 制作好的模板为一个yaml文件,里面存取了特征信息,这边就是读取该文件
    	// 将所有的模板特征读到内存
    	ClearModel();
    	string model_name = model_root_ + class_name_ + KCG_MODEL_SUFFUX;
    	FileStorage fs(model_name, FileStorage::READ);
    	assert(fs.isOpened() && "load model failed.");
    	FileNode fn = fs.root();
    	angle_range_.begin = fn["angle_range_bgin"];
    	angle_range_.end = fn["angle_range_end"];
    	angle_range_.step = fn["angle_range_step"];
    	scale_range_.begin = fn["scale_range_bgin"];
    	scale_range_.end = fn["scale_range_end"];
    	scale_range_.step = fn["scale_range_step"];
    	FileNode tps_fn = fn["templates"];
    	FileNodeIterator tps_it = tps_fn.begin(), tps_it_end = tps_fn.end();
    	for (; tps_it != tps_it_end; ++tps_it)
    	{
    		int template_id = (*tps_it)["template_id"];
    		FileNode pyrds_fn = (*tps_it)["template_pyrds"];
    		FileNodeIterator pyrd_it = pyrds_fn.begin(), pyrd_it_end = pyrds_fn.end();
    		int pl = 0;
    		for (; pyrd_it != pyrd_it_end; ++pyrd_it)
    		{
    			FileNode pyrd_fn = (*pyrd_it);
    			Template templ;
    			templ.id = pyrd_fn["id"];
    			templ.pyramid_level = pyrd_fn["pyramid_level"];
    			templ.is_valid = pyrd_fn["is_valid"];
    			templ.x = pyrd_fn["x"];
    			templ.y = pyrd_fn["y"];
    			templ.w = pyrd_fn["w"];
    			templ.h = pyrd_fn["h"];
    			templ.shape_info.scale = pyrd_fn["shape_scale"];
    			templ.shape_info.angle = pyrd_fn["shape_angle"];
    			FileNode features_fn = pyrd_fn["features"];
    			FileNodeIterator feature_it = features_fn.begin(), feature_it_end = features_fn.end();
    			for (; feature_it != feature_it_end; ++feature_it)
    			{
    				FileNode feature_fn = (*feature_it);
    				FileNodeIterator feature_info = feature_fn.begin();
    				Feature feat;
    				feature_info >> feat.x >> feat.y >> feat.lbl;
    				templ.features.push_back(feat);
    			}
    			templ_all_[pl].push_back(templ);
    			pl++;
    		}
    	}
    
    	LoadRegion8Idxes();
    }

    匹配

    vector<Match> KcgMatch::Matching(Mat source, float score_thresh, float overlap,
    	float mag_thresh, float greediness, PyramidLevel pyrd_level, int T, int top_k,
    	MatchingStrategy strategy, const Mat mask) {
    
    	// 初始化匹配参数
    	InitMatchParameter(score_thresh, overlap, mag_thresh, greediness, T, top_k, strategy);
    	// 获取搜索图像进行所有的有效金字塔层图像
    	GetAllPyramidLevelValidSource(source, pyrd_level);
    
    	vector<Match> matches; 
    	// 从最高层金字塔开始匹配量化为8方向的图像相似度矩阵
    	matches = MatchingPyrd8(sources_[pyrd_level], pyrd_level, region8_idxes_);
    	// 获取前K个最相似match
    	matches = GetTopKMatches(matches);
    
    	// 再次确认
    	matches = ReconfirmMatches(matches, pyrd_level);
    	matches = GetTopKMatches(matches);
    
    	// 最后匹配,从金字塔顶层还原到底层ROI开始匹配
    	matches = MatchingFinal(matches, pyrd_level);
    	matches = GetTopKMatches(matches);
    
    	// 返回指定的匹配
    	return matches;
    }

    完整代码

    KcgMatch.h

    /*M///
    //
    // Author	: KayChan
    // Explain	: Shape matching
    //
    //M*/
    
    #ifndef _KCG_MATCH_H_
    #define _KCG_MATCH_H_
    
    #include <opencv2/opencv.hpp>
    #include <omp.h>
     
    #ifndef ATTR_ALIGN
    #  if defined(__GNUC__)
    #    define ATTR_ALIGN(n)	__attribute__((aligned(n)))
    #  else
    #    define ATTR_ALIGN(n)	__declspec(align(n))
    #  endif
    #endif // #ifndef ATTR_ALIGN
    
    using namespace cv;
    using namespace std;
    
    namespace kcg{
    
    struct MatchRange
    {
    	float begin;
    	float end;
    	float step;
    
    	MatchRange() : begin(0.f), end(0.f), step(0.f) {}
    	MatchRange(float b, float e, float s);
    };
    inline MatchRange::MatchRange(float b, float e, float s) : begin(b), end(e), step(s) {}
    typedef struct MatchRange AngleRange;
    typedef struct MatchRange ScaleRange;
    
    typedef struct ShapeInfo_S
    {
    	float angle;
    	float scale;
    }ShapeInfo;
    
    typedef struct Feature_S
    {
    	int x;
    	int y;
    	int lbl;
    }Feature;
    
    typedef struct Candidate_S
    {
    	/// Sort candidates with high score to the front
    	bool operator<(const struct Candidate_S &rhs) const
    	{
    		return score > rhs.score;
    	}
    	float score;
    	Feature feature;
    }Candidate;
    
    typedef struct Template_S
    {
    	int id = 0;
    	int pyramid_level = 0;
    	int is_valid = 0;
    	int x = 0;
    	int y = 0;
    	int w = 0;
    	int h = 0;
    	ShapeInfo shape_info;
    	vector<Feature> features;
    }Template;
    
    typedef struct Match_S
    {
    	/// Sort matches with high similarity to the front
    	bool operator<(const struct Match_S &rhs) const
    	{
    		// Secondarily sort on template_id for the sake of duplicate removal
    		if (similarity != rhs.similarity)
    			return similarity > rhs.similarity;
    		else
    			return template_id < rhs.template_id;
    	}
    
    	bool operator==(const struct Match_S &rhs) const
    	{
    		return x == rhs.x && y == rhs.y && similarity == rhs.similarity;
    	}
    
    	int x;
    	int y;
    	float similarity;
    	int template_id;
    }Match;
    
    typedef enum PyramidLevel_E
    {
    	PyramidLevel_0 = 0,
    	PyramidLevel_1 = 1,
    	PyramidLevel_2 = 2,
    	PyramidLevel_3 = 3,
    	PyramidLevel_4 = 4,
    	PyramidLevel_5 = 5,
    	PyramidLevel_6 = 6,
    	PyramidLevel_7 = 7,
    	PyramidLevel_TabooUse = 16,
    }PyramidLevel;
    
    typedef enum MatchingStrategy_E
    {
    	Strategy_Accurate = 0,
    	Strategy_Middling = 1,
    	Strategy_Rough = 2,
    }MatchingStrategy;
    
    class KcgMatch
    {
    public:
    
    	KcgMatch(string model_root, string class_name);
    	~KcgMatch();
    	/*
    	@model: 输入图像
    	@angle_range: 角度范围
    	@scale_range: 尺度范围
    	@num_features: 特征数
    	@weak_thresh:弱阈值
    	@strong_thresh: 强阈值
    	@mask: 掩码
    	*/
    	void MakingTemplates(Mat model, AngleRange angle_range, ScaleRange scale_range,
    		int num_features, float weak_thresh = 30.0f, float strong_thresh = 60.0f,
    		Mat mask = Mat());
    	/*
    	加载模型
    	*/
    	void LoadModel();
    	/*
    	@source: 输入图像
    	@score_thresh: 匹配分数阈值
    	@overlap: 重叠阈值
    	@mag_thresh: 最小梯度阈值
    	@greediness: 贪婪度,越小匹配越快,但是可能无法匹配到目标
    	@pyrd_level: 金字塔层数,越大匹配越快,但是可能无法匹配到目标
    	@T: T参数
    	@top_k: 最多匹配多少个
    	@strategy: 精确匹配(0), 普通匹配(1), 粗略匹配(2)
    	@mask: 匹配掩码
    	*/
    	vector<Match> Matching(Mat source, float score_thresh = 0.9f, float overlap = 0.4f,
    		float mag_thresh = 30.f, float greediness = 0.8f, PyramidLevel pyrd_level = PyramidLevel_3,
    		int T = 2, int top_k = 0, MatchingStrategy strategy = Strategy_Accurate, const Mat mask = Mat());
    	void DrawMatches(Mat &image, vector<Match> matches, Scalar color);
    
    protected:
    	void PaddingModelAndMask(Mat &model, Mat &mask, float max_scale);
    	vector<ShapeInfo> ProduceShapeInfos(AngleRange angle_range, ScaleRange scale_range);
    	Mat Transform(Mat src, float angle, float scale);
    	Mat MdlOf(Mat model, ShapeInfo info);
    	Mat MskOf(Mat mask, ShapeInfo info);
    	void DrawTemplate(Mat &image, Template templ, Scalar color);
    	void QuantifyEdge(Mat image, Mat &angle, Mat &quantized_angle, Mat &mag, float mag_thresh, bool calc_180 = true);
    	void Quantify8(Mat angle, Mat &quantized_angle, Mat mag, float mag_thresh);
    	void Quantify180(Mat angle, Mat &quantized_angle, Mat mag, float mag_thresh);
    	Template ExtractTemplate(Mat angle, Mat quantized_angle, Mat mag, ShapeInfo shape_info, 
    		PyramidLevel pl, float weak_thresh, float strong_thresh, int num_features, Mat mask);
    	Template SelectScatteredFeatures(vector<Candidate> candidates, int num_features, float distance);
    	Rect CropTemplate(Template &templ);
    	void LoadRegion8Idxes();
    	void ClearModel();
    	void SaveModel();
    	void InitMatchParameter(float score_thresh, float overlap, float mag_thresh, float greediness, int T, int top_k, MatchingStrategy strategy);
    	void GetAllPyramidLevelValidSource(Mat &source, PyramidLevel pyrd_level);
    	vector<Match> GetTopKMatches(vector<Match> matches);
    	vector<Match> DoNmsMatches(vector<Match> matches, PyramidLevel pl, float overlap);
    	vector<Match> MatchingPyrd180(Mat src, PyramidLevel pl, vector<int> region_idxes = vector<int>());
    	vector<Match> MatchingPyrd8(Mat src, PyramidLevel pl, vector<int> region_idxes = vector<int>());
    	void Spread(const Mat quantized_angle, Mat &spread_angle, int T);
    	void ComputeResponseMaps(const Mat spread_angle, vector<Mat> &response_maps);
    	bool CalcPyUpRoiAndStartPoint(PyramidLevel cur_pl, PyramidLevel obj_pl, Match match,
    		Mat &r, Point &p, bool is_padding = false);
    	void CalcRegionIndexes(vector<int> &region_idxes, Match match, MatchingStrategy strategy);
    	vector<Match> ReconfirmMatches(vector<Match> matches, PyramidLevel pl);
    	vector<Match> MatchingFinal(vector<Match> matches, PyramidLevel pl);
    
    private:
    	typedef vector<Template> TemplateMatchRange;
    	TemplateMatchRange templ_all_[PyramidLevel_TabooUse];
    	vector<Mat> sources_;
    	ATTR_ALIGN(32) float score_table_[180][180];
    	ATTR_ALIGN(8) unsigned char score_table_8map_[8][256];
    	string model_root_;
    	string class_name_;
    	AngleRange angle_range_;
    	ScaleRange scale_range_;
    	vector<int> region8_idxes_;
    
    	float score_thresh_;
    	float overlap_;
    	float mag_thresh_;
    	float greediness_;
    	int T_;
    	int top_k_;
    	MatchingStrategy strategy_;
    };
    }
    
    #endif
    
    

    KcgMatch.cpp

    #include "KcgMatch.h"
    #include <math.h>
    
    #define KCG_EPS 0.00001f
    #define KCG_PI	3.1415926535897932384626433832795f
    #define KCG_MODEL_SUFFUX string(".yaml")
    
    const float AngleRegionTable[16][2] = {
    
    	0.f		, 22.5f	,
    	22.5f	, 45.f	,
    	45.f	, 67.5f	,
    	67.5f	, 90.f	,
    	90.f	, 112.5f,
    	112.5f	, 135.f	,
    	135.f	, 157.5f,
    	157.5f	, 180.f,
    	180.f	, 202.5f,
    	202.5f	, 225.f,
    	225.f	, 247.5f,
    	247.5f	, 270.f,
    	270.f	, 292.5f,
    	292.5f	, 315.f,
    	315.f	, 337.5f,
    	337.5f	, 360.f
    };
    
    namespace cv_dnn_nms {
    
    template <typename T>
    static inline bool SortScorePairDescend(const std::pair<float, T>& pair1, const std::pair<float, T>& pair2) {
    	
    	return pair1.first > pair2.first;
    }
    
    inline void GetMaxScoreIndex(const std::vector<float>& scores, const float threshold, const int top_k,
    	std::vector<std::pair<float, int> >& score_index_vec) {
    
    	for (size_t i = 0; i < scores.size(); ++i)
    	{
    		if (scores[i] > threshold)
    		{
    			//score_index_vec.push_back(std::make_pair(scores[i], i));
    			std::pair<float, int> psi;
    			psi.first = scores[i];
    			psi.second = (int)i;
    			score_index_vec.push_back(psi);
    		}
    	}
    	std::stable_sort(score_index_vec.begin(), score_index_vec.end(),
    		SortScorePairDescend<int>);
    	if (top_k > 0 && top_k < (int)score_index_vec.size())
    	{
    		score_index_vec.resize(top_k);
    	}
    }
    
    template <typename BoxType>
    inline void NMSFast_(const std::vector<BoxType>& bboxes,
    	const std::vector<float>& scores, const float score_threshold,
    	const float nms_threshold, const float eta, const int top_k,
    	std::vector<int>& indices, float(*computeOverlap)(const BoxType&, const BoxType&)) {
    
    	CV_Assert(bboxes.size() == scores.size());
    	std::vector<std::pair<float, int> > score_index_vec;
    	GetMaxScoreIndex(scores, score_threshold, top_k, score_index_vec);
    
    	float adaptive_threshold = nms_threshold;
    	indices.clear();
    	for (size_t i = 0; i < score_index_vec.size(); ++i) {
    		const int idx = score_index_vec[i].second;
    		bool keep = true;
    		for (int k = 0; k < (int)indices.size() && keep; ++k) {
    			const int kept_idx = indices[k];
    			float overlap = computeOverlap(bboxes[idx], bboxes[kept_idx]);
    			keep = overlap <= adaptive_threshold;
    		}
    		if (keep)
    			indices.push_back(idx);
    		if (keep && eta < 1 && adaptive_threshold > 0.5) {
    			adaptive_threshold *= eta;
    		}
    	}
    }
    
    template<typename _Tp> static inline
    	double jaccardDistance__(const Rect_<_Tp>& a, const Rect_<_Tp>& b) {
    	_Tp Aa = a.area();
    	_Tp Ab = b.area();
    
    	if ((Aa + Ab) <= std::numeric_limits<_Tp>::epsilon()) {
    		// jaccard_index = 1 -> distance = 0
    		return 0.0;
    	}
    
    	double Aab = (a & b).area();
    	// distance = 1 - jaccard_index
    	return 1.0 - Aab / (Aa + Ab - Aab);
    }
    
    template <typename T>
    static inline float rectOverlap(const T& a, const T& b) {
    
    	return 1.f - static_cast<float>(jaccardDistance__(a, b));
    }
    
    void NMSBoxes(const std::vector<Rect>& bboxes, const std::vector<float>& scores,
    	const float score_threshold, const float nms_threshold,
    	std::vector<int>& indices, const float eta = 1, const int top_k = 0) {
    
    	NMSFast_(bboxes, scores, score_threshold, nms_threshold, eta, top_k, indices, rectOverlap);
    }
    
    } // end namespace cv_dnn_nms
    
    namespace kcg_matching{
    
    KcgMatch::KcgMatch(string model_root, string class_name) {
    
    	assert(!model_root.empty() && "model_root should not empty.");
    	assert(!class_name.empty() && "class_name should not empty.");
    	if (model_root[model_root.length() - 1] != '/') {
    
    		model_root.push_back('/');
    	}
    	model_root_ = model_root;
    	class_name_ = class_name;
    
    	/// Create 180*180 table
    	for (int i = 0; i < 180; i++) {
    
    		for (int j = 0; j < 180; j++) {
    
    			float rad = (i - j) * KCG_PI / 180.f;
    			score_table_[i][j] = fabs(cosf(rad));
    		}
    	}
    
    	/// Create 8*8 table
    	ATTR_ALIGN(8) unsigned char score_table_8d[8][8];
    	for (int i = 0; i < 8; i++) {
    
    		for (int j = 0; j < 8; j++) {
    
    			float rad = (i - j) * (180.f / 8.f) * KCG_PI / 180.f;
    			score_table_8d[i][j] = (unsigned char)(fabs(cosf(rad))*100.f);
    		}
    	}
    
    	/// Create 8*256 table
    	for (int i = 0; i < 8; i++) {
    
    		for (int j = 0; j < 256; j++) {
    
    			unsigned char max_score = 0;
    			for (int shift_time = 0; shift_time < 8; shift_time++) {
    
    				unsigned char flg = (j >> shift_time) & 0b00000001;
    				if (flg) {
    
    					if (score_table_8d[i][shift_time] > max_score) {
    
    						max_score = score_table_8d[i][shift_time];
    					}
    				}
    			}
    			score_table_8map_[i][j] = max_score;
    		}
    	}
    }
    
    KcgMatch::~KcgMatch() {
    
    }
    
    void KcgMatch::MakingTemplates(Mat model, AngleRange angle_range, ScaleRange scale_range,
    	int num_features,float weak_thresh, float strong_thresh, Mat mask) {
    
    	ClearModel();
    	PaddingModelAndMask(model, mask, scale_range.end);
    	angle_range_ = angle_range;
    	scale_range_ = scale_range;
    	vector<ShapeInfo> shape_infos = ProduceShapeInfos(angle_range, scale_range);
    	vector<Mat> l0_mdls; l0_mdls.clear();
    	vector<Mat> l0_msks; l0_msks.clear();
    	for (int s = 0; s < shape_infos.size(); s++) {
    
    		l0_mdls.push_back(MdlOf(model, shape_infos[s]));
    		l0_msks.push_back(MskOf(mask, shape_infos[s]));
    	}
    	for (int p = 0; p <= PyramidLevel_7; p++) {
    
    		for (int s = 0; s < shape_infos.size(); s++) {
    
    			Mat mdl_pyrd = l0_mdls[s];
    			Mat msk_pyrd = l0_msks[s];
    			if (p > 0) {
    
    				Size sz = Size(l0_mdls[s].cols >> 1, l0_mdls[s].rows >> 1);
    				pyrDown(l0_mdls[s], mdl_pyrd, sz);
    				pyrDown(l0_msks[s], msk_pyrd, sz);
    			}
    			erode(msk_pyrd, msk_pyrd, Mat(), Point(-1, -1), 1, BORDER_REPLICATE);
    			l0_mdls[s] = mdl_pyrd;
    			l0_msks[s] = msk_pyrd;
    
    			int features_pyrd = (int)((num_features >> p) * shape_infos[s].scale);
    
    			Mat mag8, angle8, quantized_angle8;
    			QuantifyEdge(mdl_pyrd, angle8, quantized_angle8, mag8, weak_thresh, false);
    			Template templ = ExtractTemplate(	angle8, quantized_angle8, mag8, 
    												shape_infos[s], PyramidLevel(p),
    												weak_thresh, strong_thresh,
    												features_pyrd, msk_pyrd);
    			templ_all_[p].push_back(templ);
    
    			Mat mag180, angle180, quantized_angle180;
    			QuantifyEdge(mdl_pyrd, angle180, quantized_angle180, mag180, weak_thresh, true);
    			templ = ExtractTemplate(	angle180, quantized_angle180, mag180,
    										shape_infos[s], PyramidLevel(p),
    										weak_thresh, strong_thresh,
    										features_pyrd, msk_pyrd);
    			templ_all_[p + 8].push_back(templ);
    
    			/// draw
    			/*Mat draw_mask;
    			msk_pyrd.copyTo(draw_mask);
    			DrawTemplate(draw_mask, templ, Scalar(0));
    			imshow("draw_mask", draw_mask);
    			waitKey(1);*/
    		}
    		cout << "train pyramid level " << p << " complete." << endl;
    	}
    	SaveModel();
    }
    
    vector<Match> KcgMatch::Matching(Mat source, float score_thresh, float overlap,
    	float mag_thresh, float greediness, PyramidLevel pyrd_level, int T, int top_k,
    	MatchingStrategy strategy, const Mat mask) {
    
    	InitMatchParameter(score_thresh, overlap, mag_thresh, greediness, T, top_k, strategy);
    	GetAllPyramidLevelValidSource(source, pyrd_level);
    
    	vector<Match> matches; 
    	matches = MatchingPyrd8(sources_[pyrd_level], pyrd_level, region8_idxes_);
    	matches = GetTopKMatches(matches);
    
    	matches = ReconfirmMatches(matches, pyrd_level);
    	matches = GetTopKMatches(matches);
    
    	matches = MatchingFinal(matches, pyrd_level);
    	matches = GetTopKMatches(matches);
    
    	return matches;
    }
    
    void KcgMatch::DrawMatches(Mat &image, vector<Match> matches, Scalar color) {
    
    	//#pragma omp parallel for
    	for (int i = 0; i < matches.size(); i++) {
    
    		auto match = matches[i];
    		auto templ = templ_all_[8][match.template_id];
    		int w = match.x + templ.w;
    		int h = match.y + templ.h;
    		for (int i = 0; i < (int)templ.features.size(); i++) {
    
    			auto feature = templ.features[i];
    			//circle(image, cv::Point(match.x + feature.x, match.y + feature.y), 1, color, 1);
    			line(image,
    				Point(match.x + feature.x, match.y + feature.y),
    				Point(match.x + feature.x, match.y + feature.y),
    				color, 1);
    		}
    		cv::rectangle(image, { match.x, match.y }, { w, h }, color, 1);
    		char info[128];
    		sprintf(info,
    				"%.2f%% [%.2f, %.2f]",
    				match.similarity * 100,
    				templ.shape_info.angle,
    				templ.shape_info.scale);
    		cv::putText(image,
    					info,
    					Point(match.x, match.y), FONT_HERSHEY_PLAIN, 1.f, color, 1);
    	}
    }
    
    void KcgMatch::PaddingModelAndMask(Mat &model, Mat &mask, float max_scale) {
    
    	CV_Assert(!model.empty() && "model is empty.");
    	if (mask.empty()) 
    		mask = Mat(model.size(), CV_8UC1, { 255 });
    	else 
    		CV_Assert(model.size() == mask.size());
    	int min_side_length = std::min(model.rows, model.cols);
    	int diagonal_line_length =
    		(int)ceil(std::sqrt(model.rows*model.rows + model.cols*model.cols)*max_scale);
    	int padding = ((diagonal_line_length - min_side_length) >> 1) + 16;
    	int double_padding = (padding << 1);
    	Mat model_padded = Mat(model.rows + double_padding, model.cols + double_padding, model.type(), Scalar::all(0));
    	model.copyTo(model_padded(Rect(padding, padding, model.cols, model.rows)));
    	Mat mask_padded = Mat(mask.rows + double_padding, mask.cols + double_padding, mask.type(), Scalar::all(0));
    	mask.copyTo(mask_padded(Rect(padding, padding, mask.cols, mask.rows)));
    	model = model_padded;
    	mask = mask_padded;
    }
    
    vector<ShapeInfo> KcgMatch::ProduceShapeInfos(AngleRange angle_range, ScaleRange scale_range) {
    
    	assert(scale_range.begin > KCG_EPS && scale_range.end > KCG_EPS);
    	assert(angle_range.end >= angle_range.begin);
    	assert(scale_range.end >= scale_range.begin);
    	assert(angle_range.step > KCG_EPS);
    	assert(scale_range.step > KCG_EPS);
    	vector<ShapeInfo> shape_infos;
    	shape_infos.clear();
    	for (float scale = scale_range.begin; scale <= scale_range.end + KCG_EPS; scale += scale_range.step) {
    
    		for (float angle = angle_range.begin; angle <= angle_range.end + KCG_EPS; angle += angle_range.step) {
    
    			ShapeInfo info;
    			info.angle = angle;
    			info.scale = scale;
    			shape_infos.push_back(info);
    		}
    	}
    	return shape_infos;
    }
    
    Mat KcgMatch::Transform(Mat src, float angle, float scale) {
    
    	Mat dst;
    	Point center(src.cols / 2, src.rows / 2);
    	Mat rot_mat = cv::getRotationMatrix2D(center, angle, scale);
    	warpAffine(src, dst, rot_mat, src.size());
    	return dst;
    }
    
    Mat KcgMatch::MdlOf(Mat model, ShapeInfo info) {
    
    	return Transform(model, info.angle, info.scale);
    }
    
    Mat KcgMatch::MskOf(Mat mask, ShapeInfo info) {
    
    	return (Transform(mask, info.angle, info.scale) > 0);
    }
    
    void KcgMatch::DrawTemplate(Mat &image, Template templ, Scalar color) {
    
    	for (int i = 0; i < templ.features.size(); i++) {
    
    		auto feature = templ.features[i];
    		line(image,
    			Point(templ.x + feature.x, templ.y + feature.y),
    			Point(templ.x + feature.x, templ.y + feature.y),
    			color, 1);
    	}
    }
    
    void KcgMatch::QuantifyEdge(Mat image, Mat &angle, Mat &quantized_angle, Mat &mag, float mag_thresh, bool calc_180) {
    
    	Mat dx, dy;
    	//Sobel(image, dx, CV_32F, 1, 0, 3, 1.0, 0.0, BORDER_REPLICATE);
    	//Sobel(image, dy, CV_32F, 0, 1, 3, 1.0, 0.0, BORDER_REPLICATE);
    	float mask_x[3][3] = { { -1,0,1 },{ -2,0,2 },{ -1,0,1 } };
    	float mask_y[3][3] = { { 1,2,1 },{ 0,0,0 },{ -1,-2,-1 } };
    	Mat kernel_x = Mat(3, 3, CV_32F, mask_x);
    	Mat kernel_y = Mat(3, 3, CV_32F, mask_y);
    	filter2D(image, dx, CV_32F, kernel_x);
    	filter2D(image, dy, CV_32F, kernel_y);
    	//dx = abs(dx);
    	//dy = abs(dy);
    	mag = dx.mul(dx) + dy.mul(dy);
    	phase(dx, dy, angle, true);
    
    	if(calc_180)
    		Quantify180(angle, quantized_angle, mag, mag_thresh);
    	else
    		Quantify8(angle, quantized_angle, mag, mag_thresh);
    }
    
    void KcgMatch::Quantify8(Mat angle, Mat &quantized_angle, Mat mag, float mag_thresh) {
    
    	Mat_<unsigned char> quantized_unfiltered;
    	angle.convertTo(quantized_unfiltered, CV_8U, 16.0f / 360.0f);
    	for (int r =0 ; r < angle.rows; ++r)
    	{
    		unsigned char *quant_ptr = quantized_unfiltered.ptr<unsigned char>(r);
    		for (int c = 0; c < angle.cols; ++c)
    		{
    			quant_ptr[c] &= 7;
    		}
    	}
    	//quantized_unfiltered.copyTo(quantized_angle);
    	quantized_angle = Mat::zeros(angle.size(), CV_8U);
    	for (int r = 0; r < quantized_angle.rows; ++r) {
    
    		quantized_angle.ptr<unsigned char>(r)[0] = 255;
    		quantized_angle.ptr<unsigned char>(r)[quantized_angle.cols - 1] = 255;
    	}
    	for (int c = 0; c < quantized_angle.cols; ++c) {
    
    		quantized_angle.ptr<unsigned char>(0)[c] = 255;
    		quantized_angle.ptr<unsigned char>(quantized_angle.rows - 1)[c] = 255;
    	}
    		
    	for (int r = 1; r < angle.rows - 1; ++r)
    	{
    		float *mag_ptr= mag.ptr<float>(r);
    		for (int c = 1; c < angle.cols - 1; ++c)
    		{
    			if (mag_ptr[c] >= (mag_thresh * mag_thresh))
    			{
    				int histogram[8] = { 0, 0, 0, 0, 0, 0, 0, 0 };
    
    				unsigned char *patch3x3_row = &quantized_unfiltered(r - 1, c - 1);
    				histogram[patch3x3_row[0]]++;
    				histogram[patch3x3_row[1]]++;
    				histogram[patch3x3_row[2]]++;
    
    				patch3x3_row += quantized_unfiltered.step1();
    				histogram[patch3x3_row[0]]++;
    				histogram[patch3x3_row[1]]++;
    				histogram[patch3x3_row[2]]++;
    
    				patch3x3_row += quantized_unfiltered.step1();
    				histogram[patch3x3_row[0]]++;
    				histogram[patch3x3_row[1]]++;
    				histogram[patch3x3_row[2]]++;
    
    				// Find bin with the most votes from the patch
    				int max_votes = 0;
    				int index = -1;
    				for (int i = 0; i < 8; ++i)
    				{
    					if (max_votes < histogram[i])
    					{
    						index = i;
    						max_votes = histogram[i];
    					}
    				}
    
    				// Only accept the quantization if majority of pixels in the patch agree
    				static const int NEIGHBOR_THRESHOLD = 5;
    				if (max_votes >= NEIGHBOR_THRESHOLD)
    					quantized_angle.at<unsigned char>(r, c) = index;
    				else
    					quantized_angle.at<unsigned char>(r, c) = 255;
    			}
    			else
    			{
    				quantized_angle.at<unsigned char>(r, c) = 255;
    			}
    		}
    	}
    }
    
    void KcgMatch::Quantify180(Mat angle, Mat &quantized_angle, Mat mag, float mag_thresh) {
    
    	quantized_angle = Mat::zeros(angle.size(), CV_8U);
    	#pragma omp parallel for
    	for (int r = 0; r < angle.rows; ++r)
    	{
    		unsigned char *quantized_angle_ptr = quantized_angle.ptr<unsigned char>(r);
    		float *angle_ptr = angle.ptr<float>(r);
    		float *mag_ptr = mag.ptr<float>(r);
    		for (int c = 0; c < angle.cols; ++c)
    		{
    			if (mag_ptr[c] >= (mag_thresh * mag_thresh))
    				quantized_angle_ptr[c] = (int)round(angle_ptr[c]) % 180;
    			else
    				quantized_angle_ptr[c] = 255;
    		}
    	}
    }
    
    Template KcgMatch::ExtractTemplate(Mat angle, Mat quantized_angle, Mat mag, ShapeInfo shape_info, 
    	PyramidLevel pl, float weak_thresh, float strong_thresh, int num_features, Mat mask) {
    
    	Mat local_angle = Mat(angle.size(), angle.type());
    	for (int r = 0; r < angle.rows; ++r) {
    
    		float *angle_ptr = angle.ptr<float>(r);
    		float *local_angle_ptr = local_angle.ptr<float>(r);
    		for (int c = 0; c < angle.cols; ++c) {
    
    			float dir = angle_ptr[c];
    			if ((dir > 0. && dir < 22.5) || (dir > 157.5 && dir < 202.5) || (dir > 337.5 && dir < 360.))
    				local_angle_ptr[c] = 0.f;
    			else if ((dir > 22.5 && dir < 67.5) || (dir > 202.5 && dir < 247.5))
    				local_angle_ptr[c] = 45.f;
    			else if ((dir > 67.5 && dir < 112.5) || (dir > 247.5 && dir < 292.5))
    				local_angle_ptr[c] = 90.f;
    			else if ((dir > 112.5 && dir < 157.5) || (dir > 292.5 && dir < 337.5))
    				local_angle_ptr[c] = 135.f;
    			else
    				local_angle_ptr[c] = 0.f;
    		}
    	}
    
    	vector<Candidate> candidates;
    	candidates.clear();
    	bool no_mask = mask.empty();
    	float weak_sq = weak_thresh * weak_thresh;
    	float strong_sq = strong_thresh * strong_thresh;
    	float pre_grad, lst_grad;
    	for (int r = 1; r < mag.rows - 1; ++r)
    	{
    		const unsigned char *mask_ptr = no_mask ? NULL : mask.ptr<unsigned char>(r);
    		const float* pre_ptr = mag.ptr<float>(r - 1);
    		const float* cur_ptr = mag.ptr<float>(r);
    		const float* lst_ptr = mag.ptr<float>(r + 1);
    		float *local_angle_ptr = local_angle.ptr<float>(r);
    
    		for (int c = 1; c < mag.cols - 1; ++c)
    		{
    			if (no_mask || mask_ptr[c])
    			{
    				switch ((int)local_angle_ptr[c]) {
    
    				case 0:
    					pre_grad = cur_ptr[c - 1];
    					lst_grad = cur_ptr[c + 1];
    					break;
    				case 45:
    					pre_grad = pre_ptr[c + 1];
    					lst_grad = lst_ptr[c - 1];
    					break;
    				case 90:
    					pre_grad = pre_ptr[c];
    					lst_grad = lst_ptr[c];
    					break;
    				case 135:
    					pre_grad = pre_ptr[c - 1];
    					lst_grad = lst_ptr[c + 1];
    					break;
    				}
    				if ((cur_ptr[c] > pre_grad) && (cur_ptr[c] > lst_grad)) {
    
    					float score = cur_ptr[c];
    					bool validity = false;
    					if (score >= weak_sq) {
    
    						if (score >= strong_sq) {
    
    							validity = true;
    						}
    						else {
    
    							if (((pre_ptr[c - 1])	>= strong_sq) ||
    								((pre_ptr[c])		>= strong_sq) ||
    								((pre_ptr[c + 1])	>= strong_sq) ||
    								((cur_ptr[c - 1])	>= strong_sq) ||
    								((cur_ptr[c + 1])	>= strong_sq) ||
    								((lst_ptr[c - 1])	>= strong_sq) ||
    								((lst_ptr[c])		>= strong_sq) ||
    								((lst_ptr[c + 1])	>= strong_sq))
    							{
    								validity = true;
    							}
    						}
    					}
    					if (validity == true && 
    						quantized_angle.at<unsigned char>(r, c) != 255) {
    
    						Candidate cd;
    						cd.score = score;
    						cd.feature.x = c;
    						cd.feature.y = r;
    						cd.feature.lbl = quantized_angle.at<unsigned char>(r, c);
    						candidates.push_back(cd);
    					}
    				}
    
    			}
    		}
    	}
    
    	Template templ;
    	templ.shape_info.angle = shape_info.angle;
    	templ.shape_info.scale = shape_info.scale;
    	templ.pyramid_level = pl;
    	templ.is_valid = 0;
    	templ.features.clear();
    
    	if (candidates.size() >= num_features && num_features > 0) {
    
    		std::stable_sort(candidates.begin(), candidates.end());
    		float distance = static_cast<float>(candidates.size() / num_features + 1);
    		templ = SelectScatteredFeatures(candidates, num_features, distance);
    	}
    	else {
    
    		for (int c = 0; c < candidates.size(); c++) {
    
    			templ.features.push_back(candidates[c].feature);
    		}
    	}
    	
    	if (templ.features.size() > 0) {
    
    		templ.is_valid = 1;
    		CropTemplate(templ);
    	}
    
    	return templ;
    }
    
    Template KcgMatch::SelectScatteredFeatures(vector<Candidate> candidates, int num_features, float distance) {
    
    	Template templ;
    	templ.features.clear();
    	float distance_sq = distance * distance;
    	int i = 0;
    	while (templ.features.size() < num_features) {
    
    		Candidate c = candidates[i];
    		// Add if sufficient distance away from any previously chosen feature
    		bool keep = true;
    		for (int j = 0; (j < (int)templ.features.size()) && keep; ++j)
    		{
    			Feature f = templ.features[j];
    			keep = ((c.feature.x - f.x) * (c.feature.x - f.x) + (c.feature.y - f.y) * (c.feature.y - f.y) >= distance_sq);
    		}
    		if (keep)
    			templ.features.push_back(c.feature);
    
    		if (++i == (int)candidates.size())
    		{
    			// Start back at beginning, and relax required distance
    			i = 0;
    			distance -= 1.0f;
    			distance_sq = distance * distance;
    			// if (distance < 3)
    			// {
    			//     // we don't want two features too close
    			//     break;
    			// }
    		}
    	}
    	return templ;
    }
    
    Rect KcgMatch::CropTemplate(Template &templ) {
    
    	int min_x = std::numeric_limits<int>::max();
    	int min_y = std::numeric_limits<int>::max();
    	int max_x = std::numeric_limits<int>::min();
    	int max_y = std::numeric_limits<int>::min();
    
    	// First pass: find min/max feature x,y 
    	for (int i = 0; i < (int)templ.features.size(); ++i)
    	{
    		int x = templ.features[i].x;
    		int y = templ.features[i].y;
    		min_x = std::min(min_x, x);
    		min_y = std::min(min_y, y);
    		max_x = std::max(max_x, x);
    		max_y = std::max(max_y, y);
    	}
    	
    	/// @todo Why require even min_x, min_y?
    	if (min_x % 2 == 1)
    		--min_x;
    	if (min_y % 2 == 1)
    		--min_y;
    
    	// Second pass: set width/height and shift all feature positions
    	templ.w = (max_x - min_x);
    	templ.h = (max_y - min_y);
    	templ.x = min_x;
    	templ.y = min_y;
    
    	for (int i = 0; i < (int)templ.features.size(); ++i)
    	{
    		templ.features[i].x -= templ.x;
    		templ.features[i].y -= templ.y;
    	}
    	return Rect(min_x, min_y, max_x - min_x, max_y - min_y);
    }
    
    void KcgMatch::LoadRegion8Idxes() {
    
    	int keys[16] = { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 };
    	region8_idxes_.clear();
    	int angle_region = (int)((angle_range_.end - angle_range_.begin) / angle_range_.step) + 1;
    	int scale_region = (int)((scale_range_.end - scale_range_.begin) / scale_range_.step) + 1;
    	for (int ar = 0; ar < angle_region; ar++) {
    
    		float cur_agl = templ_all_[PyramidLevel_0][ar].shape_info.angle;
    		if (cur_agl < 0.f) cur_agl += 360.f;
    		int idx = 0;
    		for (int i = 0; i < 16; i++) {
    
    			if (cur_agl >= AngleRegionTable[i][0] &&
    				cur_agl < AngleRegionTable[i][1]) {
    
    				idx = i;
    				break;
    			}
    		}
    		if (keys[idx] == 0) {
    
    			for (int sr = 0; sr < scale_region; sr++) {
    
    				region8_idxes_.push_back(ar + sr * angle_region);
    			}
    		}
    		keys[idx] = 1;
    	}
    }
    
    void KcgMatch::SaveModel() {
    
    	int total_templ = 0;
    	for (int i = 0; i < PyramidLevel_TabooUse; i++) {
    
    		total_templ += (int)templ_all_[i].size();
    	}
    	assert((total_templ / PyramidLevel_TabooUse) == templ_all_[0].size());
    	int match_range_size = (int)templ_all_[0].size();
    	string model_name = model_root_ + class_name_ + KCG_MODEL_SUFFUX;
    	FileStorage fs(model_name, FileStorage::WRITE);
    	fs << "class_name" << class_name_;
    	fs << "total_pyramid_levels" << PyramidLevel_7;
    	fs << "angle_range_bgin" << angle_range_.begin;
    	fs << "angle_range_end" << angle_range_.end;
    	fs << "angle_range_step" << angle_range_.step;
    	fs << "scale_range_bgin" << scale_range_.begin;
    	fs << "scale_range_end" << scale_range_.end;
    	fs << "scale_range_step" << scale_range_.step;
    	fs << "templates"
    	<< "[";
    	{
    		for (int i = 0; i < match_range_size; i++) {
    
    			fs << "{";
    			fs << "template_id" << int(i);
    			fs << "template_pyrds"
    			<< "[";
    			{
    				for (int j = 0; j < PyramidLevel_TabooUse; j++) {
    
    					auto templ = templ_all_[j][i];
    					fs << "{";
    					fs << "id" << int(i);
    					fs << "pyramid_level" << templ.pyramid_level;
    					fs << "is_valid" << templ.is_valid;
    					fs << "x" << templ.x;
    					fs << "y" << templ.y;
    					fs << "w" << templ.w;
    					fs << "h" << templ.h;
    					fs << "shape_scale" << templ.shape_info.scale;
    					fs << "shape_angle" << templ.shape_info.angle;
    					fs << "feature_size" << (int)templ.features.size();
    					fs << "features"
    					<< "[";
    					{
    						for (int k = 0; k < (int)templ.features.size(); k++) {
    
    							auto feat = templ.features[k];
    							fs << "[:" << feat.x << feat.y << feat.lbl << "]";
    						}
    					}
    					fs << "]";
    					fs << "}";
    				}
    			}
    			fs << "]";
    			fs << "}";
    		}
    	}
    	fs << "]";
    }
    
    void KcgMatch::LoadModel() {
    
    	ClearModel();
    	string model_name = model_root_ + class_name_ + KCG_MODEL_SUFFUX;
    	FileStorage fs(model_name, FileStorage::READ);
    	assert(fs.isOpened() && "load model failed.");
    	FileNode fn = fs.root();
    	angle_range_.begin = fn["angle_range_bgin"];
    	angle_range_.end = fn["angle_range_end"];
    	angle_range_.step = fn["angle_range_step"];
    	scale_range_.begin = fn["scale_range_bgin"];
    	scale_range_.end = fn["scale_range_end"];
    	scale_range_.step = fn["scale_range_step"];
    	FileNode tps_fn = fn["templates"];
    	FileNodeIterator tps_it = tps_fn.begin(), tps_it_end = tps_fn.end();
    	for (; tps_it != tps_it_end; ++tps_it)
    	{
    		int template_id = (*tps_it)["template_id"];
    		FileNode pyrds_fn = (*tps_it)["template_pyrds"];
    		FileNodeIterator pyrd_it = pyrds_fn.begin(), pyrd_it_end = pyrds_fn.end();
    		int pl = 0;
    		for (; pyrd_it != pyrd_it_end; ++pyrd_it)
    		{
    			FileNode pyrd_fn = (*pyrd_it);
    			Template templ;
    			templ.id = pyrd_fn["id"];
    			templ.pyramid_level = pyrd_fn["pyramid_level"];
    			templ.is_valid = pyrd_fn["is_valid"];
    			templ.x = pyrd_fn["x"];
    			templ.y = pyrd_fn["y"];
    			templ.w = pyrd_fn["w"];
    			templ.h = pyrd_fn["h"];
    			templ.shape_info.scale = pyrd_fn["shape_scale"];
    			templ.shape_info.angle = pyrd_fn["shape_angle"];
    			FileNode features_fn = pyrd_fn["features"];
    			FileNodeIterator feature_it = features_fn.begin(), feature_it_end = features_fn.end();
    			for (; feature_it != feature_it_end; ++feature_it)
    			{
    				FileNode feature_fn = (*feature_it);
    				FileNodeIterator feature_info = feature_fn.begin();
    				Feature feat;
    				feature_info >> feat.x >> feat.y >> feat.lbl;
    				templ.features.push_back(feat);
    			}
    			templ_all_[pl].push_back(templ);
    			pl++;
    		}
    	}
    
    	LoadRegion8Idxes();
    }
    
    void KcgMatch::ClearModel() {
    
    	for (int i = 0; i < PyramidLevel_TabooUse; i++) {
    
    		templ_all_[i].clear();
    	}
    }
    
    void KcgMatch::InitMatchParameter(float score_thresh, float overlap, float mag_thresh, float greediness, int T, int top_k, MatchingStrategy strategy) {
    
    	score_thresh_ = score_thresh;
    	overlap_ = overlap;
    	mag_thresh_ = mag_thresh;
    	greediness_ = greediness;
    	T_ = T;
    	top_k_ = top_k;
    	strategy_ = strategy;
    }
    
    void KcgMatch::GetAllPyramidLevelValidSource(cv::Mat &source, PyramidLevel pyrd_level) {
    
    	sources_.clear();
    	for (int pl = 0; pl <= pyrd_level; pl++) {
    
    		Mat source_pyrd;
    		if (pl == 0) source_pyrd = source;
    		else pyrDown(source, source_pyrd, Size(source.cols >> 1, source.rows >> 1));
    		source = source_pyrd;
    		sources_.push_back(source_pyrd);
    	}
    }
    
    vector<Match> KcgMatch::GetTopKMatches(vector<Match> matches) {
    
    	vector<Match> top_k_matches;
    	top_k_matches.clear();
    	if (top_k_ > 0 && (top_k_ < matches.size()) && (matches.size() > 0)) {
    
    		int k = 0;
    		top_k_matches.push_back(matches[0]);
    		for (int m = 1; m < matches.size(); m++) {
    
    			if (matches[m].similarity < matches[m - 1].similarity) {
    
    				++k;
    				if(k >= top_k_) break;
    			}
    			top_k_matches.push_back(matches[m]);
    		}
    	}
    	else
    	{
    		top_k_matches = matches;
    	}
    	return top_k_matches;
    }
    
    vector<Match> KcgMatch::DoNmsMatches(vector<Match> matches, PyramidLevel pl, float overlap) {
    
    	vector<Rect> boxes; boxes.clear();
    	vector<float> scores; scores.clear();
    	vector<int> indices; indices.clear();
    	for (int m = 0; m < matches.size(); m++) {
    
    		auto templ = templ_all_[pl][matches[m].template_id];
    		Rect box = Rect(matches[m].x, matches[m].y, templ.w, templ.h);
    		boxes.insert(boxes.end(), box);
    		scores.insert(scores.end(), matches[m].similarity);
    	}
    	cv_dnn_nms::NMSBoxes(boxes, scores, overlap, overlap, indices);
    	vector<Match> final_matches; final_matches.clear();
    	for (auto index : indices) {
    
    		final_matches.push_back(matches[index]);
    	}
    	return final_matches;
    }
    
    vector<Match> KcgMatch::MatchingPyrd180(Mat src, PyramidLevel pl, vector<int> region_idxes) {
    
    	pl = PyramidLevel(pl + 8);
    	vector<Match> matches; matches.clear();
    	Mat angle, quantized_angle, mag;
    	QuantifyEdge(src, angle, quantized_angle, mag, mag_thresh_, true);
    	#pragma omp parallel 
    	{
    		int tlsz = region_idxes.empty() ? ((int)templ_all_[pl].size()) : ((int)region_idxes.size());
    		#pragma omp for nowait
    		for (int t = 0; t < tlsz; t++) {
    
    			Template templ = region_idxes.empty() ? (templ_all_[pl][t]) : (templ_all_[pl][region_idxes[t]]);
    			for (int r = 0; r < quantized_angle.rows - templ.h; r++) {
    
    				for (int c = 0; c < quantized_angle.cols - templ.w; c++) {
    
    					int fsz = (int)templ.features.size();
    					float partial_sum = 0.f;
    					bool valid = true;
    					for (int f = 0; f < fsz; f++) {
    
    						Feature feat = templ.features[f];
    						int sidx = quantized_angle.ptr<unsigned char>(r + feat.y)[c + feat.x];
    						int tidx = feat.lbl;
    						if (sidx != 255) {
    
    							partial_sum += score_table_[sidx][tidx];
    						}
    						if (partial_sum + (fsz - f) * greediness_ < score_thresh_ * fsz) {
    
    							valid = false;
    							break;
    						}
    					}
    					if (valid) {
    
    						float score = partial_sum / fsz;
    						if (score >= score_thresh_) {
    
    							Match match;
    							match.x = c;
    							match.y = r;
    							match.similarity = score;
    							match.template_id = templ.id;
    							#pragma omp critical
    							matches.insert(matches.end(), match);
    						}
    					}
    
    				}
    			}
    		}
    	}
    	matches = DoNmsMatches(matches, pl, overlap_);
    	return matches;
    }
    
    vector<Match> KcgMatch::MatchingPyrd8(Mat src, PyramidLevel pl, vector<int> region_idxes) {
    
    	vector<Match> matches; matches.clear();
    	Mat angle, quantized_angle, mag;
    	QuantifyEdge(src, angle, quantized_angle, mag, mag_thresh_, false);
    	Mat spread_angle;
    	Spread(quantized_angle, spread_angle, T_);
    	vector<Mat> response_maps;
    	ComputeResponseMaps(spread_angle, response_maps);
    	#pragma omp parallel 
    	{
    		int tlsz = region_idxes.empty() ? ((int)templ_all_[pl].size()) : ((int)region_idxes.size());
    		#pragma omp for nowait
    		for (int t = 0; t < tlsz; t++) {
    
    			Template templ = region_idxes.empty() ? (templ_all_[pl][t]) : (templ_all_[pl][region_idxes[t]]);
    			for (int r = 0; r < quantized_angle.rows - templ.h; r += T_) {
    
    				for (int c = 0; c < quantized_angle.cols - templ.w; c += T_) {
    
    					int fsz = (int)templ.features.size();
    					int partial_sum = 0;
    					bool valid = true;
    					for (int f = 0; f < fsz; f++) {
    
    						Feature feat = templ.features[f];
    						int label = feat.lbl;
    						partial_sum += 
    							response_maps[label].ptr<unsigned char>(r + feat.y)[c + feat.x];
    						if (partial_sum + (fsz - f) * greediness_ < score_thresh_ * fsz) {
    
    							valid = false;
    							break;
    						}
    					}
    					if (valid) {
    
    						float score = partial_sum / (100.f * fsz);
    						if (score >= score_thresh_) {
    
    							Match match;
    							match.x = c;
    							match.y = r;
    							match.similarity = score;
    							match.template_id = templ.id;
    							#pragma omp critical
    							matches.insert(matches.end(), match);
    						}
    					}
    				}
    			}
    		}
    	}
    	matches = DoNmsMatches(matches, pl, overlap_);
    	return matches;
    }
    
    void KcgMatch::Spread(const Mat quantized_angle, Mat &spread_angle, int T) {
    
    	spread_angle = Mat::zeros(quantized_angle.size(), CV_8U);
    	int cols = quantized_angle.cols;
    	int rows = quantized_angle.rows;
    	int half_T = 0;
    	if (T != 1) half_T = T / 2;
    	#pragma omp parallel for
    	for (int r = half_T; r < rows - half_T; r++) {
    
    		for (int c = half_T; c < cols - half_T; c++) {
    
    			for (int i = -half_T; i <= half_T; i++) {
    
    				for (int j = -half_T; j <= half_T; j++) {
    
    					unsigned char shift_bits =
    						quantized_angle.ptr<unsigned char>(r + i)[c + j];
    					if (shift_bits < 8) {
    
    						spread_angle.ptr<unsigned char>(r)[c] |=
    							(unsigned char)(1 << shift_bits);
    					}
    				}
    			}
    		}
    	}
    }
    
    void KcgMatch::ComputeResponseMaps(const Mat spread_angle, vector<Mat> &response_maps) {
    
    	response_maps.clear();
    	for (int i = 0; i < 8; i++) {
    
    		Mat rm;
    		rm.create(spread_angle.size(), CV_8U);
    		response_maps.push_back(rm);
    	}
    	int cols = spread_angle.cols;
    	int rows = spread_angle.rows;
    	#pragma omp parallel for
    	for (int i = 0; i < 8; i++) {
    
    		for (int r = 0; r < rows; r++) {
    
    			for (int c = 0; c < cols; c++) {
    
    				response_maps[i].ptr<unsigned char>(r)[c] =
    					score_table_8map_[i][spread_angle.ptr<unsigned char>(r)[c]];
    			}
    		}
    	}
    }
    
    bool KcgMatch::CalcPyUpRoiAndStartPoint(PyramidLevel cur_pl, PyramidLevel obj_pl, Match match,
    	Mat &r, Point &p, bool is_padding) {
    
    	auto templ = templ_all_[cur_pl][match.template_id];
    	int padding = 0;
    	if (is_padding) {
    
    		int min_side = std::min(templ.w, templ.h);
    		int diagonal_line_length = (int)ceil(sqrt(templ.w*templ.w + templ.h*templ.h));
    		padding = diagonal_line_length - min_side;
    	}
    	int err_pl = cur_pl - obj_pl;
    	int T = 2 * T_;
    	int extend_pixel = 1;
    	cv::Point bp, ep;
    	int multiple = (1 << err_pl);
    	match.x -= (T + padding) / 2;
    	match.y -= (T + padding) / 2;
    	templ.w += (T + padding);
    	templ.h += (T + padding);
    	bp.x = (match.x - extend_pixel) * multiple;
    	bp.y = (match.y - extend_pixel) * multiple;
    	ep.x = (match.x + templ.w + extend_pixel) * multiple;
    	ep.y = (match.y + templ.h + extend_pixel) * multiple;
    	if (bp.x < 0) bp.x = 0;
    	if (bp.y < 0) bp.y = 0;
    	if (ep.x < 0) ep.x = 0;
    	if (ep.y < 0) ep.y = 0;
    	if (bp.x >= sources_[obj_pl].cols) bp.x = sources_[obj_pl].cols - 1;
    	if (bp.y >= sources_[obj_pl].rows) bp.y = sources_[obj_pl].rows - 1;
    	if (ep.x >= sources_[obj_pl].cols) ep.x = sources_[obj_pl].cols - 1;
    	if (ep.y >= sources_[obj_pl].rows) ep.y = sources_[obj_pl].rows - 1;
    	if (bp.x != ep.x || bp.y != ep.y) {
    
    		Rect rect = Rect(bp, ep);
    		Mat roi(sources_[obj_pl], rect);
    		r = roi;
    		p = bp;
    		return true;
    	}
    	else
    	{
    		return false;
    	}
    }
    
    void KcgMatch::CalcRegionIndexes(vector<int> &region_idxes, Match match, MatchingStrategy strategy) {
    
    	region_idxes.clear();
    	Template templ = templ_all_[PyramidLevel_0][match.template_id];
    	float match_agl = templ.shape_info.angle;
    	float match_sal = templ.shape_info.scale;
    	int angle_region = (int)((angle_range_.end - angle_range_.begin) / angle_range_.step) + 1;
    	int scale_region = (int)((scale_range_.end - scale_range_.begin) / scale_range_.step) + 1;
    	if (strategy <= Strategy_Middling) {
    
    		if (match_agl < 0.f) match_agl += 360.f;
    		int key = (int)floor(match_agl / 22.5f);
    		float left_agl = match_agl - key * 22.5f;
    		for (int ar = 0; ar < angle_region; ar++) {
    
    			float cur_agl = templ_all_[PyramidLevel_0][ar].shape_info.angle;
    			if (cur_agl < 0.f) cur_agl += 360.f;
    			int k = key;
    			if (cur_agl >= AngleRegionTable[k][0] && cur_agl < AngleRegionTable[k][1]) {
    
    				for (int sr = 0; sr < scale_region; sr++) {
    
    					region_idxes.push_back(ar + sr * angle_region);
    				}
    			}
    			if (strategy == Strategy_Accurate) {
    
    				if (left_agl < 11.25f) {
    
    					k = key - 1;
    					if (k < 0) k = 15;
    					if (cur_agl >= AngleRegionTable[k][0] && cur_agl < AngleRegionTable[k][1]) {
    
    						for (int sr = 0; sr < scale_region; sr++) {
    
    							region_idxes.push_back(ar + sr * angle_region);
    						}
    					}
    				}
    				else
    				{
    					k = key + 1;
    					if (k > 15) k = 0;
    					if (cur_agl >= AngleRegionTable[k][0] && cur_agl < AngleRegionTable[k][1]) {
    
    						for (int sr = 0; sr < scale_region; sr++) {
    
    							region_idxes.push_back(ar + sr * angle_region);
    						}
    					}
    				}
    			}
    		}
    	}
    	else if(strategy == Strategy_Rough) {
    
    		float err_range = 3.f;
    		for (int ar = 0; ar < angle_region; ar++) {
    
    			float cur_agl = templ_all_[PyramidLevel_0][ar].shape_info.angle;
    			if (cur_agl >= (match_agl - angle_range_.step * err_range) && 
    				cur_agl <= (match_agl + angle_range_.step * err_range)) {
    
    				for (int sr = 0; sr < scale_region; sr++) {
    
    					float cur_sal = templ_all_[PyramidLevel_0][ar + sr * angle_region].shape_info.scale;
    					if (cur_sal >= (match_sal - scale_range_.step * err_range) &&
    						cur_sal <= (match_sal + scale_range_.step * err_range)) {
    
    						region_idxes.push_back(ar + sr * angle_region);
    					}
    				}
    			}
    		}
    	}
    }
    
    vector<Match> KcgMatch::ReconfirmMatches(vector<Match> matches, PyramidLevel pl) {
    
    	vector<Match> rf_matches;
    	rf_matches.clear();
    	for (int i = 0; i < matches.size(); i++) {
    
    		Mat roi;
    		Point sp;
    		CalcPyUpRoiAndStartPoint(pl, pl, matches[i], roi, sp, true);
    		vector<int> region_idxes;
    		CalcRegionIndexes(region_idxes, matches[i], Strategy_Accurate);
    		auto tmp_matches = MatchingPyrd8(roi, pl, region_idxes);
    		if (tmp_matches.size() > 0) {
    
    			tmp_matches[0].x += sp.x;
    			tmp_matches[0].y += sp.y;
    			rf_matches.push_back(tmp_matches[0]);
    		}
    	}
    	rf_matches = DoNmsMatches(rf_matches, pl, overlap_);
    	return rf_matches;
    }
    
    vector<Match> KcgMatch::MatchingFinal(vector<Match> matches, PyramidLevel pl) {
    
    	vector<Match> final_matches;
    	final_matches.clear();
    	for (int i = 0; i < matches.size(); i++) {
    
    		Mat roi;
    		Point sp;
    		CalcPyUpRoiAndStartPoint(pl, PyramidLevel_0, matches[i], roi, sp, false);
    		vector<int> region_idxes;
    		CalcRegionIndexes(region_idxes, matches[i], strategy_);
    		auto tmp_matches = MatchingPyrd180(roi, PyramidLevel_0, region_idxes);
    		if (tmp_matches.size() > 0) {
    
    			tmp_matches[0].x += sp.x;
    			tmp_matches[0].y += sp.y;
    			final_matches.push_back(tmp_matches[0]);
    		}
    	}
    	final_matches = DoNmsMatches(final_matches, pl, overlap_);
    	return final_matches;
    }
    
    } // end namespace kcg_matching
    
    
    

    调用过程

    新建C++工程,将KcgMatch.h,KcgMatch.cpp添加到工程,添加main.cpp,在main.cpp里添加以下代码即可实现以上演示的功能。

    #include "KcgMatch.h"
    
    int main(int argc, char **argv) {
    	
    	// 实例化KcgMatch 
    	// "demo/k"为存储模板的根目录 
    	// "k"为模板的名字
    	kcg_matching::KcgMatch kcg("demo/k", "k");
    	// 读取模板图像
    	Mat model = imread("demo/k/template.png");
    	// 转灰度
    	cvtColor(model, model, COLOR_BGR2GRAY);
    	// 指定要制作的模板角度,尺度范围
    	kcg_matching::AngleRange ar(-180.f, 180.f, 10.f);
    	kcg_matching::ScaleRange sr(0.70f, 1.3f, 0.05f);
    	// 开始制作模板
    	kcg.MakingTemplates(model, ar, sr, 0, 30.f, 60.f);
    	
    	// 加载模板
    	cout << "Loading model ......" << endl;
    	kcg.LoadModel();
    	cout << "Load succeed." << endl;
    	
    	// 读取搜索图像
    	Mat source = imread("demo/k/search.png");
    	Mat draw_source;
    	source.copyTo(draw_source);
    	cvtColor(source, source, COLOR_BGR2GRAY);
    
    	Timer timer;
    	// 开始匹配
    	auto matches = 
    		kcg.Matching(source, 0.80f, 0.1f, 30.f, 0.9f,
    			kcg_matching::PyramidLevel_2, 2, 12); 
    	double t = timer.out("=== Match time ===");
    	cout << "Final match size: " << matches.size() << endl << endl;
    
    	// 画出匹配结果
    	kcg.DrawMatches(draw_source, matches, Scalar(255, 0, 0));
    
    	// 画出匹配时间
    	rectangle(draw_source, Rect(Point(0, 0), Point(136, 20)), Scalar(255, 255, 255), -1);
    	cv::putText(draw_source,
    				"time: " + to_string(t) + "s",
    				Point(0, 16), FONT_HERSHEY_PLAIN, 1.f, Scalar(0, 0, 0), 1);
    				
    	// 显示结果图像
    	namedWindow("draw_source", 0);
    	imshow("draw_source", draw_source);
    	waitKey(0);
    	system("pause");
    }
    

    后话

    当然,里面有很多细节代码,暂时没时间介绍(待续···),算法原来还是很简单的,主要就是速度的问题,因为模板个数太多了,还要滑动窗口来遍历图像,实在是很耗时间。最初没有任何优化的代码匹配这样的图片完整傻乎乎的跑下来需要1分钟左右,现在经过一些骚炒作(优化加速)只需要200ms,主要的优化过程就是查表,并行,量化,金字塔等等。还可以考虑用CPU指令集加速(sse avx)等。有空的话再一起来细扣(待续···)。

    码代码很累,且行且珍惜,如果觉得对您有点帮助,求各位大侠赏赐两包拉条,谢谢。

    参考资料

    机器视觉算法与应用

    Edge Based Template Matching

    知乎

    https://github.com/meiqua/shape_based_matching

     

     

     

     

     

    展开全文
  • 该存储库包含论文”中形状匹配管道的易于阅读的C ++ / Python实现。 生成线图 为了生成线图,我入侵了软件,将3D网格文件用作命令行输入,并生成了一组图像作为输出。 该图像集合显示了使用算法从不同视图渲染的输入...
  • 用opencv编写的形状匹配算法,但不具旋转和缩放功能。 https://www.codeproject.com/KB/graphics/Edge_Based_template_match/GeoMatch_src.zip https://www.codeproject.com/Articles/99457/E...
  • Halcon形状匹配

    千次阅读 2018-09-11 10:55:36
    去年有过一段时间的集中学习,做了许多的练习和实验,并对基于HDevelop的形状匹配算法的参数优化进行了研究,写了一篇《基于HDevelop的形状匹配算法参数的优化研究》文章,总结了在形状匹配过程中哪些参数影响到模板...
  • 使用高光谱图像中的导数的光谱曲线形状匹配
  • 机器视觉(一)2D形状匹配

    千次阅读 2019-11-08 16:28:16
    常见的匹配算法有基于灰度的匹配,基于边缘的匹配,基于形状匹配等。推荐《机器视觉算法与应用》这本书,这本书是halcon开发人员撰写,对于模板匹配介绍的比较详细。下面介绍基于形状的模板匹配。 在pcl和opencv...
  • 因此,为了提高道路网络中的地图匹配精度,提出了基于Zernike形状矩的地图匹配算法。新算法引入Zernike矩描述轨迹曲线的形状,进一步修正了错误结果。通过仿真和实验表明,新算法在复杂环境下具有较强的有效性和可靠性。
  • 这是郭小元维护的论文“具有邻近点匹配和基于局部形状的强度分析的簇状核分割”的重新实现。 表中的内容 -[数据集](用于测试的图像补丁) * ( Image patches: 20, random size, 1024*1376 pixels ) -[用法](如何...
  • 本文以螺丝的分拣为例,设计一种基于形状匹配的螺丝识别算法 一、算法设计 流水线工作环境应处于明亮的室内,流水线的输送皮带颜色应该是暗色系,黑色最佳,相机摄取到的图片如图 对其进行阈值分
  • 算子: ...(1)参数1是待匹配轮廓或者灰度图像 (2)参数2同参数1 (3)比较参数1和2相似度的方法,opencv提供了三种如下: CV_CONTOURS_MATCH_I1 CV_CONTOURS_MATCH_I2 CV_CONTOURS_MATCH_I3 (4)
  • OpenCV基于形状的模板匹配

    千次阅读 2020-12-11 14:51:53
    OpenCV基于形状的模板匹配引言基于形状的匹配算法具体代码KcgMatch.hKcgMatch....而在halcon中有个基于形状匹配的算子,这个算子非常好用,随便截取一个ROI区域做模板就可以在搜索图像中匹配到相似的区域,并且能输
  • 一个Matlab项目,实现了具有自适应窗口形状的Kanade立体匹配算法。 文件索引: Benchmark1.m-关于9种以上窗口形状的初始工作 Benchmark3.m-我们算法的主要实现。 它使用以下功能: 不确定度.m:在给定窗口内计算不...
  • halcon 形状匹配

    万次阅读 多人点赞 2018-10-31 14:01:57
    如果检测对象几乎是对称的,那么就应该限制角度范围,否则模板会在不同的角度上找到多个最匹配的模板,十字形状或者正方形形状必须将角度范围限定在90°,对于长方形,应该限制在180°,圆的话应该是0°。...
  • 基于形状的模板匹配cpp,由halcon程序编写后导出
  • Opencv之利用matchshape算子实现简单的形状匹配

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 33,885
精华内容 13,554
关键字:

形状匹配代码