锥盒碰撞 [英] Cone to box collision

查看:12
本文介绍了锥盒碰撞的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我希望在一个圆锥体(底部是圆的.所以它基本上是一个球体)和一个盒子之间实现碰撞检测.我不太在意它是 AABB 还是 OBB,因为转换应该足够简单.我找到的每个解决方案都使用三角锥,但我的锥更像是一个有角度和距离的弧".

是否有一个简单的解决方案来进行这种碰撞检测?或者是做几种测试的情况?IE.比如在一个球体上获取交点,r 是我的圆锥距离,然后测试它们是否在一个角度内相交?

解决方案

我很好奇并计划用 GLSL 数学风格做一些需要的东西.所以这里有不同的方法.让我们考虑一下这个圆锥的定义:

  1. 创建一组基本几何图元

你需要支持点、线、三角形、凸三角网格、球面(圆锥).

  1. 在点和三角形、网格、圆锥之间实现内部测试

对于 triangle 任意边和点之间的交叉结果 - 边原点应指向三角形的同一侧(如正常).如果不是点在外面.

对于 convex mesh 点-面原点和面法线指向之间的点积对于所有面应为 <=0.

对于cone,该点应在球体半径内,并且锥轴与点锥原点之间的角度应为<= ang.再次点积可以用于此.

  1. 在基本基元之间实现最近线

这就像在形成一条线的每个图元上找到最近的点.它类似于垂直距离.

point-point很简单,因为它们是最近的线.

point-line 可以使用将点投影到线上(点积)来完成.但是,您需要限制结果,使其位于行内,而不是外推.

point-triangle 可以作为所有圆周线与点组合的最小值以及到表面的垂直距离(与三角形法线的点积)获得.

所有其他原语组合都可以从这些基本组合构建.

  1. 网格和圆锥之间的最近线

只需在圆锥球体中心和网格之间使用最近的线.如果线位于圆锥内,则将其缩短球体半径 R.这将考虑所有帽的相互作用.

然后在锥体表面上测试线,沿其圆周采样,从锥体球心开始,到最外圈(锥体和盖子之间的边缘)结束.如果您需要更高的精度,您也可以测试三角形.

  1. 网格和圆锥的交点

这很容易,只需计算网格和圆锥之间的最近留置权.然后测试它在网格侧的点是否在圆锥内.

检查

 `bool intersect(convex_mesh m0,spherical_sector s0);`

在下面的代码中实现.

这里是小的 C++/OpenGL 示例(使用 GLSL 风格的数学):

