OpenSimMirror/libraries/ode-0.9/contrib/TerrainAndCone/dCone.cpp

505 lines
13 KiB
C++

//Benoit CHAPEROT 2003-2004 www.jstarlab.com
//some code inspired by Magic Software
#include <ode/common.h>
#include <ode/collision.h>
#include <ode/matrix.h>
#include <ode/rotation.h>
#include <ode/odemath.h>
#include "collision_kernel.h"
#include "collision_std.h"
#include "collision_std_internal.h"
#include "collision_util.h"
#include <drawstuff/drawstuff.h>
#include "windows.h"
#include "ode\ode.h"
#define CONTACT(p,skip) ((dContactGeom*) (((char*)p) + (skip)))
const dReal fEPSILON = 1e-9f;
dxCone::dxCone (dSpaceID space, dReal _radius,dReal _length) :
dxGeom (space,1)
{
dAASSERT(_radius > 0.f);
dAASSERT(_length > 0.f);
type = dConeClass;
radius = _radius;
lz = _length;
}
dxCone::~dxCone()
{
}
void dxCone::computeAABB()
{
const dMatrix3& R = final_posr->R;
const dVector3& pos = final_posr->pos;
dReal xrange = dFabs(R[2] * lz) + radius;
dReal yrange = dFabs(R[6] * lz) + radius;
dReal zrange = dFabs(R[10] * lz) + radius;
aabb[0] = pos[0] - xrange;
aabb[1] = pos[0] + xrange;
aabb[2] = pos[1] - yrange;
aabb[3] = pos[1] + yrange;
aabb[4] = pos[2] - zrange;
aabb[5] = pos[2] + zrange;
}
dGeomID dCreateCone(dSpaceID space, dReal _radius,dReal _length)
{
return new dxCone(space,_radius,_length);
}
void dGeomConeSetParams (dGeomID g, dReal _radius, dReal _length)
{
dUASSERT (g && g->type == dConeClass,"argument not a cone");
dAASSERT (_radius > 0.f);
dAASSERT (_length > 0.f);
g->recomputePosr();
dxCone *c = (dxCone*) g;
c->radius = _radius;
c->lz = _length;
dGeomMoved (g);
}
void dGeomConeGetParams (dGeomID g, dReal *_radius, dReal *_length)
{
dUASSERT (g && g->type == dConeClass,"argument not a cone");
g->recomputePosr();
dxCone *c = (dxCone*) g;
*_radius = c->radius;
*_length = c->lz;
}
//positive inside
dReal dGeomConePointDepth(dGeomID g, dReal x, dReal y, dReal z)
{
dUASSERT (g && g->type == dConeClass,"argument not a cone");
g->recomputePosr();
dxCone *cone = (dxCone*) g;
dVector3 tmp,q;
tmp[0] = x - cone->final_posr->pos[0];
tmp[1] = y - cone->final_posr->pos[1];
tmp[2] = z - cone->final_posr->pos[2];
dMULTIPLY1_331 (q,cone->final_posr->R,tmp);
dReal r = cone->radius;
dReal h = cone->lz;
dReal d0 = (r - r*q[2]/h) - dSqrt(q[0]*q[0]+q[1]*q[1]);
dReal d1 = q[2];
dReal d2 = h-q[2];
if (d0 < d1) {
if (d0 < d2) return d0; else return d2;
}
else {
if (d1 < d2) return d1; else return d2;
}
}
//plane plane
bool FindIntersectionPlanePlane(const dReal Plane0[4], const dReal Plane1[4],
dVector3 LinePos,dVector3 LineDir)
{
// If Cross(N0,N1) is zero, then either planes are parallel and separated
// or the same plane. In both cases, 'false' is returned. Otherwise,
// the intersection line is
//
// L(t) = t*Cross(N0,N1) + c0*N0 + c1*N1
//
// for some coefficients c0 and c1 and for t any real number (the line
// parameter). Taking dot products with the normals,
//
// d0 = Dot(N0,L) = c0*Dot(N0,N0) + c1*Dot(N0,N1)
// d1 = Dot(N1,L) = c0*Dot(N0,N1) + c1*Dot(N1,N1)
//
// which are two equations in two unknowns. The solution is
//
// c0 = (Dot(N1,N1)*d0 - Dot(N0,N1)*d1)/det
// c1 = (Dot(N0,N0)*d1 - Dot(N0,N1)*d0)/det
//
// where det = Dot(N0,N0)*Dot(N1,N1)-Dot(N0,N1)^2.
/*
Real fN00 = rkPlane0.Normal().SquaredLength();
Real fN01 = rkPlane0.Normal().Dot(rkPlane1.Normal());
Real fN11 = rkPlane1.Normal().SquaredLength();
Real fDet = fN00*fN11 - fN01*fN01;
if ( Math::FAbs(fDet) < gs_fEpsilon )
return false;
Real fInvDet = 1.0f/fDet;
Real fC0 = (fN11*rkPlane0.Constant() - fN01*rkPlane1.Constant())*fInvDet;
Real fC1 = (fN00*rkPlane1.Constant() - fN01*rkPlane0.Constant())*fInvDet;
rkLine.Direction() = rkPlane0.Normal().Cross(rkPlane1.Normal());
rkLine.Origin() = fC0*rkPlane0.Normal() + fC1*rkPlane1.Normal();
return true;
*/
dReal fN00 = dLENGTHSQUARED(Plane0);
dReal fN01 = dDOT(Plane0,Plane1);
dReal fN11 = dLENGTHSQUARED(Plane1);
dReal fDet = fN00*fN11 - fN01*fN01;
if ( fabs(fDet) < fEPSILON)
return false;
dReal fInvDet = 1.0f/fDet;
dReal fC0 = (fN11*Plane0[3] - fN01*Plane1[3])*fInvDet;
dReal fC1 = (fN00*Plane1[3] - fN01*Plane0[3])*fInvDet;
dCROSS(LineDir,=,Plane0,Plane1);
dNormalize3(LineDir);
dVector3 Temp0,Temp1;
dOPC(Temp0,*,Plane0,fC0);
dOPC(Temp1,*,Plane1,fC1);
dOP(LinePos,+,Temp0,Temp1);
return true;
}
//plane ray
bool FindIntersectionPlaneRay(const dReal Plane[4],
const dVector3 &LinePos,const dVector3 &LineDir,
dReal &u,dVector3 &Pos)
{
/*
u = (A*X1 + B*Y1 + C*Z1 + D) / (A*(X1-X2) + B*(Y1-Y2)+C*(Z1-Z2))
*/
dReal fDet = -dDot(Plane,LineDir,3);
if ( fabs(fDet) < fEPSILON)
return false;
u = (dDot(Plane,LinePos,3) - Plane[3]) / fDet;
dOPC(Pos,*,LineDir,u);
dOPE(Pos,+=,LinePos);
return true;
}
int SolveQuadraticPolynomial(dReal a,dReal b,dReal c,dReal &x0,dReal &x1)
{
dReal d = b*b - 4*a*c;
int NumRoots = 0;
dReal dr;
if (d < 0.f)
return NumRoots;
if (d == 0.f)
{
NumRoots = 1;
dr = 0.f;
}
else
{
NumRoots = 2;
dr = sqrtf(d);
}
x0 = (-b -dr) / (2.f * a);
x1 = (-b +dr) / (2.f * a);
return NumRoots;
}
/*
const int VALID_INTERSECTION = 1<<0;
const int POS_TEST_FAILEDT0 = 1<<0;
const int POS_TEST_FAILEDT1 = 1<<1;
*/
int ProcessConeRayIntersectionPoint( dReal r,dReal h,
const dVector3 &q,const dVector3 &v,dReal t,
dVector3 &p,
dVector3 &n,
int &f)
{
dOPC(p,*,v,t);
dOPE(p,+=,q);
n[0] = 2*p[0];
n[1] = 2*p[1];
n[2] = -2*p[2]*r*r/(h*h);
f = 0;
if (p[2] > h) return 0;
if (p[2] < 0) return 0;
if (t > 1) return 0;
if (t < 0) return 0;
return 1;
}
//cone ray
//line in cone space (position,direction)
//distance from line position (direction normalized)(if any)
//return the number of intersection
int FindIntersectionConeRay(dReal r,dReal h,
const dVector3 &q,const dVector3 &v,dContactGeom *pContact)
{
dVector3 qp,vp;
dOPE(qp,=,q);
dOPE(vp,=,v);
qp[2] = h-q[2];
vp[2] = -v[2];
dReal ts = (r/h);
ts *= ts;
dReal a = vp[0]*vp[0] + vp[1]*vp[1] - ts*vp[2]*vp[2];
dReal b = 2.f*qp[0]*vp[0] + 2.f*qp[1]*vp[1] - 2.f*ts*qp[2]*vp[2];
dReal c = qp[0]*qp[0] + qp[1]*qp[1] - ts*qp[2]*qp[2];
/*
dReal a = v[0]*v[0] + v[1]*v[1] - (v[2]*v[2]*r*r) / (h*h);
dReal b = 2.f*q[0]*v[0] + 2.f*q[1]*v[1] + 2.f*r*r*v[2]/h - 2*r*r*q[0]*v[0]/(h*h);
dReal c = q[0]*q[0] + q[1]*q[1] + 2*r*r*q[2]/h - r*r*q[2]/(h*h) - r*r;
*/
int nNumRoots=SolveQuadraticPolynomial(a,b,c,pContact[0].depth,pContact[1].depth);
int flag = 0;
dContactGeom ValidContact[2];
int nNumValidContacts = 0;
for (int i=0;i<nNumRoots;i++)
{
if (ProcessConeRayIntersectionPoint(r,h,q,v,pContact[i].depth,pContact[i].pos,
pContact[i].normal,flag))
{
ValidContact[nNumValidContacts] = pContact[i];
nNumValidContacts++;
}
}
dOP(qp,+,q,v);
if ((nNumValidContacts < 2) && (v[2] != 0.f))
{
dReal d = (0.f-q[2]) / (v[2]);
if ((d>=0) && (d<=1))
{
dOPC(vp,*,v,d);
dOP(qp,+,q,vp);
if (qp[0]*qp[0]+qp[1]*qp[1] < r*r)
{
dOPE(ValidContact[nNumValidContacts].pos,=,qp);
ValidContact[nNumValidContacts].normal[0] = 0.f;
ValidContact[nNumValidContacts].normal[1] = 0.f;
ValidContact[nNumValidContacts].normal[2] = -1.f;
ValidContact[nNumValidContacts].depth = d;
nNumValidContacts++;
}
}
}
if (nNumValidContacts == 2)
{
if (ValidContact[0].depth > ValidContact[1].depth)
{
pContact[0] = ValidContact[1];
pContact[1] = ValidContact[0];
}
else
{
pContact[0] = ValidContact[0];
pContact[1] = ValidContact[1];
}
}
else if (nNumValidContacts == 1)
{
pContact[0] = ValidContact[0];
}
return nNumValidContacts;
}
int dCollideConePlane (dxGeom *o1, dxGeom *o2, int flags,
dContactGeom *contact, int skip)
{
dIASSERT (skip >= (int)sizeof(dContactGeom));
dIASSERT (o1->type == dConeClass);
dIASSERT (o2->type == dPlaneClass);
dxCone *cone = (dxCone*) o1;
dxPlane *plane = (dxPlane*) o2;
contact->g1 = o1;
contact->g2 = o2;
dVector3 p0,p1,pp0,pp1;
dOPE(p0,=,cone->final_posr->pos);
p1[0] = cone->final_posr->R[0*4+2] * cone->lz + p0[0];
p1[1] = cone->final_posr->R[1*4+2] * cone->lz + p0[1];
p1[2] = cone->final_posr->R[2*4+2] * cone->lz + p0[2];
dReal u;
FindIntersectionPlaneRay(plane->p,p0,plane->p,u,pp0);
FindIntersectionPlaneRay(plane->p,p1,plane->p,u,pp1);
if (dDISTANCE(pp0,pp1) < fEPSILON)
{
p1[0] = cone->final_posr->R[0*4+0] * cone->lz + p0[0];
p1[1] = cone->final_posr->R[1*4+0] * cone->lz + p0[1];
p1[2] = cone->final_posr->R[2*4+0] * cone->lz + p0[2];
FindIntersectionPlaneRay(plane->p,p1,plane->p,u,pp1);
dIASSERT(dDISTANCE(pp0,pp1) >= fEPSILON);
}
dVector3 h,r0,r1;
h[0] = cone->final_posr->R[0*4+2];
h[1] = cone->final_posr->R[1*4+2];
h[2] = cone->final_posr->R[2*4+2];
dOP(r0,-,pp0,pp1);
dCROSS(r1,=,h,r0);
dCROSS(r0,=,r1,h);
dNormalize3(r0);
dOPEC(h,*=,cone->lz);
dOPEC(r0,*=,cone->radius);
dVector3 p[3];
dOP(p[0],+,cone->final_posr->pos,h);
dOP(p[1],+,cone->final_posr->pos,r0);
dOP(p[2],-,cone->final_posr->pos,r0);
int numMaxContacts = flags & 0xffff;
if (numMaxContacts == 0)
numMaxContacts = 1;
int n=0;
for (int i=0;i<3;i++)
{
dReal d = dGeomPlanePointDepth(o2, p[i][0], p[i][1], p[i][2]);
if (d>0.f)
{
CONTACT(contact,n*skip)->g1 = o1;
CONTACT(contact,n*skip)->g2 = o2;
dOPE(CONTACT(contact,n*skip)->normal,=,plane->p);
dOPE(CONTACT(contact,n*skip)->pos,=,p[i]);
CONTACT(contact,n*skip)->depth = d;
n++;
if (n == numMaxContacts)
return n;
}
}
return n;
}
int dCollideRayCone (dxGeom *o1, dxGeom *o2, int flags,
dContactGeom *contact, int skip)
{
dIASSERT (skip >= (int)sizeof(dContactGeom));
dIASSERT (o1->type == dRayClass);
dIASSERT (o2->type == dConeClass);
dxRay *ray = (dxRay*) o1;
dxCone *cone = (dxCone*) o2;
contact->g1 = o1;
contact->g2 = o2;
dVector3 tmp,q,v;
tmp[0] = ray->final_posr->pos[0] - cone->final_posr->pos[0];
tmp[1] = ray->final_posr->pos[1] - cone->final_posr->pos[1];
tmp[2] = ray->final_posr->pos[2] - cone->final_posr->pos[2];
dMULTIPLY1_331 (q,cone->final_posr->R,tmp);
tmp[0] = ray->final_posr->R[0*4+2] * ray->length;
tmp[1] = ray->final_posr->R[1*4+2] * ray->length;
tmp[2] = ray->final_posr->R[2*4+2] * ray->length;
dMULTIPLY1_331 (v,cone->final_posr->R,tmp);
dReal r = cone->radius;
dReal h = cone->lz;
dContactGeom Contact[2];
if (FindIntersectionConeRay(r,h,q,v,Contact))
{
dMULTIPLY0_331(contact->normal,cone->final_posr->R,Contact[0].normal);
dMULTIPLY0_331(contact->pos,cone->final_posr->R,Contact[0].pos);
dOPE(contact->pos,+=,cone->final_posr->pos);
contact->depth = Contact[0].depth * dLENGTH(v);
/*
dMatrix3 RI;
dRSetIdentity (RI);
dVector3 ss;
ss[0] = 0.01f;
ss[1] = 0.01f;
ss[2] = 0.01f;
dsSetColorAlpha (1,0,0,0.8f);
dsDrawBox(contact->pos,RI,ss);
*/
return 1;
}
return 0;
}
int dCollideConeSphere(dxGeom *o1, dxGeom *o2, int flags, dContactGeom *contact, int skip)
{
dIASSERT (skip >= (int)sizeof(dContactGeom));
dIASSERT (o1->type == dConeClass);
dIASSERT (o2->type == dSphereClass);
dxCone *cone = (dxCone*) o1;
dxSphere ASphere(0,cone->radius);
dGeomSetRotation(&ASphere,cone->final_posr->R);
dGeomSetPosition(&ASphere,cone->final_posr->pos[0],cone->final_posr->pos[1],cone->final_posr->pos[2]);
return dCollideSphereSphere(&ASphere, o2, flags, contact, skip);
}
int dCollideConeBox(dxGeom *o1, dxGeom *o2, int flags, dContactGeom *contact, int skip)
{
dIASSERT (skip >= (int)sizeof(dContactGeom));
dIASSERT (o1->type == dConeClass);
dIASSERT (o2->type == dBoxClass);
dxCone *cone = (dxCone*) o1;
dxSphere ASphere(0,cone->radius);
dGeomSetRotation(&ASphere,cone->final_posr->R);
dGeomSetPosition(&ASphere,cone->final_posr->pos[0],cone->final_posr->pos[1],cone->final_posr->pos[2]);
return dCollideSphereBox(&ASphere, o2, flags, contact, skip);
}
int dCollideCCylinderCone(dxGeom *o1, dxGeom *o2, int flags, dContactGeom *contact, int skip)
{
dIASSERT (skip >= (int)sizeof(dContactGeom));
dIASSERT (o1->type == dCCylinderClass);
dIASSERT (o2->type == dConeClass);
dxCone *cone = (dxCone*) o2;
dxSphere ASphere(0,cone->radius);
dGeomSetRotation(&ASphere,cone->final_posr->R);
dGeomSetPosition(&ASphere,cone->final_posr->pos[0],cone->final_posr->pos[1],cone->final_posr->pos[2]);
return dCollideCCylinderSphere(o1, &ASphere, flags, contact, skip);
}
extern int dCollideSTL(dxGeom *o1, dxGeom *o2, int flags, dContactGeom *contact, int skip);
int dCollideTriMeshCone(dxGeom *o1, dxGeom *o2, int flags, dContactGeom *contact, int skip)
{
dIASSERT (skip >= (int)sizeof(dContactGeom));
dIASSERT (o1->type == dTriMeshClass);
dIASSERT (o2->type == dConeClass);
dxCone *cone = (dxCone*) o2;
dxSphere ASphere(0,cone->radius);
dGeomSetRotation(&ASphere,cone->final_posr->R);
dGeomSetPosition(&ASphere,cone->final_posr->pos[0],cone->final_posr->pos[1],cone->final_posr->pos[2]);
return dCollideSTL(o1, &ASphere, flags, contact, skip);
}