1. 程式人生 > >將多個凸包連線起來形成一個大凸包找到最接近的輪廓

將多個凸包連線起來形成一個大凸包找到最接近的輪廓

目的:對於很暗很暗的影象 用普通方法不好一下找到輪廓 現在需要找到其最接近的輪廓

像這樣的影象  拍的是一塊完整的石頭  藍色部分只是石頭的反光部分而已  現在是要找到這個石頭的輪廓  最終得到它的面積。

#include<opencv2/opencv.hpp>
#include<opencv2/imgproc/imgproc.hpp>
#include<opencv2/highgui/highgui.hpp>
#include <afxwin.h>
#include<fstream>
using namespace std;
using namespace cv;

Mat my_contour(IplImage* dst, IplImage bimgipl)
{
	cvZero(dst);
	CvMemStorage *storage = cvCreateMemStorage();
	CvSeq *contour = NULL, *hull = NULL;
	vector<CvPoint> allpoints;
	CvContourScanner scanner = cvStartFindContours(&bimgipl, storage);
	while ((contour = cvFindNextContour(scanner)) != NULL)
	{
		cvDrawContours(dst, contour, cv::Scalar(255), cv::Scalar(255), 0);
		hull = cvConvexHull2(contour, 0, CV_CLOCKWISE, 0);
		CvPoint pt0 = **(CvPoint**)cvGetSeqElem(hull, hull->total - 1);
		allpoints.push_back(pt0);
		for (int i = 0; i < hull->total; ++i)
		{
			CvPoint pt1 = **(CvPoint**)cvGetSeqElem(hull, i);
			allpoints.push_back(pt1);
			cvLine(dst, pt0, pt1, cv::Scalar(255));
			pt0 = pt1;
		}
	}
	//上面是第一次尋找 可能找到了多個輪廓即有多個凸包 下面儲存每個凸包的凸點到容器中 然後將這些點連線起來 這樣就有凹陷的地方了 就可以在此基礎上尋找最後的凸包了 
	Mat dstmat(dst, 0);
	std::vector<std::vector<cv::Point>> myown_contours;
	cv::findContours(dstmat, myown_contours, CV_RETR_EXTERNAL, CV_CHAIN_APPROX_NONE);
	//判斷也許本來就只有一個凸包 就不用存點畫線再尋一次凸包了
	if (myown_contours.size() > 1)
	{
		for (int i = 0; i < allpoints.size() - 1; i++)
		{
			CvPoint firstdot = allpoints[i];
			CvPoint secdot = allpoints[i + 1];
			cvLine(dst, firstdot, secdot, cv::Scalar(255), 2);
		}
		CvContourScanner scanner2 = cvStartFindContours(dst, storage);
		while ((contour = cvFindNextContour(scanner2)) != NULL)
		{
			cvDrawContours(dst, contour, cv::Scalar(255), cv::Scalar(255), 0);
			hull = cvConvexHull2(contour, 0, CV_CLOCKWISE, 0);
			CvPoint pt0 = **(CvPoint**)cvGetSeqElem(hull, hull->total - 1);
			for (int i = 0; i < hull->total; ++i)
			{
				CvPoint pt1 = **(CvPoint**)cvGetSeqElem(hull, i);
				cvLine(dst, pt0, pt1, cv::Scalar(255));
				pt0 = pt1;
			}
		}
	}
	Mat bimgdst(dst, 0);
	std::vector<std::vector<cv::Point>> contours;
	cv::findContours(bimgdst, contours, CV_RETR_EXTERNAL, CV_CHAIN_APPROX_NONE);
	Mat contoursimg(bimgdst.size(), CV_8UC1, cv::Scalar(0));
	drawContours(contoursimg, contours, -1, Scalar(255), CV_FILLED);
	return contoursimg;
}

char filename[100];
char savename[100];