//--------------------------------------------------------------------------//--- GL几何 ------------------------------------------------------------------------//--------------------------------------------------------------------------#ifndef _gl_geometry_h#define _gl_geometry_h//--------------------------------------------------------------------------常量浮点度=M_PI/180.0;常量浮动弧度=180.0/M_PI;浮动除法(float a,float b){ if (fabs(b)<1e-10) return 0.0;否则返回 a/b;}双除(双a,双b){如果(fabs(b)<1e-10)返回0.0;否则返回 a/b;}#include GLSL_math.h"#include List.h"//--------------------------------------------------------------------------上课点{上市://配置文件vec3 p0;观点()     {}点(点& a){ *this=a;}〜点(){}点*运算符 = (const point *a) { *this=*a;返回这个;}//point* operator = (const point &a) { ...copy... return this;}点(vec3 _p0){p0=_p0;计算();}无效计算(){};无效绘制(){glBegin(GL_POINTS);glVertex3fv(p0.dat);glEnd();}};//--------------------------------------------------------------------------类轴{上市://配置文件vec3 p0,dp;轴(){}轴(轴& a){ *this=a;}〜轴(){}轴* 运算符 = (常量轴 *a) { *this=*a;返回这个;}//axis* operator = (const axis &a) { ...copy... return this;}轴(vec3 _p0,vec3 _dp){p0=_p0;dp=_dp;计算();}无效计算(){dp=标准化(dp);}无效绘制(){vec3 p;p=p0+100.0*dp;glBegin(GL_LINES);glVertex3fv(p0.dat);glVertex3fv(p .dat);glEnd();}};//--------------------------------------------------------------------------班线{上市://配置文件vec3 p0,p1;//计算浮动 l;vec3 dp;线()  {}line(line&a) { *this=a;}〜线(){}line* operator = (const line *a) { *this=*a;返回这个;}//line* operator = (const line &a) { ...copy... return this;}线(vec3 _p0,vec3 _p1){p0=_p0;p1=_p1;计算();}无效交换(){vec3 p=p0;p0=p1;p1=p;}无效计算(){dp=p1-p0;l=长度(dp);}无效绘制(){glBegin(GL_LINES);glVertex3fv(p0.dat);glVertex3fv(p1.dat);glEnd();}};//--------------------------------------------------------------------------类三角形{上市://配置文件vec3 p0,p1,p2;//计算vec3 n;三角形()  {}三角形(三角形& a) { *this=a;}〜三角形(){}三角形* 运算符 = (const triangle *a) { *this=*a;返回这个;}//triangle* operator = (const triangle &a) { ...copy... return this;}三角形(vec3 _p0,vec3 _p1,vec3 _p2){p0=_p0;p1=_p1;p2=_p2;计算();}无效交换(){vec3 p=p1;p1=p2;p2=p;n=-n;}无效计算(){n=标准化(交叉(p1-p0,p2-p1));}无效绘制(){glBegin(GL_TRIANGLES);glNormal3fv(n.dat);glVertex3fv(p0.dat);glVertex3fv(p1.dat);glVertex3fv(p2.dat);glEnd();}};//--------------------------------------------------------------------------类凸网格{上市://配置文件列表<三角形>三;//计算vec3 p0;//中央凸网格(){ tri.num = 0;}凸网格(凸网格& a){*this=a;}~convex_mesh() {}凸网格* 运算符 = (常量凸网格 *a) { *this=*a;返回这个;}//convex_mesh* operator = (const convex_mesh &a) { ...copy... return this;}void init_box(vec3 _p0,vec3 _u,vec3 _v,vec3 _w)//中心,一半大小{常量 vec3 p[8]={_p0-_u+_v-_w,_p0+_u+_v-_w,_p0+_u-_v-_w,_p0-_u-_v-_w,_p0-_u-_v+_w,_p0+_u-_v+_w,_p0+_u+_v+_w,_p0-_u+_v+_w,};常量 int ix[36]={0,1,2,0,2,3,4,5,6,4,6,7,3,2,5,3,5,4,2,1,6,2,6,5,1,0,7,1,7,6,0,3,4,0,4,7,};tri.num=0;for (int i=0;i<36;i+=3) tri.add(triangle(p[ix[i+0]],p[ix[i+1]],p[ix[i+2]]));计算();}无效计算(){诠释我,n;p0=vec3(0.0,0.0,0.0);如果(!tri.num)返回;对于 (i=0,n=0;i0.75) u=vec3(0.0,1.0,0.0);v=交叉(u,w);u=交叉(v,w);u=标准化(u);v=标准化(v);//计算表da=2.0*M_PI/float(N-1);db=ang/float(M-1);对于 (a=0.0,i=0;i1e-6) 返回 false;浮动 d0,d1,d2;d0=点(t0.n,交叉(p0.p0-t0.p0,t0.p1-t0.p0));d1=点(t0.n,交叉(p0.p0-t0.p1,t0.p2-t0.p1));d2=点(t0.n,交叉(p0.p0-t0.p2,t0.p0-t0.p2));如果 (d0*d1<-1e-6) 返回 false;如果 (d0*d2<-1e-6) 返回假;如果 (d1*d2<-1e-6) 返回假;返回真;}bool inside(point p0,convex_mesh m0){for (int i=0;i0.0)返回假;返回真;}bool inside(point p0,spherical_sector s0){浮动 t,l;vec3你;u=p0.p0-s0.p0;l=长度(u);如果(l>s0.R)返回假;t=除法(点(u,s0.dp),(l*s0.R));if (t1e-6)&&(inside(p,t0))){ ll=line(p0.p0,p.p0);如果 (cl.l>ll.l) cl=ll;}ll=最近的(p0,线(t0.p0,t0.p1));如果 (cl.l>ll.l) cl=ll;ll=最近的(p0,线(t0.p1,t0.p2));如果 (cl.l>ll.l) cl=ll;ll=最近的(p0,线(t0.p2,t0.p0));如果 (cl.l>ll.l) cl=ll;返回 cl;}最接近的线(点 p0,convex_mesh m0){诠释我;线 cl,ll;cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0));cl.l=1e300;对于 (i=0;ill.l) cl=ll;}返回 cl;}//--------------------------------------------------------------------------最接近的线(轴 a0,点 p0){ 线 cl;cl=最近的(p0,a0);cl.swap();返回 cl;}最近的线(轴 a0,轴 a1){vec3 u=a0.dp;vec3 v=a1.dp;vec3 w=a0.p0-a1.p0;浮动一个=点(u,u);//总是 >= 0浮动 b=dot(u,v);浮动 c=点(v,v);//总是 >= 0浮动 d=dot(u,w);浮动 e=dot(v,w);浮动 D=a*c-b*b;//总是 >= 0浮动t0,t1;//计算两个最近点的线参数if (D<1e-6)//线几乎平行{t0=0.0;t1=(b>c?d/b:e/c);//使用最大的分母}别的{t0=(b*e-c*d)/D;t1=(a*e-b*d)/D;}返回线(a0.p0+(a0.dp*t0),a1.p0+(a1.dp*t1));}最近的线(轴 a0,线 l1){vec3 u=a0.dp;vec3 v=l1.dp;vec3 w=a0.p0-l1.p0;浮动一个=点(u,u);//总是 >= 0浮动 b=dot(u,v);浮动 c=点(v,v);//总是 >= 0浮动 d=dot(u,w);浮动 e=dot(v,w);浮动 D=a*c-b*b;//总是 >= 0浮动t0,t1;//计算两个最近点的线参数if (D<1e-6)//线几乎平行{t0=0.0;t1=(b>c?d/b:e/c);//使用最大的分母}别的{t0=(b*e-c*d)/D;t1=(a*e-b*d)/D;}如果(t1<0.0)t1=0.0;如果 (t1>1.0) t1=1.0;返回线(a0.p0+(a0.dp*t0),l1.p0+(l1.dp*t1));}最近的线(轴 a0,三角形 t0){线 cl,ll;cl=最近的(a0,线(t0.p0,t0.p1));ll=最近的(a0,线(t0.p1,t0.p2));如果 (cl.l>ll.l) cl=ll;ll=最近的(a0,线(t0.p2,t0.p0));如果 (cl.l>ll.l) cl=ll;返回 cl;}最近的线(轴 a0,凸网格 m0){诠释我;线 cl,ll;cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0));cl.l=1e300;对于 (i=0;ill.l) cl=ll;}返回 cl;}//--------------------------------------------------------------------------线最近(线 l0,点 p0){ 线 cl;cl=最近的(p0,l0);cl.swap();返回 cl;}线最近(线 l0,轴 a0){ 线 cl;cl=最近的(a0,l0);cl.swap();返回 cl;}最近的线路(线路 l0,线路 l1){vec3 u=l0.p1-l0.p0;vec3 v=l1.p1-l1.p0;vec3 w=l0.p0-l1.p0;浮动一个=点(u,u);//总是 >= 0浮动 b=dot(u,v);浮动 c=点(v,v);//总是 >= 0浮动 d=dot(u,w);浮动 e=dot(v,w);浮动 D=a*c-b*b;//总是 >= 0浮动t0,t1;//计算两个最近点的线参数if (D<1e-6)//线几乎平行{t0=0.0;t1=(b>c?d/b:e/c);//使用最大的分母}别的{t0=(b*e-c*d)/D;t1=(a*e-b*d)/D;}如果(t0<0.0)t0=0.0;如果 (t0>1.0) t0=1.0;如果(t1<0.0)t1=0.0;如果 (t1>1.0) t1=1.0;返回线(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));}最近的线(线 l0,三角形 t0){浮动 t;p点;线 cl,ll;cl.l=1e300;t=点(l0.p0-t0.p0,t0.n);p=l0.p0-t*t0.n;if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(l0.p0,p.p0);如果 (cl.l>ll.l) cl=ll;}t=点(l0.p1-t0.p0,t0.n);p=l0.p1-t*t0.n;if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(l0.p1,p.p0);如果 (cl.l>ll.l) cl=ll;}ll=最近的(l0,线(t0.p0,t0.p1));如果 (cl.l>ll.l) cl=ll;ll=最近的(l0,线(t0.p1,t0.p2));如果 (cl.l>ll.l) cl=ll;ll=最近的(l0,线(t0.p2,t0.p0));如果 (cl.l>ll.l) cl=ll;返回 cl;}最近的线(线 l0,convex_mesh m0){诠释我;线 cl,ll;cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0));cl.l=1e300;对于 (i=0;ill.l) cl=ll;}返回 cl;}//--------------------------------------------------------------------------线最近(三角形 t0,点 p0){ 线 cl;cl=最近的(p0,t0);cl.swap();返回 cl;}最近的线(三角形 t0,轴 a0){ 线 cl;cl=最近的(a0,t0);cl.swap();返回 cl;}线最近(三角形 t0,线 l0){ 线 cl;cl=最近的(l0,t0);cl.swap();返回 cl;}线最近(三角形 t0,三角形 t1){浮动 t;p点;第 l0,l1,l2,l3,l4,l5,cl,ll 行;l0=线(t0.p0,t0.p1);l3=线(t1.p0,t1.p1);l1=线(t0.p1,t0.p2);l4=线(t1.p1,t1.p2);l2=线(t0.p2,t0.p0);l5=线(t1.p2,t1.p0);cl.l=1e300;t=点(t0.p0-t1.p0,t1.n);p=t0.p0-t*t1.n;if ((fabs(t)>1e-6)&&(inside(p,t1))){ ll=line(t0.p0,p.p0);如果 (cl.l>ll.l) cl=ll;}t=点(t0.p1-t1.p0,t1.n);p=t0.p1-t*t1.n;if ((fabs(t)>1e-6)&&(inside(p,t1))){ ll=line(t0.p1,p.p0);如果 (cl.l>ll.l) cl=ll;}t=点(t0.p2-t1.p0,t1.n);p=t0.p2-t*t1.n;if ((fabs(t)>1e-6)&&(inside(p,t1))){ ll=line(t0.p2,p.p0);如果 (cl.l>ll.l) cl=ll;}t=点(t1.p0-t0.p0,t0.n);p=t1.p0-t*t0.n;if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p.p0,t1.p0);如果 (cl.l>ll.l) cl=ll;}t=点(t1.p1-t0.p0,t0.n);p=t1.p1-t*t0.n;if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p.p0,t1.p1);如果 (cl.l>ll.l) cl=ll;}t=点(t1.p2-t0.p0,t0.n);p=t1.p2-t*t0.n;if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p.p0,t1.p2);如果 (cl.l>ll.l) cl=ll;}ll=最近的(l0,l3);如果 (cl.l>ll.l) cl=ll;ll=最近的(l0,l4);如果 (cl.l>ll.l) cl=ll;ll=最近的(l0,l5);如果 (cl.l>ll.l) cl=ll;ll=最近的(l1,l3);如果 (cl.l>ll.l) cl=ll;ll=最近的(l1,l4);如果 (cl.l>ll.l) cl=ll;ll=最近的(l1,l5);如果 (cl.l>ll.l) cl=ll;ll=最近的(l2,l3);如果 (cl.l>ll.l) cl=ll;ll=最近的(l2,l4);如果 (cl.l>ll.l) cl=ll;ll=最近的(l2,l5);如果 (cl.l>ll.l) cl=ll;返回 cl;}线最近(三角形 t0,convex_mesh m0){诠释我;线 cl,ll;cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0));cl.l=1e300;对于 (i=0;ill.l) cl=ll;}返回 cl;}//--------------------------------------------------------------------------线最近(凸网格 m0,点 p0){ 线 cl;cl=最近的(p0,m0);cl.swap();返回 cl;}最接近的线(凸网格 m0,轴 a0){ 线 cl;cl=最近的(a0,m0);cl.swap();返回 cl;}最接近的线(凸网格 m0,线 l0){ 线 cl;cl=最近的(l0,m0);cl.swap();返回 cl;}最接近的线(凸网格 m0,三角形 t0){ 线 cl;cl=最近的(t0,m0);cl.swap();返回 cl;}线最近(凸网格 m0,凸网格 m1){诠释 i0, i1;线 cl,ll;cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0));cl.l=1e300;对于 (i0=0;i0ll.l) cl=ll;}返回 cl;}线最近(凸网格 m0,球形扇区 s0){诠释我,N = 18;浮动 a,da,ca,sa,cb,sb;vec3 u,v,w,q;线 cl,ll;//上限ll=最近的(m0,点(s0.p0));//球体if (dot(ll.dp,s0.dp)/(ll.l*s0.R)>=cos(s0.ang))//上限ll=线(ll.p0,ll.p1+(ll.dp*s0.R/ll.l));cl=ll;//圆锥w=标准化(s0.dp);u=vec3(1.0,0.0,0.0);if (fabs(dot(u,w))>0.75) u=vec3(0.0,1.0,0.0);v=交叉(u,w);u=交叉(v,w);u=标准化(u)*s0.r;v=标准化(v)*s0.r;da=2.0*M_PI/float(N-1);cb=cos(s0.ang);sb=sin(s0.ang);对于 (a=0.0,i=0;ill.l) cl=ll;}返回 cl;}//--------------------------------------------------------------------------bool intersect(convex_mesh m0,spherical_sector s0){线 cl;cl=最近的(m0,s0);if (cl.l<=1e-6) 返回真;如果(内部(cl.p0,s0))返回真;返回假;}//--------------------------------------------------------------------------#万一//--------------------------------------------------------------------------

