一功能要求:
1.给定一堆三维点,能自动生成其最小凸多面体。
2.随机生成一些空间点,可以判断随机点是否在凸多面体内。
二代码实现:
#define _CRT_SECURE_NO_WARNINGS
#include <random>
#include <ctime>
#include <iostream>
#include <vector>
#include<algorithm>
#include <math.h>
#include <glut.h>
#define M 600
#define eps 1e-10
using namespace std;
//定义点类
class node {
public:
double x, y, z, dis; //三维坐标
node() {};
node(double xx, double yy, double zz) :x(xx), y(yy), z(zz) {};
//运算符重载
node operator + (const node n) { //向量+加法重载
node temp;
temp.x = this->x + n.x;
temp.y = this->y + n.y;
temp.z = this->z + n.z;
return temp;
}
node operator -(const node n) { //向量-减法重载
node temp;
temp.x = this->x - n.x;
temp.y = this->y - n.y;
temp.z = this->z - n.z;
return temp;
}
node operator *(const double m) { //向量数乘
return node(m*this->x, m*this->y, m*this->z);
}
node operator /(const double m) { //向量数除
return node(this->x / m, this->y / m, this->z / m);
}
node operator *(const node n) {//向量叉乘
return node(n.y *this->z - this->y*n.z, n.z*this->x - n.x*this->z, n.x*this->y - n.y*this->x);
}
double operator ^(const node n) { //向量点乘
return (this->x*n.x) + this->y*n.y + this->z*n.z;
}
};
//定义面类
class face {
public:
int a, b, c; //面上三个顶点编号
int ok; //标志位,判断该面是否是凸包中的面
};
//定义全局变量
node p[M]; //对象数组定义初始点
node pj[M]; //对象数组 定义需要判断的点
node pji[M];
node pjo[M];
face f[M * 8]; //对象数组 节点构成的凸包三角形面
vector<int> v;
class Polygin { //定义三维凸包类
public:
int n; //初始点数
int m; //判断点数
int cnt; //凸包三角形数
int to[M][M]; // 边ij属于哪个面
//功能函数
double len(node p) { //求向量长度
return sqrt(p.x *p.x + p.y*p.y + p.z*p.z);
}
double area(node a, node b, node c) { //求三个点围成三角形的面积
double pow;
pow = 0.5*len(b - a)*len(c - b)*len(a - c);
return sqrt(pow*(pow - len(b - a))*(pow - len(c - b))*(pow - len(a - c)));
}
double volume(node a, node b, node c, node d) {
return (b - a)*(c - a) ^ (d - a); //求四个点围成的体积
}
double ptof(node q, face f) { //判断面外一点与面的内外关系,可见性
node m = p[f.b] - p[f.a]; //面内两个向量
node n = p[f.c] - p[f.a];
node t = q - p[f.a]; //该点与面内一点的向量
double p;
p = m * n^t; //向叉乘再点乘 通过判断pd的正负判断
return p;
}
//凸包更新函数
void deal(int q, int a, int b) { //a,b 表示边的顶点编号 q表示:
int fa = to[a][b]; //与当前面cnt共边的另一个面
face add;
if (f[fa].ok) //f为对象数组,f[fa]是一个面对象,f[fa].ok是该面是否是凸包表面的标志位
{
if (ptof(p[q], f[fa]) > eps)//q能看到f[fa]面
delet(q, fa); //删除该面
else // q点能看到cnt 面,但是看不到f[fa]面,此时q与a,b组成一个新的面add
{
add.a = b; //将b,a,q三点在add面上重新分配
add.b = a;
add.c = q;
add.ok = 1;
to[b][a] = to[a][q] = to[q][b] = cnt; //q,a,b三个点构成了新的面
f[cnt++] = add; //当前面存入对象数组中
}
}
}
//删除面的函数
void delet(int q, int cur) //cur为当前面 q表示面外该点
{
f[cur].ok = 0;
/*删除该面后,还需判断与该面相邻的三个面与q点之间关系,需要三次调用deal函数*/
/*点是以逆时针顺序存放的,相邻的面中点的顺序与当前面刚好相反*/
deal(q, f[cur].b, f[cur].a);
deal(q, f[cur].c, f[cur].b);
deal(q, f[cur].a, f[cur].c); //三条边,判断三次
}
//判断两个三角形是否共面
int same(int s, int t)
{
node a = p[f[s].a];
node b = p[f[s].b];
node c = p[f[s].c];
if (fabs(volume(a, b, c, p[f[t].a])) < eps //通过四点体积判断
&&fabs(volume(a, b, c, p[f[t].b])) < eps
&&fabs(volume(a, b, c, p[f[t].c])) < eps)
return 1;
return 0;
}
void make() //凸包构建函数
{
//part-1
cnt = 0;
if (n < 4)return;
int sb = 1; //立flag
for (int i = 1; i < n; i++) //前两个点不共点
{
if (len(p[0] - p[i])) //若0点在1点前
{
swap(p[1], p[i]); // 交换两个点
sb = 0;
break;
}
}
if (sb)return;
sb = 1;
for (int i = 2; i < n; i++) //前三个点不共线
{
if (len((p[1] - p[0]) * (p[i] - p[0])) > eps)
{
swap(p[2], p[i]);
sb = 0;
break;
}
}
if (sb)return;
sb = 1;
for (int i = 3; i < n; i++) //前四个点不共面
{
if (fabs(volume(p[0], p[1], p[2], p[i])) > eps)
{
swap(p[3], p[i]);
sb = 0;
break;
}
}
if (sb) return;
//part-2
face add;
for (int i = 0; i < 4; i++)//构建初始四面体(4个点为p[0],p[1],[2],p[3])
{
add.a = (i + 1) % 4;
add.b = (i + 2) % 4;
add.c = (i + 3) % 4;
add.ok = 1;
if (ptof(p[i], add) > eps) swap(add.c, add.b); //保证逆时针,即法向量朝外,这样新点才可看到
to[add.a][add.b] = to[add.b][add.c] = to[add.c][add.a] = cnt; //按照逆时针方向保存边到某个面中 //每条边上存有面的编号
f[cnt++] = add; //
}
//part-3
for (int i = 4; i < n; i++) //判断每个凸包外的点与凸包每个面的关系
{
for (int j = 0; j < cnt; j++) //逐个面判断
{
if (f[j].ok&&ptof(p[i], f[j]) > eps) //判断点能否看到当前面
{
delet(i, j); //删除当前面 并进行凸包更新
break; //当前面更新完成,跳出
}
}
}
//part-4
int tmp = cnt;
cnt = 0;
for (int i = 0; i < tmp; i++)
if (f[i].ok)
{
f[cnt++] = f[i];
//将每个面上的三个点进行储存
v.push_back(f[i].a);
v.push_back(f[i].b);
v.push_back(f[i].c);
}
//去掉相同的点并按顺序排列
sort(v.begin(), v.end());
v.erase(unique(v.begin(), v.end()), v.end());
//输出
cout << "*********************************" << endl;
cout << "凸多面体顶点为:" << endl;
for (int i = 0; i < v.size(); i++)
{
cout << v[i] << ":";
cout << "(" << p[v[i]].x << "," << p[v[i]].y << "," << p[v[i]].z << ")" << endl;
}
}
}; Polygin hull;
//判断点在凸包内外的函数
int mmm = 0, nnn = 0;
void judge(node p)
{
cout <<'('<< p.x <<','<< p.y <<','<< p.z<<')';
int i;
for ( i = 0; i < hull.cnt; i++)
{
if (hull.ptof(p, f[i]) > eps) //调用对象函数ptof
{
pjo[mmm].x = p.x;
pjo[mmm].y = p.y;
pjo[mmm].z = p.z;
mmm++; //外部点计数器
cout << "该点在凸多面体外部" << endl;
break;
}
}
if (i == hull.cnt)
{
pji[nnn].x = p.x;
pji[nnn].y = p.y;
pji[nnn].z = p.z;
nnn++; //内部点计数器
cout << "该点在凸多面体内部" << endl;
}
}
//坐标轴设置
void axis(double length) //建立坐标系
{
glColor3f(1.0f, 1.0f, 1.0f); //使用颜色设置白色
glPushMatrix(); //堆栈变换,压入堆栈
glBegin(GL_LINES); //指定线单元
glVertex3d(0.0, 0.0, 0.0); //指定点(0,0,0)
glVertex3d(0.0, 0.0, length); //指定点
glEnd();//将当前操作点移到指定位置
glTranslated(0.0, 0.0, length - 0.2); //沿着Z轴平移相应的距离
glColor3f(1.0, 1.0, 1.0); //颜色设置 白色
glutSolidCone(0.04, 0.3, 8, 8); //生成坐标轴箭头
glPopMatrix(); //弹出堆栈
}
void paint(void)
{
//Part1
glClear(GL_COLOR_BUFFER_BIT); //利用之前设置好的颜色作背景颜色
glMatrixMode(GL_PROJECTION); //对投影矩阵应用随后的矩阵操作
glLoadIdentity(); //将当前的用户坐标系的原点移到了屏幕中心:类似于一个复位操作
glOrtho(-2.0, 2.0, -2.0, 2.0, -100, 100); //将当前的可视空间设置为正投影空间
glPointSize(1); //指定栅格化点的直径
glMatrixMode(GL_MODELVIEW); //对模型视景矩阵堆栈应用随后的矩阵操作
glLoadIdentity();
gluLookAt(1.2, 1.4, 2.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0); //坐标系的大体方向1.3, 1.6, 2.0
//Part2
//画坐标系
axis(2); //z轴
glPushMatrix();
glRotated(90.0, 0, 1.0, 0); //绕y轴正方向旋转90度 //x轴
axis(2);
glPopMatrix();
glPushMatrix();
glRotated(-90.0, 1.0, 0.0, 0.0); //绕x轴负方向旋转 //y轴
axis(2);
glPopMatrix();
//Part3
//画P点
for (int i = 0; i < M; i++)
{
glColor3f(0.0, 1.0, 1.0); //画笔蓝色
glPointSize(4);
glBegin(GL_POINTS);
glVertex3f(p[i].x, p[i].y, p[i].z);
glEnd();
}
//画凸多面体顶点
for (int i = 0; i < v.size(); i++)
{
glColor3f(1.0, 0.0, 0.0); //画笔红色
glPointSize(6);
glBegin(GL_POINTS);
glVertex3f(p[v[i]].x, p[v[i]].y, p[v[i]].z);
glEnd();
}
//画凸多面体面
for (int i = 0; i < hull.cnt; i++)
{
glColor3f(0.0, 0.0, 1.0); //画笔蓝色
glPointSize(6);
glBegin(GL_LINE_LOOP);
glVertex3f(p[f[i].a].x, p[f[i].a].y, p[f[i].a].z);
glVertex3f(p[f[i].b].x, p[f[i].b].y, p[f[i].b].z);
glVertex3f(p[f[i].c].x, p[f[i].c].y, p[f[i].c].z);
glEnd();
}
//画所需判断Pm点
for (int i = 0; i < mmm; i++)
{
glColor3f(1, 1, 1); //外部点画笔白色
glPointSize(6);
glBegin(GL_POINTS);
glVertex3f(pjo[i].x, pjo[i].y, pjo[i].z);
glEnd();
}
for (int i = 0; i < nnn; i++)
{
glColor3f(0, 1, 0); //内部点画笔绿色
glPointSize(6);
glBegin(GL_POINTS);
glVertex3f(pji[i].x, pji[i].y, pji[i].z);
glEnd();
}
glPopMatrix();//用来整体绕y轴旋转
glutSwapBuffers();
}
void Init()
{
glClearColor(0.5, 0.5, 0.5, 1.0); //设置渲染背景颜色
}
int main(int argv, char *argc[])
{
cout << "输入三维点的数量:" ;
scanf("%d", &hull.n);
while (hull.n!= EOF) //输入点的个数
{
default_random_engine e(time(0)); //引擎,随机数生成引擎
uniform_real_distribution<double> u(0, 2.0); //分布,u内为生成随机数范围
for (int i = 0; i < hull.n; i++) {
p[i].x = static_cast<int>((1000 * u(e) + 5) / 10) / 100.0;
p[i].y = static_cast<int>((1000 * u(e) + 5) / 10) / 100.0;
p[i].z = static_cast<int>((1000 * u(e) + 5) / 10) / 100.0;
cout<<"输入第"<<i<<"个点坐标为: "<< "(" << p[i].x << "," << p[i].y << "," << p[i].z << ")" << endl;
}
hull.make();
cout << endl<<endl<< endl;
cout << "输入待判断的三维点个数:";
scanf("%d", &hull.m);
default_random_engine r(time(0));
uniform_real_distribution<double> n(0.1, 2);
for (int i = 0; i < hull.m; i++)
{
pj[i].x = static_cast<int>((1000 * n(e) + 5) / 10) / 100.0;
pj[i].y = static_cast<int>((1000 * n(e) + 5) / 10) / 100.0;
pj[i].z = static_cast<int>((1000 * n(e) + 5) / 10) / 100.0;
judge(pj[i]);
}
cout << "*********************************" << endl;
printf("三维凸包显示");
glutInit(&argv, argc);
glutInitDisplayMode(GLUT_RGB | GLUT_DOUBLE); //设定显示模式
glutInitWindowSize(800, 600); //设置初始窗口的大小
glutInitWindowPosition(0, 0); //设置初始窗口的位置(相对屏幕左上角)
glutCreateWindow("最小凸多面体的生成"); //生成一个分辨率为600*600名称为“最小凸多面体的窗口”
Init(); //OpenGL初始化代码
glutDisplayFunc(paint); //注册回调函数paint
glutMainLoop();
}
return 0;
}
三.运行测试