void main()
{
	TickMeter tm;
	tm.start();

	vector<Mat> bgrimgs;
	for (int i = 1; i <= 1657; i++)
	{
		sprintf(filename, "正常圖_壓縮\\%d.bmp", i);
		sprintf(savename, "my_result3\\%d.bmp", i);

		Mat img = imread(filename);
		//
		split(img, bgrimgs);
		Mat bimg(img.size(), CV_8UC1, cv::Scalar(0));
		bgrimgs[0].copyTo(bimg);
		//     
		threshold(bimg, bimg, 12, 255, CV_THRESH_BINARY/*|CV_THRESH_OTSU*/);
		//imshow("thresh result", bimg);
		//
		Mat element(4, 4, CV_8U, cv::Scalar(1));
		morphologyEx(bimg, bimg, cv::MORPH_OPEN, element);
		//imshow("dilate result", bimg);
		//前20行置0 去除上部分的水滴噪聲 後20行也置0 去除下面的噪聲
		/*  //這個置0沒有下面寫的置0好 用下面的
		IplImage bimgipl = bimg;
		cvSetImageROI(&bimgipl, cvRect(0, 0, bimg.size().width,20));
		IplImage* bimg_front = cvCreateImage(cvSize(bimg.size().width, 20), IPL_DEPTH_8U, 1);
		cvCopy(&bimgipl, bimg_front, 0);
		cvShowImage("front", bimg_front);

		cvSetImageROI(&bimgipl, cvRect(0, bimg.size().height-20, bimg.size().width, 20));
		IplImage* bimg_end = cvCreateImage(cvSize(bimg.size().width, 20), IPL_DEPTH_8U, 1);
		cvCopy(&bimgipl, bimg_end, 0);
		cvShowImage("end", bimg_end);

		cvZero(bimg_front);
		cvZero(bimg_end);
		*/
		for (int i = 0; i < 20; i++)
		{
			uchar* data = bimg.ptr<uchar>(i);
			for (int j = 0; j < bimg.size().width; j++)
				data[j] = 0;
		}
		for (int i = bimg.size().height - 19; i < bimg.size().height; i++)
		{
			uchar* data = bimg.ptr<uchar>(i);
			for (int j = 0; j < bimg.size().width; j++)
				data[j] = 0;
		}
		//凸包找輪廓 (我用了2次凸包 中間用線連線 方便第二次凸包 找到完整輪廓)
		IplImage mybimgipl = bimg;
		IplImage *mydst = cvCreateImage(cvGetSize(&mybimgipl), 8, 1);
		Mat myresult(mydst, 0);
		myresult = my_contour(mydst, mybimgipl);
		imwrite(savename, myresult);
	}
	tm.stop();
	cout << "count=" << tm.getCounter() << ",process time=" << tm.getTimeMilli() << " ms" << endl;

}
拿一種一張圖的中間結果來說:

分別是閾值化 開操作 第一次凸包 第一次凸包填充後的結果  

左圖就是對於第一次凸包連線後 再凸包後的結果   右邊就是對左圖填充後的結果   也是最終結果!
done!!!!今天星期五!好開心!!!!!!!!!!!!

關於opencv中找輪廓的 : http://www.cnblogs.com/nktblog/p/4027137.html   寫得很棒

////////////////////////////////

給學校我自習室的電腦上裝了opencv2.4.13 按照http://blog.csdn.net/dcrmg/article/details/51809614   

看到巴衛和奈奈生就開心。。。。。。

//////////////////////////////////////////////////////////////////////////////////////////////////

我公司的同事說我這個比較慢  他在網上 https://github.com/abatilo/QuickHull 找到一個快速找凸包的叫我應用進來:於是我重新又寫了個:

#include <algorithm>
#include <cmath>
#include <vector>
#include <chrono>
#include <cstdio>
#include <random>
#include <vector>

#include<opencv2/opencv.hpp>
#include<opencv2/imgproc/imgproc.hpp>
#include<opencv2/highgui/highgui.hpp>
using namespace std;
using namespace cv;

struct myPoint {
	int x;
	int y;
	bool operator<(const myPoint& pt) {

		if (x == pt.x) 
		{
			return (y < pt.y);
		}
		return (x < pt.x);
	}
};