GLSL 数学可以通过

锥体根据相交测试的结果旋转并改变颜色.黄线是最近的线结果.

我在这个周末为了好玩而破坏了这个,所以它还没有经过广泛的测试,并且可能仍然存在未处理的边缘情况.

我希望代码尽可能具有可读性,因此根本没有优化.我也没有过多评论(因为低级原语和基本向量数学应该足够明显,如果不是你应该在实现这样的东西之前先学习)

看起来我弄乱了一些 closest 测试,忽略了一些边缘情况.它需要进行重大返工(对所有相关功能应用修复),我现在没有足够的时间(将更新代码一旦完成)所以现在这里只对行与行最近的测试进行快速修复:

line3D 最接近(line3D l0,line3D l1){vec3 u=l0.p1-l0.p0;vec3 v=l1.p1-l1.p0;vec3 w=l0.p0-l1.p0;浮动一个=点(u,u);//总是 >= 0浮动 b=dot(u,v);浮动 c=点(v,v);//总是 >= 0浮动 d=dot(u,w);浮动 e=dot(v,w);浮动 D=a*c-b*b;//总是 >= 0浮动t0,t1;点3D p;线3D r,rr;诠释 f;//检查垂直于: 1: l0, 2: l1 的距离f=0;r.l=-1.0;//计算两个最近点的 line3D 参数if (D1.0){ f|=1;t0=1.0;}如果 (t1<0.0){ f|=2;t1=0.0;}如果 (t1>1.0){ f|=2;t1=1.0;}r=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));}//检查与 f 中标记的每个端点的垂直距离if (int(f&1)){t0=0.0;t1=除法(点(l0.p0-l1.p0,l1.dp),l1.l*l1.l);如果(t1<0.0)t1=0.0;如果 (t1>1.0) t1=1.0;rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));如果 ((r.l<0.0)||(r.l>rr.l)) r=rr;t0=1.0;t1=除法(点(l0.p1-l1.p0,l1.dp),l1.l*l1.l);如果(t1<0.0)t1=0.0;如果 (t1>1.0) t1=1.0;rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));如果 ((r.l<0.0)||(r.l>rr.l)) r=rr;}如果 (int(f&2)){t1=0.0;t0=除法(点(l1.p0-l0.p0,l0.dp),l0.l*l0.l);如果(t0<0.0)t0=0.0;如果 (t0>1.0) t0=1.0;rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));如果 ((r.l<0.0)||(r.l>rr.l)) r=rr;t1=1.0;t0=除(点(l1.p1-l0.p0,l0.dp),l0.l*l0.l);如果(t0<0.0)t0=0.0;如果 (t0>1.0) t0=1.0;rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));如果 ((r.l<0.0)||(r.l>rr.l)) r=rr;}返回 r;}

I'm looking to implement collision detection between a cone (With a round bottom. So it's basically a slice of a sphere) and a box. I'm not too fussed about it being AABB or OBB because transforming should be simple enough. Every solution I find uses a triangular cone but my cone is more of an "arc" that has an angle and distance.

