00001 #include "ConvexQuad.h"
00002 #include <limits>
00003
00004 using namespace std;
00005
00006 namespace rcrt
00007 {
00008
00009 ConvexQuad::ConvexQuad(SolidObject* parent, const Point3D& pa, const Point3D& pb, const Point3D& pc,const Point3D& pd)
00010 :Primitive(parent), va(pa), vb(pb), vc(pc), vd(pd), center((pa+(pb+(pc+pd)))/4), normal((pb-pa).crossP(pd-pa)), box(pa, pb)
00011 {
00012 box.extend(vc);
00013 box.extend(vd);
00014 normal.normalize();
00015
00016 }
00017
00018 ConvexQuad::~ConvexQuad()
00019 {
00020 }
00021
00022 bool ConvexQuad::isInside(Intersection& is) const
00023 {
00024 char LargestComponent;
00025 Vec2D A, B, C, D;
00026 Vec2D M;
00027 Vec2D AB, BC, CD, AD, AM, AE;
00028 float u, v;
00029 float a, b, c, SqrtDelta;
00030 bool isValid = false;
00031 Vec2D Vector;
00032
00033
00034
00035
00036
00037 Vec3D absNormal = absNormal.abs();
00038
00039 LargestComponent = (fabs(normal.x()) > fabs(normal.y()) ? (fabs(normal.x()) > fabs(normal.z()) ? 'x' : 'z') : (fabs(normal.y()) > fabs(normal.z()) ? 'y' : 'z'));
00040
00041 if (LargestComponent == 'x') {
00042 A[0] = va[1]; B[0] = vb[1]; C[0] = vc[1]; D[0] = vd[1];
00043 M[0] = is.getPosition()[1];
00044 }
00045 else {
00046 A[0] = va[0]; B[0] = vb[0]; C[0] = vc[0]; D[0] = vd[0];
00047 M[0] = is.getPosition()[0];
00048 }
00049
00050 if (LargestComponent == 'z') {
00051 A[1] = va[1]; B[1] = vb[1]; C[1] = vc[1]; D[1] = vd[1];
00052 M[1] = is.getPosition()[1];
00053 }
00054 else {
00055 A[1] = va[2]; B[1] = vb[2]; C[1] = vc[2]; D[1] = vd[2];
00056 M[1] = is.getPosition()[2];
00057 }
00058
00059 AB = B - A; BC = C - B;
00060 CD = D - C; AD = D - A;
00061 AE = CD + AB; AE = AE * -1; AM = M - A;
00062
00063 if (fabs(AB.det(CD)) < numeric_limits<float>::epsilon())
00064 {
00065 Vector = AB - CD;
00066 v = AM.det(Vector) / AD.det(Vector);
00067 if ((v >= 0.0) && (v <= 1.0)) {
00068 b = AB.det(AD) - AM.det(AE);
00069 c = AM.det(AD);
00070 u = fabs(b) < numeric_limits<float>::epsilon() ? -1.0 : c/b;
00071 isValid = ((u >= 0.0) && (u <= 1.0));
00072 }
00073 }
00074 else if (fabs(BC.det(AD)) < numeric_limits<float>::epsilon())
00075 {
00076 Vector = AD + BC;
00077 u = AM.det(Vector) / AB.det(Vector);
00078 if ((u >= 0.0) && (u <= 1.0)) {
00079 b = AD.det(AB) - AM.det(AE);
00080 c = AM.det(AB);
00081 v = fabs(b) < numeric_limits<float>::epsilon() ? -1.0 : c/b;
00082 isValid = ((v >= 0.0) && (v <= 1.0));
00083 }
00084 }
00085 else
00086 {
00087 a = AB.det(AE); c = - AM.det(AD);
00088 b = AB.det(AD) - AM.det(AE);
00089 a = -0.5/a; b *= a; c *= (a + a); SqrtDelta = b*b + c;
00090 if (SqrtDelta >= 0.0) {
00091 SqrtDelta = sqrt(SqrtDelta);
00092 u = b - SqrtDelta;
00093 if ((u < 0.0) || (u > 1.0))
00094 u = b + SqrtDelta;
00095 if ((u >= 0.0) && (u <= 1.0)) {
00096 v = AD[0] + u * AE[0];
00097 if (fabs(v) < numeric_limits<float>::epsilon())
00098 v = (AM[1] - u * AB[1]) / (AD[1] + u * AE[1]);
00099 else
00100 v = (AM[0] - u * AB[0]) / v;
00101 isValid = ((v >= 0.0) && (v <= 1.0));
00102 }
00103 }
00104 }
00105
00106 if (isValid) {
00107 is.setUV(Point2D(u,v));
00108 }
00109 return (isValid);
00110 }
00111
00112 const AABB& ConvexQuad::getBoundingBox() const
00113 {
00114 return box;
00115 }
00116
00117 const Point3D& ConvexQuad::getCentroid() const
00118 {
00119 return center;
00120 }
00121
00122 Intersection ConvexQuad::intersect(Ray& r) const
00123 {
00124
00125 float dist = r.dir() * normal;
00126 if (fabs(dist) < numeric_limits<float>::epsilon())
00127 return Intersection();
00128
00129
00130 Vec3D v = va - r.org();
00131 dist = (v*normal) / dist;
00132
00133 Intersection is(dist, this, r.atDistance(dist), 0, 0);
00134
00135 if (isInside(is)) {
00136 is.setSNormalL(normal.abs());
00137 return is;
00138 }
00139 return Intersection();
00140 }
00141
00142 Vec3D ConvexQuad::getSNormal(float a, float b) const
00143 {
00144 return normal;
00145 }
00146
00147 Vec3D ConvexQuad::getGNormal(float a, float b) const
00148 {
00149 return normal;
00150 }
00151
00152 Point2D ConvexQuad::getUV(float a, float b) const
00153 {
00154
00155 return Point2D(a,b);
00156 }
00157 }