// http://stackoverflow.com/questions/1560492/how-to-tell-whether-a-myPoint-is-to-the-right-or-left-side-of-a-line
// A and B are used to define a line
// C is the myPoint whose side we're trying to determine
int SideOfLine(const myPoint &P1, const myPoint &P2, const myPoint &P3) {
	return (P2.x - P1.x) * (P3.y - P1.y) - (P2.y - P1.y) * (P3.x - P1.x);
}

// https://en.wikipedia.org/wiki/Distance_from_a_myPoint_to_a_line#Line_defined_by_two_myPoints
float DistanceFromLine(const myPoint &P1, const myPoint &P2, const myPoint &P3) {
	return (std::abs((P2.y - P1.y) * P3.x - (P2.x - P1.x) * P3.y + P2.x * P1.y - P2.y * P1.x)
		/ std::sqrt((P2.y - P1.y) * (P2.y - P1.y) + (P2.x - P1.x) * (P2.x - P1.x)));
}

// http://stackoverflow.com/questions/13300904/determine-whether-myPoint-lies-inside-triangle
bool myPointInTriangle(const myPoint &p, const myPoint &p1, const myPoint &p2, const myPoint &p3) {
	float a = ((p2.y - p3.y) * (p.x - p3.x) + (p3.x - p2.x) * (p.y - p3.y)) / ((p2.y - p3.y) * (p1.x - p3.x) + (p3.x - p2.x) * (p1.y - p3.y));
	float b = ((p3.y - p1.y) * (p.x - p3.x) + (p1.x - p3.x) * (p.y - p3.y)) / ((p2.y - p3.y) * (p1.x - p3.x) + (p3.x - p2.x) * (p1.y - p3.y));
	float c = 1.0f - a - b;
	return (0.0f < a && 0.0f < b && 0.0f < c);
}

// http://www.cse.yorku.ca/~aaw/Hang/quick_hull/Algorithm.html
void FindHull(const std::vector<myPoint> &Sk, const myPoint P, const myPoint Q, std::vector<myPoint> &hullmyPoints) {
	// If Sk has no myPoint, then  return
	if (Sk.size() == 0) return;

	std::vector<myPoint> S0;
	std::vector<myPoint> S1;
	std::vector<myPoint> S2;

	// From the given set of myPoints in Sk, find farthest myPoint, say C, from segment PQ 
	float furthestDistance = 0.0f;
	myPoint C;
	for (const auto &pt : Sk) {
		float distance = DistanceFromLine(P, Q, pt);
		if (distance > furthestDistance) {
			furthestDistance = distance;
			C = pt;
		}
	}
	// Add myPoint C to convex hull at the location between P and Q 
	hullmyPoints.push_back(C);

	/*
	* Three myPoints P, Q, and C partition the remaining myPoints of Sk into 3 subsets: S0, S1, and S2
	* where S0 are myPoints inside triangle PCQ, S1are myPoints on the right side of the oriented
	* line from  P to C, and S2 are myPoints on the right side of the oriented line from C to Q.
	*/
	for (const auto &pt : Sk) {
		if (myPointInTriangle(pt, P, C, Q)) {
			S0.push_back(pt);
		}
		else if (0 < SideOfLine(P, C, pt)) {
			S1.push_back(pt);
		}
		else if (0 < SideOfLine(C, Q, pt)) {
			S2.push_back(pt);
		}
	}

	FindHull(S1, P, C, hullmyPoints);
	FindHull(S2, C, Q, hullmyPoints);
}

// http://www.cse.yorku.ca/~aaw/Hang/quick_hull/Algorithm.html
void QuickHull(const std::vector<myPoint> &s, std::vector<myPoint> &hullmyPoints) {
	// Find left and right most myPoints, say A & B, and add A & B to convex hull 
	myPoint A = s[0];
	myPoint B = s[s.size() - 1];
	hullmyPoints.push_back(A);
	hullmyPoints.push_back(B);

	std::vector<myPoint> S1;
	std::vector<myPoint> S2;

	/*
	* Segment AB divides the remaining (n-2) myPoints into 2 groups S1 and S2
	* where S1 are myPoints in S that are on the right side of the oriented line from A to B,
	* and S2 are myPoints in S that are on the right side of the oriented line from B to A
	*/
	for (auto it = s.begin() + 1; it != s.end() - 1; ++it) {
		const myPoint pt = *it;
		const int s1 = SideOfLine(A, B, pt);
		const int s2 = SideOfLine(B, A, pt);
		if (0 < s1) {
			S1.push_back(pt);
		}
		else if (0 < s2) {
			S2.push_back(pt);
		}
	}

	FindHull(S1, A, B, hullmyPoints);
	FindHull(S2, B, A, hullmyPoints);
}