Is there a simple solution to doing this collision detection? Or is it a case of doing several types of tests? ie. something like getting intersection points on a sphere with r being my cone distance then testing if they intersect within an angle or something?

解决方案

I was curious and planned to do stuff needed for this in GLSL math style anyway. So here a different approach. Let consider this definition of your cone:

  1. create a set of basic geometry primitives

You need to support points,lines,triangles,convex triangulated mesh,spherical sector (cone).

  1. implement inside test between point and triangle,mesh,cone

for triangle the results of cross between any side and point - side origin should point on the same side of triangle (like normal). If not point is outside.

for convex mesh dot product between point-face origin and face normal pointing out should be <=0 for all faces.

for cone the point should be inside sphere radius and angle between cone axis and point-cone origin should be <= ang. again dot product can be used for this.

  1. implement closest line between basic primitives

this is like finding closest points on each primitive that forms a line. Its similar to perpendicular distance.

point-point its easy as they are the closest line.

point-line can be done using projection of point onto the line (dot product). However you need to bound the result so it is inside line and not extrapolated beond it.

point-triangle can be obtained as minimum of all the circumference lines vs point combinations and perpendicular distance to surface (dot product with triangle normal).

All the other combinations of primitives can be build from these basic ones.

  1. closest line between mesh and cone

simply use closest line between cone sphere center and mesh. If the line lies inside cone shorten it by sphere radius R. This will account all cap interactions.

Then test lines on surface of the cone so sample along its circumference starting on cone sphere center and ending on the outermost circle (edge between cone and cap). You van test also triangles instead if you need better precision.

  1. intersection between mesh and cone

this one is easy just compute closest lien between mesh and cone. And then test if its point on mesh side is inside cone or not.

check the

    `bool intersect(convex_mesh m0,spherical_sector s0);`

implementation in the code below.

Here small C++/OpenGL example (using GLSL style math):

//---------------------------------------------------------------------------
//--- GL geometry -----------------------------------------------------------
//---------------------------------------------------------------------------
#ifndef _gl_geometry_h
#define _gl_geometry_h
//---------------------------------------------------------------------------
const float deg=M_PI/180.0;
const float rad=180.0/M_PI;
float divide(float a,float b){ if (fabs(b)<1e-10) return 0.0; else return a/b; }
double divide(double a,double b){ if (fabs(b)<1e-10) return 0.0; else return a/b; }
#include "GLSL_math.h"
#include "List.h"
//---------------------------------------------------------------------------
class point
    {
public:
    // cfg
    vec3 p0;

    point()     {}
    point(point& a) { *this=a; }
    ~point()    {}
    point* operator = (const point *a) { *this=*a; return this; }
    //point* operator = (const point &a) { ...copy... return this; }

    point(vec3 _p0)
        {
        p0=_p0;
        compute();
        }
    void compute(){};
    void draw()
        {
        glBegin(GL_POINTS);
        glVertex3fv(p0.dat);
        glEnd();
        }
    };
//---------------------------------------------------------------------------
class axis
    {
public:
    // cfg
    vec3 p0,dp;

    axis()      {}
    axis(axis& a)   { *this=a; }
    ~axis() {}
    axis* operator = (const axis *a) { *this=*a; return this; }
    //axis* operator = (const axis &a) { ...copy... return this; }

    axis(vec3 _p0,vec3 _dp)
        {
        p0=_p0;
        dp=_dp;
        compute();
        }
    void compute()
        {
        dp=normalize(dp);
        }
    void draw()
        {
        vec3 p; p=p0+100.0*dp;
        glBegin(GL_LINES);
        glVertex3fv(p0.dat);
        glVertex3fv(p .dat);
        glEnd();
        }
    };
//---------------------------------------------------------------------------
class line
    {
public:
    // cfg
    vec3 p0,p1;
    // computed
    float l;
    vec3 dp;

    line()  {}
    line(line& a)   { *this=a; }
    ~line() {}
    line* operator = (const line *a) { *this=*a; return this; }
    //line* operator = (const line &a) { ...copy... return this; }

    line(vec3 _p0,vec3 _p1)
        {
        p0=_p0;
        p1=_p1;
        compute();
        }
    void swap()
        {
        vec3 p=p0; p0=p1; p1=p;
        }
    void compute()
        {
        dp=p1-p0;
        l=length(dp);
        }
    void draw()
        {
        glBegin(GL_LINES);
        glVertex3fv(p0.dat);
        glVertex3fv(p1.dat);
        glEnd();
        }
    };
