一个有趣的数学现象

转自：http://www.bfcat.com/index.php/2012/11/random-polygon-ellipse/

今天在一个博客上看到的。

在一张纸上随机画一些点，然后将他们连起来形成随机的多边形。找到多边形所有边的中点并将这些中点连成一个新的多边形，如此重复下去，最终的结果会是一个越来越小的形状，而且会趋向与一个椭圆！

这个过程如下图，也可以进入作者提供的一个app网页自己动手操作一下：

matthen 博客里还给了这个过程的代码。这个博客还介绍了其他很多有意思的数学和科学问题，值得订阅。

==================     以下内容为原创 ==================

本文地址：http://blog.csdn.net/houston11235/article/details/8498938转载请注明

这个现象看起来挺有意思的，于是用 OpenCV 模拟了一下，遇到了几个问题

1、坐标归一化，由于这么取中点多边形肯定是越来越小的，不做归一化一会就收到到一点了

2、坐标的数据类型，开始用的整型，发现点多了很难收敛，可能是在近似的过程中达到稳态了，于是改成浮点型，这下收敛了，点数多也能收敛到椭圆

  

  

把代码贴出来，有想一探究竟的朋友可以拿去看看

// Poly2Ellipse.h

#include <opencv.hpp>
#include <vector>
#include <iostream>

using namespace std;
using namespace cv;

#define img_width 400
#define img_height 400
#define n_pts 100

// Poly2Ellipse.cpp

#include "Poly2Ellipse.h"

void display(Mat& img, string win_name, int wait_time = 1)
{
	namedWindow(win_name);
	imshow(win_name, img);
	waitKey(wait_time);
}

void drawPts(Mat& img, vector<Point2f>& pts)
{
	for (size_t i = 0; i < pts.size(); ++i)
	{
		circle(img, pts[i], 3, Scalar(0,255,0));
	}

	for (size_t i = 0; i < pts.size() - 1; ++i)
	{
		line(img, pts[i], pts[i+1], Scalar(0,0,255));
	}
	line(img, pts.back(), pts.front(), Scalar(0,0,255));
}

void interpPt(const vector<Point2f>& pts1, vector<Point2f>& pts2)
{
	float x_max = -1;
	float y_max = -1;
	float x_min = 10000;
	float y_min = 10000;
	for (size_t i = 0; i < pts1.size() - 1; ++i)
	{
		pts2[i].x = (pts1[i].x + pts1[i + 1].x) / 2;
		pts2[i].y = (pts1[i].y + pts1[i + 1].y) / 2;
		if (pts2[i].x > x_max)
		{
			x_max = pts2[i].x;
		}
		if (pts2[i].y > y_max)
		{
			y_max = pts2[i].y;
		}
		if (pts2[i].x < x_min)
		{
			x_min = pts2[i].x;
		}
		if (pts2[i].y < y_min)
		{
			y_min = pts2[i].y;
		}
	}
	pts2.back().x = (pts1.front().x + pts1.back().x) / 2;
	pts2.back().y = (pts1.front().y + pts1.back().y) / 2;
	if (pts2.back().x > x_max)
	{
		x_max = pts2.back().x;
	}
	if (pts2.back().y > y_max)
	{
		y_max = pts2.back().y;
	}
	if (pts2.back().x < x_min)
	{
		x_min = pts2.back().x;
	}
	if (pts2.back().y < y_min)
	{
		y_min = pts2.back().y;
	}
	float x_ratio = 0.8 * img_width / (x_max - x_min);
	float y_ratio = 0.8 * img_height / (y_max - y_min);
	int x_shift = 0.1 * img_width;
	int y_shift = 0.1 * img_height;

	for (size_t i = 0; i < pts2.size(); ++i)
	{
		pts2[i].x = x_ratio * (pts2[i].x - x_min) + x_shift;
		pts2[i].y = y_ratio * (pts2[i].y - y_min) + y_shift;
	}
}

int main(void)
{
	Mat img_show(img_height, img_width, CV_8UC3, Scalar(0));
	Mat img_black = img_show.clone();

	//int n_pts = 10;
	RNG rng(getTickCount());
	Mat tmp(1, 2, CV_32S);
	vector<Point2f> pts1(n_pts);
	vector<Point2f> pts2(n_pts);
	for (int i = 0; i < n_pts; ++i)
	{
		pts1[i].x = rng.uniform(0, img_width);
		pts1[i].y = rng.uniform(0, img_height);
	}
	drawPts(img_show, pts1);
	display(img_show, "img");

	for (;;)
	{
		interpPt(pts1, pts2);
		img_black.copyTo(img_show);
		pts1 = pts2;
		drawPts(img_show, pts1);
		display(img_show, "interp");
		if (waitKey(50)=='e')
		{
			break;
		}
	}

	waitKey(0);
	return 0;
}