Mat mypredone(IplImage* srcipl, vector<Mat> bgrimgs)
{
	////預處理
	Mat img(srcipl, 0);
	split(img, bgrimgs);
	Mat bimg(img.size(), CV_8UC1, cv::Scalar(0));
	bgrimgs[0].copyTo(bimg);
	threshold(bimg, bimg, 12, 255, CV_THRESH_BINARY/*|CV_THRESH_OTSU*/);
	Mat element(4, 4, CV_8U, cv::Scalar(1));
	morphologyEx(bimg, bimg, cv::MORPH_OPEN, element);
	for (int i = 0; i < 20; i++)
	{
		uchar* data = bimg.ptr<uchar>(i);
		for (int j = 0; j < bimg.size().width; j++)
			data[j] = 0;
	}
	//for (int i = bimg.size().height - 19; i < bimg.size().height; i++)
	//{
	//	uchar* data = bimg.ptr<uchar>(i);
	//	for (int j = 0; j < bimg.size().width; j++)
	//		data[j] = 0;
	//}
	return bimg;
}

IplImage* luobin_findcontour(Mat bimg)
{
	IplImage mybimgipl = bimg;
	
	CvMemStorage *storage = cvCreateMemStorage();
	CvSeq *contour = NULL, *hull = NULL;
	vector<CvPoint> allpoints;
	vector<myPoint> myallpoints;
	CvContourScanner scanner = cvStartFindContours(&mybimgipl, storage);
	while ((contour = cvFindNextContour(scanner)) != NULL)
	{
		hull = cvConvexHull2(contour, 0, CV_CLOCKWISE, 0);
		CvPoint pt0 = **(CvPoint**)cvGetSeqElem(hull, hull->total - 1);
		myPoint mypt0;
		mypt0.x = (float)pt0.x;
		mypt0.y = (float)pt0.y;
		////mypt0.operator<(mypt0);  
		allpoints.push_back(pt0);
		myallpoints.push_back(mypt0);
		for (int i = 0; i < hull->total; ++i)
		{
			CvPoint pt1 = **(CvPoint**)cvGetSeqElem(hull, i);
			allpoints.push_back(pt1);
			//cvLine(mydst, pt0, pt1, cv::Scalar(255));
			//pt0 = pt1;
			myPoint mypt1;
			mypt1.x = (float)pt1.x;
			mypt1.y = (float)pt1.y;
			////mypt1.operator<(mypt1);  //////////
			myallpoints.push_back(mypt1);
		}
	}
	std::sort(myallpoints.begin(), myallpoints.end());
	std::vector<myPoint> myhullPoints;
	QuickHull(myallpoints, myhullPoints);

	CvSeq *hull2 = NULL;
	CvSeq* ptseq = cvCreateSeq(CV_SEQ_KIND_GENERIC | CV_32SC2, sizeof(CvContour), sizeof(CvPoint), storage);
	for (int i = 0; i < myhullPoints.size(); ++i)
	{
		myPoint pt1 = myhullPoints[i];
		CvPoint myfinalpot1;
		myfinalpot1.x = pt1.x;
		myfinalpot1.y = pt1.y;
		cvSeqPush(ptseq, &myfinalpot1);
	}
	hull2 = cvConvexHull2(ptseq, 0, CV_CLOCKWISE, 1);
	IplImage *mydst3 = cvCreateImage(cvGetSize(&mybimgipl), 8, 1);
	cvZero(mydst3);
	cvDrawContours(mydst3, hull2, CV_RGB(255, 255, 255), CV_RGB(255, 255, 255), -1, CV_FILLED, 8);
	return mydst3;
}


char filename[100];
char newcontour[100];