//---------------------------------------------------------------------------
class triangle
    {
public:
    // cfg
    vec3 p0,p1,p2;
    // computed
    vec3 n;

    triangle()  {}
    triangle(triangle& a)   { *this=a; }
    ~triangle() {}
    triangle* operator = (const triangle *a) { *this=*a; return this; }
    //triangle* operator = (const triangle &a) { ...copy... return this; }

    triangle(vec3 _p0,vec3 _p1,vec3 _p2)
        {
        p0=_p0;
        p1=_p1;
        p2=_p2;
        compute();
        }
    void swap()
        {
        vec3 p=p1; p1=p2; p2=p;
        n=-n;
        }
    void compute()
        {
        n=normalize(cross(p1-p0,p2-p1));
        }
    void draw()
        {
        glBegin(GL_TRIANGLES);
        glNormal3fv(n.dat);
        glVertex3fv(p0.dat);
        glVertex3fv(p1.dat);
        glVertex3fv(p2.dat);
        glEnd();
        }
    };
//---------------------------------------------------------------------------
class convex_mesh
    {
public:
    // cfg
    List<triangle> tri;
    // computed
    vec3 p0;            // center

    convex_mesh()   { tri.num=0; }
    convex_mesh(convex_mesh& a) { *this=a; }
    ~convex_mesh()  {}
    convex_mesh* operator = (const convex_mesh *a) { *this=*a; return this; }
    //convex_mesh* operator = (const convex_mesh &a) { ...copy... return this; }

    void init_box(vec3 _p0,vec3 _u,vec3 _v,vec3 _w) // center, half sizes
        {
        const vec3 p[8]=
            {
            _p0-_u+_v-_w,
            _p0+_u+_v-_w,
            _p0+_u-_v-_w,
            _p0-_u-_v-_w,
            _p0-_u-_v+_w,
            _p0+_u-_v+_w,
            _p0+_u+_v+_w,
            _p0-_u+_v+_w,
            };
        const int ix[36]=
            {
            0,1,2,0,2,3,
            4,5,6,4,6,7,
            3,2,5,3,5,4,
            2,1,6,2,6,5,
            1,0,7,1,7,6,
            0,3,4,0,4,7,
            };
        tri.num=0;
        for (int i=0;i<36;i+=3) tri.add(triangle(p[ix[i+0]],p[ix[i+1]],p[ix[i+2]]));
        compute();
        }
    void compute()
        {
        int i,n;
        p0=vec3(0.0,0.0,0.0);
        if (!tri.num) return;
        for (i=0,n=0;i<tri.num;i++,n+=3)
            {
            p0+=tri.dat[i].p0;
            p0+=tri.dat[i].p1;
            p0+=tri.dat[i].p2;
            } p0/=float(n);
        for (i=0;i<tri.num;i++)
         if (dot(tri.dat[i].p0-p0,tri.dat[i].n)<0.0)
          tri.dat[i].swap();
        }
    void draw()
        {
        int i;
        glBegin(GL_TRIANGLES);
        for (i=0;i<tri.num;i++) tri.dat[i].draw();
        glEnd();
        }
    };
//---------------------------------------------------------------------------
class spherical_sector
    {
public:
    // cfg
    vec3 p0,p1;
    float ang;
    // computed
    vec3 dp;
    float r,R;

    spherical_sector()  {}
    spherical_sector(spherical_sector& a)   { *this=a; }
    ~spherical_sector() {}
    spherical_sector* operator = (const spherical_sector *a) { *this=*a; return this; }
    //spherical_sector* operator = (const spherical_sector &a) { ...copy... return this; }

    spherical_sector(vec3 _p0,vec3 _p1,float _ang)
        {
        p0=_p0;
        p1=_p1;
        ang=_ang;
        compute();
        }
    void compute()
        {
        dp=p1-p0;
        R=length(dp);
        r=R*tan(ang);
        }
    void draw()
        {
        const int N=32;
        const int M=16;
        vec3 pnt[M][N]; // points
        vec3 n0[N];     // normals for cine
        vec3 n1[M][N];  // normals for cap
        int i,j;
        float a,b,da,db,ca,sa,cb,sb;
        vec3 q,u,v,w;
        // basis vectors
        w=normalize(dp);         u=vec3(1.0,0.0,0.0);
        if (fabs(dot(u,w))>0.75) u=vec3(0.0,1.0,0.0);
        v=cross(u,w);
        u=cross(v,w);
        u=normalize(u);
        v=normalize(v);
        // compute tables
        da=2.0*M_PI/float(N-1);
        db=ang/float(M-1);
        for (a=0.0,i=0;i<N;i++,a+=da)
            {
            ca=cos(a);
            sa=sin(a);
            n0[i]=u*ca+v*sa;
            for (b=0.0,j=0;j<M;j++,b+=db)
                {
                cb=cos(b);
                sb=sin(b);
                q=vec3(ca*sb,sa*sb,cb);
                pnt[j][i]=p0+((q.x*u+q.y*v+q.z*w)*R);
                n1[j][i]=normalize(pnt[j][i]);
                }
            }
        // render
        glBegin(GL_TRIANGLES);
        for (i=1,j=M-1;i<N;i++)
            {
            glNormal3fv(n0[i].dat);         // p0 should have average 0.5*(n0[i]+n0[i-1]) as nomal
            glVertex3fv(p0.dat);
            glVertex3fv(pnt[j][i+0].dat);
            glNormal3fv(n0[i-1].dat);
            glVertex3fv(pnt[j][i-1].dat);
            glNormal3fv( n1[0][0].dat);
            glVertex3fv(pnt[0][0].dat);
            glNormal3fv( n1[1][i-1].dat);
            glVertex3fv(pnt[1][i-1].dat);
            glNormal3fv( n1[1][i+0].dat);
            glVertex3fv(pnt[1][i+0].dat);
            }
        glEnd();
        glBegin(GL_QUADS);
        for (i=0;i<N;i++)
         for (j=2;j<M;j++)
            {
            glNormal3fv( n1[j-1][i+0].dat);
            glVertex3fv(pnt[j-1][i+0].dat);
            glNormal3fv( n1[j-1][i-1].dat);
            glVertex3fv(pnt[j-1][i-1].dat);
            glNormal3fv( n1[j+0][i-1].dat);
            glVertex3fv(pnt[j+0][i-1].dat);
            glNormal3fv( n1[j+0][i+0].dat);
            glVertex3fv(pnt[j+0][i+0].dat);
            }
        glEnd();
        }
    };