void main()
{
	IplImage* src_ipl = cvLoadImage("12.bmp");
	cvShowImage("srcimage", src_ipl);
	vector<Mat> bgrimgs;
	Mat bimg = mypredone(src_ipl, bgrimgs);

	IplImage* luobin_result = luobin_findcontour(bimg);
	cvShowImage("luobin result2", luobin_result);
	waitKey(0);

	/*
		TickMeter tm;
		tm.start();
		for (int i = 55124; i <= 56460; i++)
		{
			int newi = i - 55123;
			sprintf(filename, "背面圖\\%d.bmp", i);
			sprintf(newcontour, "羅彬背面輪廓\\%d.bmp", newi);

			IplImage* src_ipl = cvLoadImage(filename);
			vector<Mat> bgrimgs;
			Mat bimg = mypredone(src_ipl, bgrimgs);

			IplImage* luobin_result = luobin_findcontour(bimg);
			cvSaveImage(newcontour, luobin_result);
		}

		tm.stop();
		cout << "count=" << tm.getCounter() << ",process time=" << tm.getTimeMilli() << endl;
	*/
}
和上面我自己的演算法是同一張原圖  的確比我的快了將近一倍!

但對帶噪聲的圖,直接用上面這個快速凸包尋找輪廓的辦法就不行:

我把原圖、灰度圖、上面這個快速凸包找輪廓的圖、我後一篇自己寫的結果圖  對比了一下:

char filename[100];
char newcontour[100];
char grayimg[100];
char myback[100];
char contrast[100];

void main()
{
	TickMeter tm;
	tm.start();

	CvSize mysize = cvSize(484 * 4, 364);
	for (int i = 1; i <= 1337; i++)
	{
		sprintf(filename, "背面圖\\%d.bmp", i);
		sprintf(newcontour, "羅彬背面輪廓\\%d.bmp", i);
		sprintf(grayimg, "背面灰度\\%d.bmp", i);
		sprintf(myback, "我的背面輪廓\\%d.bmp", i); 
		sprintf(contrast, "背面對比圖\\%d.bmp", i);
	
		IplImage* left = cvLoadImage(filename);
		IplImage* combine = cvCreateImage(mysize, left->depth, left->nChannels);
		cvZero(combine); 
		cvSetImageROI(combine, cvRect(0, 0, left->width, left->height)); 
		cvCopy(left, combine);
		cvResetImageROI(combine);
		IplImage* right = cvLoadImage(grayimg);
		cvSetImageROI(combine, cvRect(483, 0, right->width, right->height));
		cvCopy(right, combine);
		cvResetImageROI(combine);
		IplImage* right2 = cvLoadImage(newcontour);
		cvSetImageROI(combine, cvRect(484*2-1, 0, right2->width, right2->height));
		cvCopy(right2, combine);
		cvResetImageROI(combine);
		IplImage* right3 = cvLoadImage(myback);
		cvSetImageROI(combine, cvRect(484 * 3 - 1, 0, right3->width, right3->height));
		cvCopy(right3, combine);
		cvResetImageROI(combine);
		cvSaveImage(contrast, combine);
		

	}

	tm.stop();
	cout << "count=" << tm.getCounter() << ",process time=" << tm.getTimeMilli() << endl;
}

所以 對有雜質的還是不能夠直接用  還需要修改。
最終改成了:

#include <algorithm>
#include <cmath>
#include <vector>
#include <chrono>
#include <cstdio>
#include <random>
#include <vector>
#include<opencv2/opencv.hpp>
#include<opencv2/imgproc/imgproc.hpp>
#include<opencv2/highgui/highgui.hpp>
using namespace std;
using namespace cv;

struct myPoint {
	int x;
	int y;
	bool operator<(const myPoint& pt) {

		if (x == pt.x)
		{
			return (y < pt.y);
		}
		return (x < pt.x);
	}
};

// http://stackoverflow.com/questions/1560492/how-to-tell-whether-a-myPoint-is-to-the-right-or-left-side-of-a-line
int SideOfLine(const myPoint &P1, const myPoint &P2, const myPoint &P3) {
	return (P2.x - P1.x) * (P3.y - P1.y) - (P2.y - P1.y) * (P3.x - P1.x);
}