//---------------------------------------------------------------------------
//---------------------------------------------------------------------------
bool inside(point p0,triangle t0);
bool inside(point p0,convex_mesh m0);
bool inside(point p0,spherical_sector s0);
//---------------------------------------------------------------------------
line closest(point p0,axis a0);
line closest(point p0,line l0);
line closest(point p0,triangle t0);
line closest(point p0,convex_mesh m0);
//---------------------------------------------------------------------------
line closest(axis a0,point p0);
line closest(axis a0,axis  a1);
line closest(axis a0,line  l1);
line closest(axis a0,triangle t0);
line closest(axis a0,convex_mesh m0);
//---------------------------------------------------------------------------
line closest(line l0,point p0);
line closest(line l0,axis  a0);
line closest(line l0,line  l1);
line closest(line l0,triangle t0);
line closest(line l0,convex_mesh m0);
//---------------------------------------------------------------------------
line closest(triangle t0,point p0);
line closest(triangle t0,axis  a0);
line closest(triangle t0,line  l0);
line closest(triangle t0,triangle t1);
line closest(triangle t0,convex_mesh m0);
//---------------------------------------------------------------------------
line closest(convex_mesh m0,point p0);
line closest(convex_mesh m0,axis  a0);
line closest(convex_mesh m0,line  l0);
line closest(convex_mesh m0,triangle t0);
line closest(convex_mesh m0,spherical_sector s0);
//---------------------------------------------------------------------------
bool intersect(convex_mesh m0,spherical_sector s0);
//---------------------------------------------------------------------------
//---------------------------------------------------------------------------
bool inside(point p0,triangle t0)
    {
    if (fabs(dot(p0.p0-t0.p0,t0.n))>1e-6) return false;
    float d0,d1,d2;
    d0=dot(t0.n,cross(p0.p0-t0.p0,t0.p1-t0.p0));
    d1=dot(t0.n,cross(p0.p0-t0.p1,t0.p2-t0.p1));
    d2=dot(t0.n,cross(p0.p0-t0.p2,t0.p0-t0.p2));
    if (d0*d1<-1e-6) return false;
    if (d0*d2<-1e-6) return false;
    if (d1*d2<-1e-6) return false;
    return true;
    }
bool inside(point p0,convex_mesh m0)
    {
    for (int i=0;i<m0.tri.num;i++)
     if (dot(p0.p0-m0.tri.dat[i].p0,m0.tri.dat[i].n)>0.0)
      return false;
    return true;
    }
bool inside(point p0,spherical_sector s0)
    {
    float t,l;
    vec3 u;
    u=p0.p0-s0.p0;
    l=length(u);
    if (l>s0.R) return false;
    t=divide(dot(u,s0.dp),(l*s0.R));
    if (t<cos(s0.ang)) return false;
    return true;
    }
//---------------------------------------------------------------------------
line closest(point p0,axis a0){ return line(p0.p0,a0.p0+(a0.dp*dot(p0.p0-a0.p0,a0.dp))); }
line closest(point p0,line l0)
    {
    float t=dot(p0.p0-l0.p0,l0.dp);
    if (t<0.0) t=0.0;
    if (t>1.0) t=1.0;
    return line(p0.p0,l0.p0+(l0.dp*t));
    }
line closest(point p0,triangle t0)
    {
    float t;
    point p;
    line cl,ll;
    cl.l=1e300;
    t=dot(p0.p0-t0.p0,t0.n); p=p0.p0-t*t0.n; if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p0.p0,p.p0); if (cl.l>ll.l) cl=ll; }
    ll=closest(p0,line(t0.p0,t0.p1)); if (cl.l>ll.l) cl=ll;
    ll=closest(p0,line(t0.p1,t0.p2)); if (cl.l>ll.l) cl=ll;
    ll=closest(p0,line(t0.p2,t0.p0)); if (cl.l>ll.l) cl=ll;
    return cl;
    }
line closest(point p0,convex_mesh m0)
    {
    int i;
    line cl,ll;
    cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0)); cl.l=1e300;
    for (i=0;i<m0.tri.num;i++)
        {
        ll=closest(p0,m0.tri.dat[i]);
        if (cl.l>ll.l) cl=ll;
        }
    return cl;
    }
//---------------------------------------------------------------------------
line closest(axis a0,point p0){ line cl; cl=closest(p0,a0); cl.swap(); return cl; }
line closest(axis a0,axis a1)
    {
    vec3 u=a0.dp;
    vec3 v=a1.dp;
    vec3 w=a0.p0-a1.p0;
    float a=dot(u,u);       // always >= 0
    float b=dot(u,v);
    float c=dot(v,v);       // always >= 0
    float d=dot(u,w);
    float e=dot(v,w);
    float D=a*c-b*b;        // always >= 0
    float t0,t1;
    // compute the line parameters of the two closest points
    if (D<1e-6)            // the lines are almost parallel
        {
        t0=0.0;
        t1=(b>c ? d/b : e/c); // use the largest denominator
        }
    else{
        t0=(b*e-c*d)/D;
        t1=(a*e-b*d)/D;
        }
    return line(a0.p0+(a0.dp*t0),a1.p0+(a1.dp*t1));
    }
line closest(axis a0,line l1)
    {
    vec3 u=a0.dp;
    vec3 v=l1.dp;
    vec3 w=a0.p0-l1.p0;
    float a=dot(u,u);       // always >= 0
    float b=dot(u,v);
    float c=dot(v,v);       // always >= 0
    float d=dot(u,w);
    float e=dot(v,w);
    float D=a*c-b*b;        // always >= 0
    float t0,t1;
    // compute the line parameters of the two closest points
    if (D<1e-6)            // the lines are almost parallel
        {
        t0=0.0;
        t1=(b>c ? d/b : e/c); // use the largest denominator
        }
    else{
        t0=(b*e-c*d)/D;
        t1=(a*e-b*d)/D;
        }
    if (t1<0.0) t1=0.0;
    if (t1>1.0) t1=1.0;
    return line(a0.p0+(a0.dp*t0),l1.p0+(l1.dp*t1));
    }
line closest(axis a0,triangle t0)
    {
    line cl,ll;
    cl=closest(a0,line(t0.p0,t0.p1));
    ll=closest(a0,line(t0.p1,t0.p2)); if (cl.l>ll.l) cl=ll;
    ll=closest(a0,line(t0.p2,t0.p0)); if (cl.l>ll.l) cl=ll;
    return cl;
    }
line closest(axis a0,convex_mesh m0)
    {
    int i;
    line cl,ll;
    cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0)); cl.l=1e300;
    for (i=0;i<m0.tri.num;i++)
        {
        ll=closest(a0,m0.tri.dat[i]);
        if (cl.l>ll.l) cl=ll;
        }
    return cl;
    }
//---------------------------------------------------------------------------
line closest(line l0,point p0){ line cl; cl=closest(p0,l0); cl.swap(); return cl; }
line closest(line l0,axis a0) { line cl; cl=closest(a0,l0); cl.swap(); return cl; }
line closest(line l0,line l1)
    {
    vec3 u=l0.p1-l0.p0;
    vec3 v=l1.p1-l1.p0;
    vec3 w=l0.p0-l1.p0;
    float a=dot(u,u);       // always >= 0
    float b=dot(u,v);
    float c=dot(v,v);       // always >= 0
    float d=dot(u,w);
    float e=dot(v,w);
    float D=a*c-b*b;        // always >= 0
    float t0,t1;
    // compute the line parameters of the two closest points
    if (D<1e-6)            // the lines are almost parallel
        {
        t0=0.0;
        t1=(b>c ? d/b : e/c); // use the largest denominator
        }
    else{
        t0=(b*e-c*d)/D;
        t1=(a*e-b*d)/D;
        }
    if (t0<0.0) t0=0.0;
    if (t0>1.0) t0=1.0;
    if (t1<0.0) t1=0.0;
    if (t1>1.0) t1=1.0;
    return line(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));
    }
line closest(line l0,triangle t0)
    {
    float t;
    point p;
    line cl,ll;
    cl.l=1e300;
    t=dot(l0.p0-t0.p0,t0.n); p=l0.p0-t*t0.n; if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(l0.p0,p.p0); if (cl.l>ll.l) cl=ll; }
    t=dot(l0.p1-t0.p0,t0.n); p=l0.p1-t*t0.n; if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(l0.p1,p.p0); if (cl.l>ll.l) cl=ll; }
    ll=closest(l0,line(t0.p0,t0.p1)); if (cl.l>ll.l) cl=ll;
    ll=closest(l0,line(t0.p1,t0.p2)); if (cl.l>ll.l) cl=ll;
    ll=closest(l0,line(t0.p2,t0.p0)); if (cl.l>ll.l) cl=ll;
    return cl;
    }
line closest(line l0,convex_mesh m0)
    {
    int i;
    line cl,ll;
    cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0)); cl.l=1e300;
    for (i=0;i<m0.tri.num;i++)
        {
        ll=closest(l0,m0.tri.dat[i]);
        if (cl.l>ll.l) cl=ll;
        }
    return cl;
    }
//---------------------------------------------------------------------------
line closest(triangle t0,point p0){ line cl; cl=closest(p0,t0); cl.swap(); return cl; }
line closest(triangle t0,axis a0) { line cl; cl=closest(a0,t0); cl.swap(); return cl; }
line closest(triangle t0,line l0) { line cl; cl=closest(l0,t0); cl.swap(); return cl; }
line closest(triangle t0,triangle t1)
    {
    float t;
    point p;
    line l0,l1,l2,l3,l4,l5,cl,ll;
    l0=line(t0.p0,t0.p1); l3=line(t1.p0,t1.p1);
    l1=line(t0.p1,t0.p2); l4=line(t1.p1,t1.p2);
    l2=line(t0.p2,t0.p0); l5=line(t1.p2,t1.p0);
    cl.l=1e300;
    t=dot(t0.p0-t1.p0,t1.n); p=t0.p0-t*t1.n; if ((fabs(t)>1e-6)&&(inside(p,t1))){ ll=line(t0.p0,p.p0); if (cl.l>ll.l) cl=ll; }
    t=dot(t0.p1-t1.p0,t1.n); p=t0.p1-t*t1.n; if ((fabs(t)>1e-6)&&(inside(p,t1))){ ll=line(t0.p1,p.p0); if (cl.l>ll.l) cl=ll; }
    t=dot(t0.p2-t1.p0,t1.n); p=t0.p2-t*t1.n; if ((fabs(t)>1e-6)&&(inside(p,t1))){ ll=line(t0.p2,p.p0); if (cl.l>ll.l) cl=ll; }
    t=dot(t1.p0-t0.p0,t0.n); p=t1.p0-t*t0.n; if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p.p0,t1.p0); if (cl.l>ll.l) cl=ll; }
    t=dot(t1.p1-t0.p0,t0.n); p=t1.p1-t*t0.n; if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p.p0,t1.p1); if (cl.l>ll.l) cl=ll; }
    t=dot(t1.p2-t0.p0,t0.n); p=t1.p2-t*t0.n; if ((fabs(t)>1e-6)&&(inside(p,t0))){ ll=line(p.p0,t1.p2); if (cl.l>ll.l) cl=ll; }
    ll=closest(l0,l3); if (cl.l>ll.l) cl=ll;
    ll=closest(l0,l4); if (cl.l>ll.l) cl=ll;
    ll=closest(l0,l5); if (cl.l>ll.l) cl=ll;
    ll=closest(l1,l3); if (cl.l>ll.l) cl=ll;
    ll=closest(l1,l4); if (cl.l>ll.l) cl=ll;
    ll=closest(l1,l5); if (cl.l>ll.l) cl=ll;
    ll=closest(l2,l3); if (cl.l>ll.l) cl=ll;
    ll=closest(l2,l4); if (cl.l>ll.l) cl=ll;
    ll=closest(l2,l5); if (cl.l>ll.l) cl=ll;
    return cl;
    }
line closest(triangle t0,convex_mesh m0)
    {
    int i;
    line cl,ll;
    cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0)); cl.l=1e300;
    for (i=0;i<m0.tri.num;i++)
        {
        ll=closest(m0.tri.dat[i],t0);
        if (cl.l>ll.l) cl=ll;
        }
    return cl;
    }