// https://en.wikipedia.org/wiki/Distance_from_a_myPoint_to_a_line#Line_defined_by_two_myPoints
float DistanceFromLine(const myPoint &P1, const myPoint &P2, const myPoint &P3) {
	return (std::abs((P2.y - P1.y) * P3.x - (P2.x - P1.x) * P3.y + P2.x * P1.y - P2.y * P1.x)
		/ std::sqrt((P2.y - P1.y) * (P2.y - P1.y) + (P2.x - P1.x) * (P2.x - P1.x)));
}

// http://stackoverflow.com/questions/13300904/determine-whether-myPoint-lies-inside-triangle
bool myPointInTriangle(const myPoint &p, const myPoint &p1, const myPoint &p2, const myPoint &p3) {
	float a = ((p2.y - p3.y) * (p.x - p3.x) + (p3.x - p2.x) * (p.y - p3.y)) / ((p2.y - p3.y) * (p1.x - p3.x) + (p3.x - p2.x) * (p1.y - p3.y));
	float b = ((p3.y - p1.y) * (p.x - p3.x) + (p1.x - p3.x) * (p.y - p3.y)) / ((p2.y - p3.y) * (p1.x - p3.x) + (p3.x - p2.x) * (p1.y - p3.y));
	float c = 1.0f - a - b;
	return (0.0f < a && 0.0f < b && 0.0f < c);
}

// http://www.cse.yorku.ca/~aaw/Hang/quick_hull/Algorithm.html
void FindHull(const std::vector<myPoint> &Sk, const myPoint P, const myPoint Q, std::vector<myPoint> &hullmyPoints) {

	if (Sk.size() == 0) return;

	std::vector<myPoint> S0;
	std::vector<myPoint> S1;
	std::vector<myPoint> S2;

	float furthestDistance = 0.0f;
	myPoint C;
	for (const auto &pt : Sk) {
		float distance = DistanceFromLine(P, Q, pt);
		if (distance > furthestDistance) {
			furthestDistance = distance;
			C = pt;
		}
	}

	hullmyPoints.push_back(C);

	for (const auto &pt : Sk) {
		if (myPointInTriangle(pt, P, C, Q)) {
			S0.push_back(pt);
		}
		else if (0 < SideOfLine(P, C, pt)) {
			S1.push_back(pt);
		}
		else if (0 < SideOfLine(C, Q, pt)) {
			S2.push_back(pt);
		}
	}

	FindHull(S1, P, C, hullmyPoints);
	FindHull(S2, C, Q, hullmyPoints);
}

// http://www.cse.yorku.ca/~aaw/Hang/quick_hull/Algorithm.html
void QuickHull(const std::vector<myPoint> &s, std::vector<myPoint> &hullmyPoints) {

	myPoint A = s[0];
	myPoint B = s[s.size() - 1];
	hullmyPoints.push_back(A);
	hullmyPoints.push_back(B);

	std::vector<myPoint> S1;
	std::vector<myPoint> S2;

	for (auto it = s.begin() + 1; it != s.end() - 1; ++it) {
		const myPoint pt = *it;
		const int s1 = SideOfLine(A, B, pt);
		const int s2 = SideOfLine(B, A, pt);
		if (0 < s1) {
			S1.push_back(pt);
		}
		else if (0 < s2) {
			S2.push_back(pt);
		}
	}

	FindHull(S1, A, B, hullmyPoints);
	FindHull(S2, B, A, hullmyPoints);
}

Mat mypredone(IplImage* srcipl, vector<Mat> bgrimgs)
{
	////預處理
	Mat img(srcipl, 0);
	split(img, bgrimgs);
	Mat bimg(img.size(), CV_8UC1, cv::Scalar(0));
	bgrimgs[0].copyTo(bimg);
	threshold(bimg, bimg, 10, 255, CV_THRESH_BINARY/*|CV_THRESH_OTSU*/);
	Mat element(4, 4, CV_8U, cv::Scalar(1));
	morphologyEx(bimg, bimg, cv::MORPH_OPEN, element);
	//for (int i = 0; i < 20; i++)
	//{
	//	uchar* data = bimg.ptr<uchar>(i);
	//	for (int j = 0; j < bimg.size().width; j++)
	//		data[j] = 0;
	//}
	//for (int i = bimg.size().height - 19; i < bimg.size().height; i++)
	//{
	//	uchar* data = bimg.ptr<uchar>(i);
	//	for (int j = 0; j < bimg.size().width; j++)
	//		data[j] = 0;
	//}
	return bimg;
}
//背面quickhull   去除溜槽走光點等噪聲
IplImage* luobin_backcontour(Mat bimg)
{
	IplImage mybimgipl = bimg;

	CvMemStorage *storage = cvCreateMemStorage();
	CvSeq *contour = NULL, *hull = NULL;
	vector<CvPoint> allpoints;
	vector<myPoint> myallpoints;
	CvContourScanner scanner = cvStartFindContours(&mybimgipl, storage);
	while ((contour = cvFindNextContour(scanner)) != NULL)
	{
		/////////////////////////////////去掉溜槽走光點的 想法
		CvMoments *moments = (CvMoments*)malloc(sizeof(CvMoments));
		cvMoments(contour, moments, 1);
		double moment10 = cvGetSpatialMoment(moments, 1, 0);
		double moment01 = cvGetSpatialMoment(moments, 0, 1);
		double area = cvGetSpatialMoment(moments, 0, 0);
		double posX = moment10 / area;
		double posY = moment01 / area;
		double mylenth = cvArcLength(contour);
		double myarea = cvContourArea(contour);
		//cout << "length=" << mylenth << "   area=" << myarea << " 重心: " << "(" << posX << " , " << posY << ")" << " " << endl;
		//center light
		if (posX > double(236.0) && posX<double(243.0) && posY>double(171.0) && posY<double(197.0) && myarea<double(141) && mylenth < double(71))
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		//right light
		if (posX > double(384.0) && posX<double(395.0) && posY>double(90.0) && posY<double(113.0) && myarea<double(316) && mylenth < double(90))
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		//below
		if (posX >160.0 && posX<331.0 && posY>295.0 && posY < 362.0 && myarea < 170 && mylenth < 100)
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		if (posX >212.0 && posX<251.0 && posY>248.0 && posY < 269.0 && myarea < 21 && mylenth < 20)
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		if (posX >256.0 && posX<271.0 && posY>296.0 && posY < 299.5 && myarea < 25 && mylenth < 20)
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		if (posX >210.0 && posX<215.0 && posY>350.0 && posY < 358.0 && myarea < 250 && mylenth < 100)
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		// special 
		if (myarea == 0)
		{
			cvSubstituteContour(scanner, NULL);
			continue;
		}
		///////////////////////////////////////////////////////////
		hull = cvConvexHull2(contour, 0, CV_CLOCKWISE, 0);
		CvPoint pt0 = **(CvPoint**)cvGetSeqElem(hull, hull->total - 1);
		myPoint mypt0;
		mypt0.x = (float)pt0.x;
		mypt0.y = (float)pt0.y;
		allpoints.push_back(pt0);
		myallpoints.push_back(mypt0);
		for (int i = 0; i < hull->total; ++i)
		{
			CvPoint pt1 = **(CvPoint**)cvGetSeqElem(hull, i);
			allpoints.push_back(pt1);
			myPoint mypt1;
			mypt1.x = (float)pt1.x;
			mypt1.y = (float)pt1.y;
			myallpoints.push_back(mypt1);
		}
	}
	std::sort(myallpoints.begin(), myallpoints.end());
	std::vector<myPoint> myhullPoints;
	QuickHull(myallpoints, myhullPoints);

	CvSeq *hull2 = NULL;
	CvSeq* ptseq = cvCreateSeq(CV_SEQ_KIND_GENERIC | CV_32SC2, sizeof(CvContour), sizeof(CvPoint), storage);
	for (int i = 0; i < myhullPoints.size(); ++i)
	{
		myPoint pt1 = myhullPoints[i];
		CvPoint myfinalpot1;
		myfinalpot1.x = pt1.x;
		myfinalpot1.y = pt1.y;
		cvSeqPush(ptseq, &myfinalpot1);
	}
	hull2 = cvConvexHull2(ptseq, 0, CV_CLOCKWISE, 1);
	IplImage *mydst3 = cvCreateImage(cvGetSize(&mybimgipl), 8, 1);
	cvZero(mydst3);
	cvDrawContours(mydst3, hull2, CV_RGB(255, 255, 255), CV_RGB(255, 255, 255), -1, CV_FILLED, 8);
	return mydst3;
}