//---------------------------------------------------------------------------
line closest(convex_mesh m0,point p0)   { line cl; cl=closest(p0,m0); cl.swap(); return cl; }
line closest(convex_mesh m0,axis a0)    { line cl; cl=closest(a0,m0); cl.swap(); return cl; }
line closest(convex_mesh m0,line l0)    { line cl; cl=closest(l0,m0); cl.swap(); return cl; }
line closest(convex_mesh m0,triangle t0){ line cl; cl=closest(t0,m0); cl.swap(); return cl; }
line closest(convex_mesh m0,convex_mesh m1)
    {
    int i0,i1;
    line cl,ll;
    cl=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0)); cl.l=1e300;
    for (i0=0;i0<m0.tri.num;i0++)
     for (i1=0;i1<m1.tri.num;i1++)
        {
        ll=closest(m0.tri.dat[i0],m1.tri.dat[i1]);
        if (cl.l>ll.l) cl=ll;
        }
    return cl;
    }
line closest(convex_mesh m0,spherical_sector s0)
    {
    int i,N=18;
    float a,da,ca,sa,cb,sb;
    vec3 u,v,w,q;
    line cl,ll;
    // cap
    ll=closest(m0,point(s0.p0));                    // sphere
    if (dot(ll.dp,s0.dp)/(ll.l*s0.R)>=cos(s0.ang))  // cap
     ll=line(ll.p0,ll.p1+(ll.dp*s0.R/ll.l));
    cl=ll;
    // cone
    w=normalize(s0.dp);      u=vec3(1.0,0.0,0.0);
    if (fabs(dot(u,w))>0.75) u=vec3(0.0,1.0,0.0);
    v=cross(u,w);
    u=cross(v,w);
    u=normalize(u)*s0.r;
    v=normalize(v)*s0.r;
    da=2.0*M_PI/float(N-1);
    cb=cos(s0.ang);
    sb=sin(s0.ang);
    for (a=0.0,i=0;i<N;i++)
        {
        ca=cos(a);
        sa=sin(a);
        q=vec3(ca*sb,sa*sb,cb);
        q=s0.p0+((q.x*u+q.y*v+q.z*w)*s0.R);
        ll=line(s0.p0,q);
        ll=closest(m0,ll);
        if (cl.l>ll.l) cl=ll;
        }
    return cl;
    }
//---------------------------------------------------------------------------
bool intersect(convex_mesh m0,spherical_sector s0)
    {
    line cl;
    cl=closest(m0,s0);
    if (cl.l<=1e-6) return true;
    if (inside(cl.p0,s0)) return true;
    return false;
    }
//---------------------------------------------------------------------------
#endif
//---------------------------------------------------------------------------

The GLSL math can be created by this or use GLM or whatever else instead.

I also used mine dynamic list template (just to store the list of triangles in mesh) so:


List<double> xxx; is the same as double xxx[];
xxx.add(5); adds 5 to end of the list
xxx[7] access array element (safe)
xxx.dat[7] access array element (unsafe but fast direct access)
xxx.num is the actual used size of the array
xxx.reset() clears the array and set xxx.num=0
xxx.allocate(100) preallocate space for 100 items

You can use whatever list you have at disposal.

And here test preview testing correctness of this:

The cone is rotating and changing color according to result of intersect test. The yellow line is the closest line result.

I busted this for fun during this weekend so its not yet extensively tested and there still might be unhandled edge cases.

I wanted the code to be as readable as I could so its not optimized at all. Also I did not comment much (as the low level primitives and basic vector math should be obvious enough if not you should learn first before implementing stuff like this)

[edit1]

Looks like I mess some closest test ignoring some edge cases.. It need major rework (apply the fix for all related functions) which I do not have enough time for right now (will update the code once finished) so for now here a quick fix for line vs. line closest test only:

line3D closest(line3D l0,line3D l1)
    {
    vec3 u=l0.p1-l0.p0;
    vec3 v=l1.p1-l1.p0;
    vec3 w=l0.p0-l1.p0;
    float a=dot(u,u);       // always >= 0
    float b=dot(u,v);
    float c=dot(v,v);       // always >= 0
    float d=dot(u,w);
    float e=dot(v,w);
    float D=a*c-b*b;        // always >= 0
    float t0,t1;
    point3D p;
    line3D r,rr;
    int f;                  // check distance perpendicular to: 1: l0, 2: l1
    f=0; r.l=-1.0;
    // compute the line3D parameters of the two closest points
    if (D<acc_dot) f=3;    // the lines are almost parallel
    else{
        t0=(b*e-c*d)/D;
        t1=(a*e-b*d)/D;
        if (t0<0.0){ f|=1; t0=0.0; }
        if (t0>1.0){ f|=1; t0=1.0; }
        if (t1<0.0){ f|=2; t1=0.0; }
        if (t1>1.0){ f|=2; t1=1.0; }
        r=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));
        }
    // check perpendicular distance from each endpoint marked in f
    if (int(f&1))
        {
        t0=0.0;
        t1=divide(dot(l0.p0-l1.p0,l1.dp),l1.l*l1.l);
        if (t1<0.0) t1=0.0;
        if (t1>1.0) t1=1.0;
        rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));
        if ((r.l<0.0)||(r.l>rr.l)) r=rr;
        t0=1.0;
        t1=divide(dot(l0.p1-l1.p0,l1.dp),l1.l*l1.l);
        if (t1<0.0) t1=0.0;
        if (t1>1.0) t1=1.0;
        rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));
        if ((r.l<0.0)||(r.l>rr.l)) r=rr;
        }
    if (int(f&2))
        {
        t1=0.0;
        t0=divide(dot(l1.p0-l0.p0,l0.dp),l0.l*l0.l);
        if (t0<0.0) t0=0.0;
        if (t0>1.0) t0=1.0;
        rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));
        if ((r.l<0.0)||(r.l>rr.l)) r=rr;
        t1=1.0;
        t0=divide(dot(l1.p1-l0.p0,l0.dp),l0.l*l0.l);
        if (t0<0.0) t0=0.0;
        if (t0>1.0) t0=1.0;
        rr=line3D(l0.p0+(l0.dp*t0),l1.p0+(l1.dp*t1));
        if ((r.l<0.0)||(r.l>rr.l)) r=rr;
        }
    return r;
    }

这篇关于锥盒碰撞的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