//正面quickhull
IplImage* luobin_frontcontour(Mat bimg)
{
	IplImage mybimgipl = bimg;

	CvMemStorage *storage = cvCreateMemStorage();
	CvSeq *contour = NULL, *hull = NULL;
	vector<CvPoint> allpoints;
	vector<myPoint> myallpoints;
	CvContourScanner scanner = cvStartFindContours(&mybimgipl, storage);
	while ((contour = cvFindNextContour(scanner)) != NULL)
	{
		hull = cvConvexHull2(contour, 0, CV_CLOCKWISE, 0);
		CvPoint pt0 = **(CvPoint**)cvGetSeqElem(hull, hull->total - 1);
		myPoint mypt0;
		mypt0.x = (float)pt0.x;
		mypt0.y = (float)pt0.y;
		allpoints.push_back(pt0);
		myallpoints.push_back(mypt0);
		for (int i = 0; i < hull->total; ++i)
		{
			CvPoint pt1 = **(CvPoint**)cvGetSeqElem(hull, i);
			allpoints.push_back(pt1);
			myPoint mypt1;
			mypt1.x = (float)pt1.x;
			mypt1.y = (float)pt1.y;
			myallpoints.push_back(mypt1);
		}
	}
	std::sort(myallpoints.begin(), myallpoints.end());
	std::vector<myPoint> myhullPoints;
	QuickHull(myallpoints, myhullPoints);

	CvSeq *hull2 = NULL;
	CvSeq* ptseq = cvCreateSeq(CV_SEQ_KIND_GENERIC | CV_32SC2, sizeof(CvContour), sizeof(CvPoint), storage);
	for (int i = 0; i < myhullPoints.size(); ++i)
	{
		myPoint pt1 = myhullPoints[i];
		CvPoint myfinalpot1;
		myfinalpot1.x = pt1.x;
		myfinalpot1.y = pt1.y;
		cvSeqPush(ptseq, &myfinalpot1);
	}
	hull2 = cvConvexHull2(ptseq, 0, CV_CLOCKWISE, 1);
	IplImage *mydst3 = cvCreateImage(cvGetSize(&mybimgipl), 8, 1);
	cvZero(mydst3);
	cvDrawContours(mydst3, hull2, CV_RGB(255, 255, 255), CV_RGB(255, 255, 255), -1, CV_FILLED, 8);
	return mydst3;
}

char filename[100];
char newcontour[100];
char grayimg[100];
char myback[100];
char contrast[100];

void main()
{
	TickMeter tm;
	tm.start();

	for (int i = 1; i <= 1337; i++)
	{
		  //////////////背面找輪廓
		sprintf(filename, "背面圖\\%d.bmp", i);
		sprintf(newcontour, "羅彬背面輪廓3\\%d.bmp", i);
		IplImage* src_ipl = cvLoadImage(filename);
		vector<Mat> bgrimgs;
		Mat bimg = mypredone(src_ipl, bgrimgs);
		IplImage* luobin_result = luobin_backcontour(bimg);
		cvSaveImage(newcontour, luobin_result);
		

		    ////////////////////正面找輪廓
		sprintf(filename, "正面圖\\%d.bmp", i);
		sprintf(newcontour, "羅彬正面輪廓\\%d.bmp", i);
		IplImage* src_ipl = cvLoadImage(filename);
		vector<Mat> bgrimgs;
		Mat bimg = mypredone(src_ipl, bgrimgs);
		IplImage* luobin_result = luobin_frontcontour(bimg);
		
	}

	tm.stop();
	cout << "count=" << tm.getCounter() << ",process time=" << tm.getTimeMilli() << endl;
}
這樣就行了!