Book Corrections for Geometric Tools for Computer Graphics 
Book Corrections Organized by Date of Change 
30 August 2017, page 771.
In the SegmentInCode listing, the case for the "vertex is convex"
has the inequalities reversed. The line of code should be
return Kross(diff, edgeR) < 0 and Kross(diff, edgeL) > 0 
02 July 2016, page 941.
All the product formulas are incorrect. They should be
sin(a)*sin(b) = (cos(ab)  cos(a+b))/2 sin(a)*cos(b) = (sin(a+b) + sin(ab))/2 cos(a)*cos(b) = (cos(ab) + cos(a+b))/2The last formula in the book appears to be an attempt to cover all the cases, but the first three equations are all the cases. 
07 July 2014, pages 157158. For simplicity of notation here, define N = hat(n), N^T is the transpose of N, and V = hat(v). The last row, first column entry of the general shearing matrix should be (tan(theta))*Dot(N,Q)*V. The first line on page 158 makes an appeal to applying the point transformation T to a vector N, which makes no sense mathematically. Even if it did, it is not clear why the author claims T causes N to change at all. Let P = P1 + h*N, where h = Dot(N,PQ). The shearing is P" = P + t*V, where from trigonometry, t = h*tan(theta), so P" = P + (tan(theta))*Dot(N,PQ)*V and P"Q = [I + (tan(theta))*V*N^T](PQ). Add Q to both sides and factor to P" = [I + (tan(theta))*V*N^T]*P  (tan(theta))*V*N^T*Q = [I + (tan(theta))*V*N^T]*P  (tan(theta))*Dot(N,Q)*V. 
14 December 2013, pages 208210. The algorithm uses barycentric coordinates to determine closest edge, but this does not work for a triangle with an obtuse angle. Thus, the SquaredDistance pseudocode is incorrect. The online source code does not use this algorithm; it uses the previously discussed (correct) algorithm. 
07 October 2013, page 229. Equation (6.24) has 5 cases (top to bottom). The final inequality in case 2 should be Dot(d0,delta+T0*d0) ≤ 0. The final inequality in case 4 should be Dot(d0,delta+T0*d0T1*d1) ≤ 0. 
21 September 2013, page 15. The "Matrix" row of Table 2.1 has pseudocode "Matrix3x3 m, m;" that should be "Matrix3x3 m, n;". 
01 March 2013, page 225. The denominator of Equation (6.17) should be dot(perp(d0),d1) rather than dot(perp(d1),d0); this has the effect of negating the righthand side of the equation. 
03 November 2012, page 78.
Figure 3.19 is confusing and has one typographical error. It is an attempt to
describe the linear transformation T(x) = 2*x. In terms of basis vectors u
and v: T(u) = 2*u and T(v) = 2*v. The vector s in the figure should be more
clearly written:
s = T(2*v) = 2*T(v) // because T is linear = 2*(2*v) // by definition of TCall the top vector of the figure r, in which case r = T(3*u) = 3*T(u) // because T is linear = 3*(2*u) // by definition of TThus, the current figure should have "3T(u) = 2(3u)". 
03 November 2012, page 79.
Figure 3.20 is confusing and has two typographical errors. It is an attempt to
describe the linear transformation that is defined on the basis vectors u and v:
T(u) = 2*u and T(v) = 1.5*v. Given a vector x = a*u+b*v, the transformation is
T(x) = T(a*u+b*v) = a*T(u)+b*T(v) = a*(2*u)+b*(1.5*v). The vector s in the
figure is incorrect and should be
s = T(2*v) = 2*T(v) // because T is linear = 2*(1.5*v) // by definition of TCall the top vector of the figure r, in which case r = T(3*u) = 3*T(u) // because T is linear = 3*(2*u) // by definition of TThus, the current figure should have "3T(u) = 2(3u)". 
10 September 2012, page 511.
Lines 4 and 5 are missing right parentheses. The pseudocode should be
t[0] = (b + sqrt(discrm))/(2*a); t[1] = (b  sqrt(discrm))/(2*a); 
10 September 2012, page 518. The page has 7 references to array elements "[3]". These should be "[2]". 
10 September 2012, page 511.
Lines 4 and 5 are missing right parentheses. The pseudocode should be
t[0] = (b + sqrt(discrm))/(2*a); t[1] = (b  sqrt(discrm))/(2*a); 
08 July 2012, page 548. In matrix M, entry 6 of row 2 says "2(u_{x} v_{y} + v_{x} v_{y})". This should be "2(u_{x} v_{y} + v_{x} u_{y})". 
24 June 2012, page 83. The sentence starting with "We've just done two forbidden operations on pointsscaling by a vector and directly adding them." This should say "We've just done two forbidden operations on pointsmultiplying points by scalars and directly adding them." 
24 June 2012, page 84. The first displayed equation has a term "alpha3*(P3P2)". This term should be "alpha3*(P3P1)". 
20 February 2010, pages 564569.
The first paragraph of Intersection Detection mentions a base point B
and refers you to Figure 11.14, which shows the base point labeled as C.
Later in Intersection Detection, there is a section Finite Cone
that refers to the base point using name B. This point is as shown in Figure
11.14 and is not the same as the apex A. On page 567, a vector hat(w) is defined. The author used the hat notation to denote unitlength vectors. However, the vector Cross(hat(d),Cross(hat(n),hat(d))) is not necessarily unit length. For example, if hat(n) is nearly parallel to hat(d), Cross(hat(n),hat(d)) is nearly the zero vector, in which case hat(w) is nearly the zero vector. The correct equation is hat(w) = Cross(hat(d),Cross(hat(n),hat(d)))/Length(Cross(hat(n),hat(d)))On pages 567, 568, and 569, the author has occurrences of arrow(w), but these should be hat(w). The variable a is the distance between points B and Ia, so the last equation on page 567 should be a = Length(Ia  B)That is, the subtraction of h is incorrect. Just to be clear, the variable b is the distance between points B and Ic, and the variable c is the distance between points Ia and Ic. On pages 568 and 569, any occurrence of Length(IaB)h should be replaced by Length(IaB). Observe that with the correct hat(w) equation, it must be that Dot(hat(n),hat(w)) = Length(Cross(hat(n),hat(d)), so it is not necessary to explicitly compute hat(w) in the source code. The equation for b^{2} is b^{2} = a^{2} [1/Length(Cross(hat(n),hat(d)))^{2}  1]However, notice that this equation has numerical problems when hat(n) and hat(d) are nearly parallel. The length of the cross product is nearly zero, so computing the reciprocal length is a problem. Instead, you want to code up the test as v = Cross(hat(n),hat(d)); lensqr = Dot(v,v); a*a*(1  lensqr) ≤ lensqr*r*r;Here is an alternative for the intersection test of plane and finite cone. Let the plane have unitlength normal N and contain point P. The plane equation is Dot(N,XP) = 0. The idea is to compute the interval of projection of the finite cone onto the line containing P with direction N. This inteval has endpoints P + t0*N and P + t1*N, where t0 < t1. If the interval [t0,t1] contains 0, then the cone and plane intersect. The projection interval endpoints are determined by one of the following.
X = C + h*D + h*tan(A)*(cos(B)*U + sin(B)*V)for B in [0,2*pi). The vectors U and V are unit length and {U,V,D} is an orthonormal set (vectors are unit length and mutually perpendicular). The projection of X onto the normal line is P + t(B)*N, where t(B) = Dot(N,XP) = Dot(N,CP) + h*Dot(N,D) + h*tan(A)*(cos(B)*Dot(N,U) + sin(B)*Dot(N,V))I have indicated that t is a function of B (different circle points project to different normal line points). The projection of C onto the normal line is Dot(N,CP). The minimum value of t(B) is t0 and the maximum value of t(B) is t1. Define d = Dot(N,CP) and r = sqrt(Dot(N,U)^2 + Dot(N,V)^2) = sqrt(1  Dot(N,D)^2)t0 is the minimum, and t0 is the maximum, of the set {d, d + h*Dot(N,D)  h*tan(A)*r, d + h*Dot(N,D) + h*tan(A)*r}The pseudocode is bool ConeIntersectsPlane (Cone cone, Plane plane) { Vector CmP = cone.C  plane.P; float d = Dot(plane.N, CmP); float NdotD = Dot(plane.N, cone.D); float r = sqrt(1  NdotD*NdotD); float tmp0 = d + cone.h*NdotD; float tmp1 = cone.h*cone.tanAngle*r; float vm = tmp0  tmp1; float vp = tmp0 + tmp1; float t0 = min(d, min(vm, vp)); float t1 = max(d, max(vm, vp)); return t0 <= 0 && 0 <= t1; }This algorithm does not require you to explicitly determine which of the bookmentioned cases 1, 2, 3a, or 3b applies. 
03 June 2009, page 454.
The pseudocode after the "Closest point on the box" comment should be
qPrime.x = box.extents.x; qPrime.y = line.origin.y; qPrime.z = line.origin.z; 
29 January 2009, page 43. Bullet 'v' in Section 2.6.3 is a repeat of bullet 'ii' and should be removed. 
29 January 2009, page 44.
The example combination of vectors in a displayed equation is
2v_{1} + 3v_{3} +1v_{3}but should be 2v_{1} + 3v_{2} +1v_{3}That is, the subscript on the middle term needed to be fixed. 
29 January 2009, page 46. The sentence: "We can map real numbers to real numbers: ..." should say "We can map nonnegative real numbers fo nonnegative real numbers: ..." 
29 January 2009, page 53. The last sentence should say "... into Equations 2.4 and 2.5, ...". 
29 January 2009, page 777. The first sentence of the last paragraph has a subjectverb mismatch. It should be "An edge s_{k} has ...". In the same paragraph, the end of the second line has "s_{1}.y" but should be "s_{1}.y_{0}". 
06 December 2008, page 872.
The first displayed equation is
max{a_0, 1+a_1, ..., 1+a_{n1}} = 1 + max{a_0, ..., a_{n1}}but the last equaltity should be an inequality max{a_0, 1+a_1, ..., 1+a_{n1}} ≤ 1 + max{a_0, ..., a_{n1}} 
14 November 2008, page 428. The paragraph about the case of region 2 has "the gradient at the corner (1,1) cannot point inside the unit square". This should say "outside" rather than "inside". 
17 August 2008, page 562.
The pseudocode has a few issues. The first line of code that computes d
should be a signed distance calculation. The first comparison of
cos(theta) should have the abs about the cosine term, not
about the epsilon term. If an epsilon comparison is made for when
cos(theta) is nearly zero, to be consistent an epsilon comparison
should also be made for when cos(theta) is nearly one (the circle
case). Here is my own rewrite of the pseudocode without epsilon tests.
The following data structures are assumed.
struct Plane { Point origin; Vector normal; }; struct Cylinder { Point origin; Vector direction; float radius; }; struct Line { Point origin; Vector direction; }; struct Circle { Point center; Vector normal; float radius; }; struct Ellipse { Point center; Vector normal, major, minor; float majorLength, minorLength; }; enum Type { EMPTY_SET, ONE_LINE, TWO_LINES, CIRCLE, ELLIPSE }; struct Intersection { Type type; Line line; // valid when type = ONE_LINE Line line1, line2; // valid when type = TWO_LINES Circle circle; // valid when type = CIRCLE Ellipse ellipse; // valid when type = ELLIPSE };The plane normal, cylinder direction, line direction, circle normal, ellipse normal, ellipse major direction, and ellipse minor direction are all required to be unitlength vectors. The output is intr. The returned bool value is true when there is an intersection, in which case intr is valid (of the type specified in its Type member). The returned value is false when there is no intersection, in which case the Type member has value EMPTY_SET. bool PlaneCylinderIntersection (Plane plane, Cylinder cylinder, Intersection& intr) { float d = Dot(plane.normal, cylinder.origin  plane.origin); Point projectedCylinderOrigin = cylinder.origin  d * plane.normal; float cosTheta = Dot(cylinder.direction, plane.normal); if (abs(cosTheta) > 0) // cylinder axis intersects plane in a unique point { if (abs(cosTheta) < 1) // intersection is an ellipse { intr.type = ELLIPSE; intr.ellipse.normal = plane.normal; intr.ellipse.center = projectedCylinderOrigin  (d/cosTheta) * cylinder.direction; intr.ellipse.major = cylinder.direction  cosTheta * plane.normal; intr.ellipse.minor = Cross(plane.normal, intr.ellipse.major); intr.ellipse.majorLength = cylinder.radius / abs(cosTheta); intr.ellipse.minorLength = cylinder.radius; Normalize(intr.ellipse.major); Normalize(intr.ellipse.minor); return true; } else // abs(cosTheta) == 1, intersection is a circle { intr.type = CIRCLE; intr.circle.normal = plane.normal; intr.circle.center = projectedCylinderOrigin; intr.circle.radius = cylinder.radius; return true; } } else // abs(cosTheta) == 0, cylinder is parallel to plane { if (abs(d) < cylinder.radius) // cylinder intersects plane in two lines { Vector offset = Cross(cylinder.direction, plane.normal); float e = sqrt(cylinder.radius * cylinder.radius  d * d); intr.type = TWO_LINES; intr.line1.origin = projectedCylinderOrigin  e * offset; intr.line1.direction = cylinder.direction; intr.line2.origin = projectedCylinderOrigin + e * offset; intr.line2.direction = cylinder.direction; return true; } else if (abs(d) == cylinder.radius) // cylinder intersects plane in one line { intr.type = ONE_LINE; intr.line.origin = projectedCylinderOrigin; intr.line.direction = cylinder.direction; return true; } else // abs(d) > cylinder.radius, cylinder does not intersect plane { intr.type = EMPTY_SET; return false; } } } 
17 August 2008, page 550.
The pseudocode computes b as an absolute value, but the absolute value should
only be used in the comparison to the sphere radius on the line after
computing b. The pseudocode on this page has had other corrections, so the
revision is listed next. I do not like the book's formulationthe revision
is more like how I write the actual code. This is the mathematically correct
formulation, without concern for "epsilon" arithmetic to deal with numerical
roundoff errors. The plane normal is assumed to be a unitlength vector.
The Circle3 data structure is for a circle in three dimensions.
bool PlaneSphereIntersection (Plane plane, Sphere sphere, Circle3& circle) { float b = dotProd(plane.normal, sphere.center  plane.pointOnPlane); float absB = fabs(b); circle.center = sphere.center  b * plane.normal; circle.normal = plane.normal; if (absB <= sphere.radius) { // The sphere intersects the plane in a circle. The circle is // degenerate when fabs(b) is equal to sphere.radius, in which // case the circle radius is zero. circle.radius = sqrt(sphere.radius * sphere.radius  b * b); return true; } return false; }The circle parameter is an output of the function. It is valid only when the function's return value is true. The plane's normal vector is copied to the circle's normal so that you can actually reconstruct the circle. If the circle center is K, the normal is N, and the radius is R, then choose vectors U and V so that {U,V,N} is an orthonormal set (vectors are unit length and mutually perpendicular). A circle parameterization is X(t) = K + R*(cos(t)*U + sin(t)*V), 0 <= t < 2*piA numerically robust way to compute U and V is if (fabs(N.x) >= fabs(N.y)) { // N.x or N.z is the largest magnitude component, swap them. invLength = 1/sqrt(N.x * N.x + N.z * N.z); U = Vector(N.z * invLength, 0, N.x * invLength); V = Vector(N.y * U.z, N.z * U.x  N.x * U.z, N.y * U.x); // = Cross(N,U) } else { // N.y or N.z is the largest magnitude component, swap them. invLength = 1/sqrt(N.y * N.y + N.z * N.z); U = Vector(0, N.z * invLength, N.y * invLength); V = Vector(N.y * U.z  N.z * U.y, N.x * U.z, N.x * U.y); // = Cross(N,U) }If you want to allow for numerical roundoff error, say, by requiring the intersection to be a single point when b is within epsilon of the radius, use bool PlaneSphereIntersection (Plane plane, Sphere sphere, Circle3& circle) { float b = dotProd(plane.normal, sphere.center  plane.pointOnPlane); float absB = fabs(b); circle.center = sphere.center  b * plane.normal; circle.normal = plane.normal; if (absB <= sphere.radius  epsilon) { // The sphere intersects the plane in a circle. The circle is // degenerate when fabs(b) is equal to sphere.radius, in which // case the circle radius is zero. circle.radius = sqrt(sphere.radius * sphere.radius  b * b); return true; } else if (absB <= sphere.radius + epsilon) { // The sphere is _nearly_ tangent to the plane. Choose the // intersection to be the point of tangency. circle.radius = 0; return true; } return false; }The magnitude of the nonnegative number epsilon is up to the programmer. Observe that if you choose epsilon to be zero, you effectively have the pseudocode for the mathematically correct case. 
14 August 2008, page 571. The end of the last sentence of the first paragraph of the subsection Nondegenerate PlaneCone Intersections should be: cos(theta) = Dot(hat(a),hat(n)). 
14 August 2008, page 575.
The style of the pseudocode differs from other pages in this section. In the
first block starting with ``float d = '', remove the 'float' types (to be
consistent). In the section about the parabola, consistent notation would be
parabola.axis = Normalize(cone.axis  cosTheta * plane.normal); parabola.perpendicularAxis = Cross(parabola.axis, plane.normal); parabola.vertex = cone.vertex + d* plane.normal + e*parabola.axis; parabola.focalLength = d/2 * cosTheta;Also on this page, the first paragraph has a condition on the dot product of the plane normal hat(n) and the cone axis hat(a) in order that the two vectors be parallel. The condition should be: Dot(hat(n),hat(a) = 1 (for the true mathematical classification) or Dot(hat(n),hat(a)) >= 1epsilon (if you want to worry about numerical roundoff errors). 
14 August 2008, page 580.
The pseudocode has a line that computes the minor axis. It should be
minorAxis = Cross(plane.normal, majorAxis);because nowhere is there defined ``ellipse.u''. 
28 October 2007, page 888.
The last displayed equation is
E(hat(c)) = (sum_{i=0}^n hat(c) * vector(v))^2 = hat(c)^T M hat(c)which is the square of a sum, but it should be a sum of squares, E(hat(c)) = sum_{i=0}^n (hat(c) * vector(v))^2 = hat(c)^T M hat(c) 
28 October 2007, page 889.
Line 5 of the first paragraph of section A.7.8 has
M = sum_{i=0}^2 vector(v) vector(v)^TThe upper limit should be n, M = sum_{i=0}^n vector(v) vector(v)^T 
22 September 2007, page 550.
The last line of pseudocode on this page should be
float radius = sqrt(sphere.radius^2  b^2);The pseudocode also rejects the "just touching" case as "not an intersection"; that is, "b < sphere.radius" should be "b ≤ sphere.radius". In fact, since the pseudocode allows a small amount of error ("radius < epsilon"), it might even be good to support an error term in the first test: "if (b ≤ sphere.radius + epsilon)". 
05 February 2007, page 228.
The symbols in Equation (6.23) on page 228 are not typeset correctly. In that equation,
let us say that the first row of the equation is region R0 and the last row of the
equation is R8 (the parameter plane is partitioned into 9 regions). The corrections are
R2: 0 < tilde(t0) < T0 and bar(t1) >= T1 R4: 0 < tilde(t1) < T1 and bar(t0) >= T0 R6: hat(t0) >= T0 and tilde(t1) <= 0 R7: tilde(t0) <= 0 and hat(t1) >= T1 R8: tilde(t0) >= T0 and tilde(t1) >= T1 
28 January 2007, page 633. The pseudocode for the rayOBB intersection is absurd. The first line of code is missing two right parentheses after "box.axis[i]". The inequality "r + box.halfLength[i] > 0" should be "r + box.halfLength[i] < 0". The rest of the pseudocode is worthless. The value s can very well be zero, leading to a division by zero when computing t0 and t1. If you find worthless pseudocode in the book, first let me know and, second, look at the code at my website. I wrote the real code (not the pseudocode), so if you find errors in the real code, please let me know. 
23 January 2007, page 612.
The last line of pseudocode is
if (value < min) min = value; else max = value;This should be if (value < min) min = value; else if (value > max) max = value; 
31 October 2006, page 550. To be consistent with the other pseudocode, "if (b < r)" should be written as "if (b < sphere.radius)". 
28 October 2006, page 548. Equation (11.20) has "r" when it should be "r^2" (r squared). 
08 September 2006, page 384.
The bookcorrection entry for page 635 (24 February 2004) was
nx*(x1+x0)+ny*(y1+y0)+nz*(z1+z0)+2*d <= nx*(x1x0)+ny*(y1y0)+nx*(z1z0)but the lefthand side is missing absolute values. It should be nx*(x1+x0)+ny*(y1+y0)+nz*(z1+z0)+2*d <= nx*(x1x0)+ny*(y1y0)+nx*(z1z0)I have modified the older entries here. 
06 September 2006, page 384.
The first displayed equation has a Q that should be a P:
Q' = P + s*e0 + t*e1 
30 March 2006, page 387. The top line of pseudocode has "Dot(p,n)d", but should be "Dot(p,n)+d". 
28 January 2006, pages 417426.
These pages cover distance calculations between linear components. The
pseudocode has major problems. Apparently, one of these blocks of code was
written and then cutpastemodified for the other blocks. They all exhibit
the same ills. The difference between closest points needs to be of the
form
v = (base0 + s*direction0)  (base1 + t*direction1)That is, the book parentheses are incorrect. Worse is the handling of the distinction between parallel and nonparallel components. If "det" is ever zero, the pseudocode will do a division by it. Also, the epsilon test is problematic when the segment direction vectors have large length. If D0 and D1 are the line segment lengths, each computed as a difference of the line segment end points, and if A is the angle between the directions, it is the case that det = Cross(D0,D1)^2 = D0^2 * D1^2 * (sin(A))^2If the angle A is very small, but the lengths of the directions are very large, "det" can be quite large, causing you to classify the segments are nonparallel when they are nearly parallel. My own implementation uses a line segment represention of C+s*D, where C is the midpoint of the segment, D is a unitlength direction, and s <= e where e is half the length of the segment. This leads to a test det = (sin(A))^2so that the epsilon test is a relative error measurement, not an absolute error measurement. Section 10.8.3 of the book describes how I think of distance measurements for linear components, and I much prefer this approach based on actually writing code and seeing how it performs on real data. 
28 January 2006, page 686. The pseudocode has two sets named CoPos and CoNeg. These should be CoSame and CoDiff, respectively. 
13 January 2006, pages 767768. The visibility graph does not contain an arc (i,i), where i is an index of a polygon vertex. The text should also state that adjacent polygon vertices are not considered to be visible to each other. Therefore, the visibility graph does not contain arcs (i,(i+1) modulo n), where n is the number of polygon vertices and 0 <= i <= n1. On page 768, the text states that the adjacency matrix for a convex polygon has 3 zero diagonals, the rest of the entries are ones. However, the entries at (row,column) of (0,n1) and (n1,0) are also zero since these correspond to adjacent polygon vertices. 
13 December 2005, page 120. In the middle of the page, I do not understand why the statement "with, as before, Vector(u) retaining its usual matrix representation" is made. This seems to imply that Vector(u) will be represented as a row matrix as shown at the top of the page. When Cross(Vector(v),Vector(u)) is represented using a skewsymmetric matrix for Vector(v), Vector(u) must be written as a column matrix. Moreover, the skewsymmetric matrix Tilde(v) is the negative of what it should be. Listing the rows from top to bottom, it should be Tilde(v) = ((0,v3,+v2),(+v3,0,v1),(v2,+v1,0)). 
22 November 2005, page 419.
The pseudocode for c should be
c = Dot(ray.direction, ray.direction); 
22 November 2005, page 531.
The pseudocode for b at the bottom of the page should be
b = (s1 * n1n2dot  s2 * n1normsqr) / (n1n2dot^2  n1normsqr * n2normsqr);That is, the first occurrence of n2normsqr in the book should have been n1normsqr. 
07 November 2005, page 356.
The first equation for the torus is incorrect. It should be the following
equation. I have factored it in a manner that leads to more efficient
numerical computation:
(x^2 + y^2 + z^2)^2  2*(r0^2 + r1^2)*(x^2 + y^2 + z^2) + 4*r0^2*z^2 + (r0^2  r1^2)^2 = 0My torus class (in Wm3Torus files) has the correct formula listed. 
07 November 2005, page 660.
The equation for c0 is missing a power 2. It should be
c0 = (Dot(P,P)  (R^2+r^2))^2  4*R^2*(r^2  P_z^2) 
02 October 2005, pages 628629.
The pseudocode for rayAABB intersection has a loop that illustrates clipping
only against the planes x=xmin and x=xmax. This makes it unclear how the
the full set of tests fit in with the rest of the loop logic. The following
pseudocode is the replacement. The Point and Vector data structures are
assumed to have operator[] functions to access the coordinates. Index 0
corresponds to x (P[0] and P.x are equivalent), index 1 corresponds to y,
and index 2 corresponds to z. It is my opinion that this function should
also return the interval of intersection (first point of intersection
and last point of intersection). This change shows up in the last input to
the function. It is guaranteed that tIntersect[0] <= tIntersect[1]. If the
ray just grazes the AABB, we will return the same value in both parameters.
boolean RayIntersectAABB (Point P, Vector d, AABB box, float tIntersect[2]) { // The first part of this code tests for lineAABB intersection. tNear = INFINITY; tFar = INFINITY; for (i = 0; i < 3; i++) { // Check if the line is parallel to the min[i] and max[i] planes. if (abs(d[i]) < epsilon) { // The line is parallel to the min[i] and max[i] planes. if (P[i] < box.min[i]  P[i] > box.max[i]) { // The line is outside the slab defined by min[i] and max[i]. return false; } // else: The line is inside the slab defined by min[i] and max[i]. // Just continue testing with the next slab. continue; } // The line is not parallel to the slab planes, so compute the // parameters of intersection. The line segment of intersection is // P+t*d, where t0 <= t <= t1. t0 = (box.min[i]  P[i])/d[i]; t1 = (box.max[i]  P[i])/d[i]; if (t0 > t1) { tmp = t1; t1 = t0; t0 = tmp; } // Compare with current values. The current slab may have increased // tNear and/or decreased tFar. if (t0 > tNear) { tNear = t0; } if (t1 < tFar) { tFar = t1; } // Check if the line misses the AABB entirely. if (tNear > tFar) { return false; } } // At this place in the code, the line intersects the AABB in the line // segment P+t*d, tNear <= t <= tFar. If all you care about is lineAABB // intersection, then use the following three lines of code and omit the // remaining pseudocode. // // tIntersect[0] = tNear; // tIntersect[1] = tFar; // return true; // The ray is (P+t*d, t >= 0). We need to compute the intersection of // the interval tNear <= t <= tFar with the interval t >= 0. if (tFar < 0) { return false; } if (tNear < 0) { tNear = 0; } tIntersect[0] = tNear; tIntersect[1] = tFar; return true; } 
28 June 2005, page 516. The last line of pseudocode has Matrix4x4 when it should be Matrix3x3. 
17 June 2005, page 861.
A bunch of errors are in "Quaternion to Matrix" and "Matrix to Quaternion".
This material was from "3D Game Engine Design", where the same errors occur,
but I had posted the correction for that book in February 2001. If R = [rij],
then equation (A.8) should have
r00 = 12y^22z^2 r01 = 2xy2wz r02 = 2xz+2wy r10 = 2xy+2wz r11 = 12x^22z^2 r12 = 2yz2wx r20 = 2xz2wy r21 = 2yz+2wx r22 = 12x^22y^2The last sentence of the first paragraph of "Matrix to Quaternion" has y = (r20r02)/(4w)but it should be y = (r02r20)/(4w)In the second paragraph of that section. When r00 is the maximum, the book has w = (r12r21)/(4x)but it should be w = (r21r12)/(4x)When r11 is the maximum, the book has w = (r20r02)/(4y)but it should be w = (r02r20)/(4y)When r22 is the maximum, the book has w = (r01r10)/(4z)but it should be w = (r10r01)/(4z) 
16 May 2005, page 377. The equation at the bottom of the page for bar(t) has the wrong numerator. The equation should be "bar(t) = (b*da*e)/(a*cb*b)". The pseudocode on page 380 does use the correct equation. 
15 December 2004, pages 844 and 846. The second line from the bottom on page 844 should have "e2 = 98". The fifth line from the bottom of page 846 should have "(1.11596,1.14945,1.63258)". 
14 December 2004, page 437.
The displayed matrix equation is incorrect. Setting the gradient to zero
leads to the system of equations
a00*u + a01*v + a02*t + b0 = 0 a01*u + a11*v + a12*t + b1 = 0 a02*u + a12*v + a22*t + b2 = 0The pseudocode on that page essentially solves this system using Cramer's rule (computes cofactors and determinant). 
13 October 2004, page 534.
The equations for Cramer's rule have errors. The inverse determinant should
be
invDet = 1/(a0*BC  b0*AC + c0*AB)The Z value should be Z = (b0*AD + a0*DB  d0*AB)*invDetI do not have a copy of Bowyer and Woodwark, the reference from which the equations are supposedly taken, so I do not know if these errors occur in that reference or were introduced in transcribing. 
14 August 2004, page 92.
Item "iv. Normal:" should be
v_perp = v  Dot(u,v)*u/Dot(u,u)The book has u as the first vector in the subtraction. 
13 July 2004, page 367. In the case t0 > 0, the distance from a point Q to a ray P+t*D should be Q(P+t0*D). The book has Q(Q+t0*D), a typographical error. 
14 June 2004, page 290.
The last line of code is missing a left parenthesis. It should be
v0.y = ((r2 * u.y) + rad) / denom; 
14 June 2004, page 291.
The first line of pseudocode with solution[0] has v.x when it should be v0.x:
solution[0].setDir(v0.y, v0.x); 
19 May 2004, page 832.
In the elimination step, the line of code
a_{k,i} = a_{k,i} + m * a_{k,j}should be a_{k,i} = a_{k,i} + m * a_{j,i} 
18 May 2004, page 645.
The upper bounds on the displayed inequalities are wrong. They should be
P_{min,x} <= Q_x <= P_{max,x} P_{min,y} <= Q_y <= P_{max,y} P_{min,z} <= Q_z <= P_{max,z}The pseudocode on that same page is inconsistent about its reference to the sphere center. I changed "sphere.c" to "sphere.center" in two places. 
12 May 2004, page 366.
This entry had been listed in the corrections organized by page number,
but it was not duplicated in the section organized by date of change. The pseudocode has a line t = Dot(l.direction, VectorSubtract(q, l.direction));This should be t = Dot(l.direction, VectorSubtract(q, l.origin)); 
09 April 2004, page 187. The displayed equation for the recursion formulas has "B_{i,0}(t)" but should be "B_{i,1}(t)". The displayed equation after that one is the recursion, but should be for "2 <= j <= n" (the text has 1, not 2). 
08 March 2004, page 328. The left image of Figure 9.2 has "r = d" and the right image has "r = d". To be consistent with the paragraph on this page "More significantly...", these should be reversed. The origin in the left image has a negative signed distance from the plane (it is "behind" the plane). The origin in the right image has a positive signed distance from the plane (it is "in front of" the plane). 
08 March 2004, page 490.
The discussion in this section is about linear components (lines, rays, or
line segments) and polygons. The text sometimes talks about lines, sometimes
about rays. The pseudocode has a test "if (t < 0 ) { return false; }" which
indicates it is about raypolygon intersection, but is incorrectly named
LinePolygonIntersection. Also, there is a line
p = line.origin + t * r.d;which should be p = line.origin + t * line.direction;another indication that the author had a ray in mind. Near the end of the pseudocode is intersection = line.origin + t * line.direction;However, the intersection p was calculated earlier, so just use intersection = p; 
08 March 2004, page 493.
The pseudocode has the declaration
Point3D intersectionThis is not needed. 
08 March 2004, page 510.
The solution for the quadratic formula is missing a minus sign on the "b"
term. It should be
t = (b +/ sqrt(b^24*a*c))/(2*a)The pseudocode has a = (tline.direction.x + tline.direction.y)^2;This should be a = tline.direction.x^2 + tline.direction.y^2; 
24 February 2004, page 635.
The return statement in the pseudocode is incorrect. It appears to be
designed for culling an AABB against the plane rather than indicating if the
AABB is intersected by the plane. The return statement also indicates that
the plane and AABB are considered not to be intersecting if the box
is just touching the plane at a single vertex or single edge or single face.
If the plane intersects the box at any point, I suggest replacing the return
statement by
return Dot(plane.normal,dMin)+plane.d <= 0 && Dot(plane.normal,dMax)+plane.d >= 0;I much prefer the "test intersection" that I use for planeOBB, which specializes easily when you have an AABB. No branching, just a few arithmetic calculations. Let N = (nx,ny,nz) be the plane normal and d be the plane constant. The plane is Dot(N,X)+d = 0. Let the AABB minima be x0, y0, and z0. Let the maxima be x1, y1, and z1. The test for intersection is nx*(x1+x0)+ny*(y1+y0)+nz*(z1+z0)+2*d <= nx*(x1x0)+ny*(y1y0)+nz*(z1z0) 
30 January 2004, page 59.
The third displayed equation should be
D^2 = (b + mx1  y1)^2 + (b + mx2  y2)^2 + (b + mx3  y3)^2The book has yn in the last term of the expression. 
16 January 2004, page 769. The sentence that starts with "If < V_{i1}, V_{i+2} > is a diagonal" should have V_{i+1}, not V_{i+2}. 
16 January 2004, page 758759.
The input N is number of points, but the first mesh insert call returns a
new triangle also named N. These are, of course, different things. For
clarity, use Nj as the return value of the two "mesh.Insert" calls, one on
page 758 and one on page 759. Use Nj as the argument to the two
"stack.Push" calls, one on page 758 and one on page 759.
In the bottom block of pseudocode on page 759, the figure reference should be 13.46, not 13.45. After that line the pseudocode should be compute i0, i1, i2, i3 wwith T.v(i1) = A.v(i1) and T.v(i2) = A.v(i2);In the "if" clause, two lines must be added to the end of the clause mesh.Remove(T); mesh.Remove(A);After the two new triangles N0 and N1 are inserted into the mesh, the old triangles T and A must be removed. 
14 January 2004, page 322. In the first paragraph of Section 8.16, the equation of the circle is listed as "(xx_c)+(yy_c)r=0". It should be "(xx_c)^2+(yy_c)^2r^2=0". 
14 January 2004, page 385. The equation of the plane is listed "ax+by+cz+d=". It should be "ax+by+cz+d=0". (The zero is missing.) 
14 January 2004, page 392. The last paragraph has an inequality "n_i * (P  V_{3i} < 0". It should be "n_i * (P  V_{3i}) < 0". (The right parenthesis is missing.) 
14 January 2004, page 403. Equation (10.9) has "gradient q(X)" on the lefthand side. This should be "gradient q(x)". (The function argument is upper case but needs to be lower case.) 
14 January 2004, page 409.
Item 4 for the stopping criteria of Newton's method has
(u_{i+1}  u_i S_u(u_i,v_i)) + (v_{i+1}  v_i S_u(u_i,v_i)) <= epsilon_2This should be (u_{i+1}  u_i) S_u(u_i,v_i) + (v_{i+1}  v_i) S_v(u_i,v_i) <= epsilon_2(Two misplaced parentheses. The subscript on the second partial derivative was u, but needs to be v.) 
14 January 2004, page 831. The "09 January 2003" entry mentions that the first displayed equation "(1) x13*x2+2*x3=6" should be "(1) x13*x22x3=6". The same correction must be made in the next occurrence of this equation near the bottom of the page (after "which gives us a new system"). 
14 January 2004, page 832.
The pseudocode at the top of the page has a loop "(for k = 1 to n)". This
should be "for (k = j+1 to n)". In that loop after the multiplier "m" is
computed, you could set a_{k,j} = 0 to make it clear you are eliminating
the column entries. If an implementation allows for the rowreduced matrix
to be used in further calculations, the assignment of zero should occur in
the code. But if the only goal is to solve the system, an implementation need
not assign zero since a_{k,j} is never used again in the code.
The first line of pseudocode for back substitution has "x_n = c_n/a_{n,h}". This should be "x_n = c_n/a_{n,n}". 
08 January 2004, page 305. The fourth sentence of the first paragraph in Section 8.9 has "r0 + 2*r" and "r0  2*r". These should be "r0 + 2*r1" and "r0  2*r1". 
08 January 2004, page 319. Figure 8.32 has a vector d1 for the top line. This should be d0 since the lines are parallel and are using the same direction vector. The first and second displayed equations are identical. I think the first one was intended to be "L1(t) = P1 + t*d0". The third displayed equation should be "L1(t) = P0 + d*hat(d0)^perp + t*hat(d0)". Following this displayed equation is the statement "if L0 is normalized". The term "normalized line" is not a standard one. The intent is that d0 is a normalized vector, meaning it has length 1. 
25 November 2003, page 377. The displayed equation for "e = e1*(BP)" should not have the leading minus sign. It should be "e = e1*(BP)". 
25 November 2003, page 380.
The pseudocode at the bottom of the page is missing some conditions that
lead to processing regions 0 and 5. The first part of the pseudocode should
be
if (s+t <= det) { if (s < 0) { if (t < 0) { region 4 } else { region 3 } } else if (t < 0) { region 5 } else { region 0 } }The document from which that was taken Distance Between Point and Triangle in 3D has the correct information. The algorithm is implemented and available at the source code page at this website (see the Distance web page, files WmlDistVec3Tri3.{h,cpp}). 
17 November 2003, page 838. The last set of equations with the d[i][j] should have alpha[2][2] replaced by alpha[2][0] and beta[2][2] replaced by beta[2][0]. 
All entries below this point are corrected in the second printing of the book. 
21 October 2003, page 504.
Section 11.3.2 starts out with a construction to find the points of
intersection of a line and a sphere. The line is parameterized by
P+t*D where D is a unitlength vector and t is a real number. The sphere
has center C and radius R. A point X on the sphere satisfies XC^2 = R^2.
The intersection calculation is reduced to computing the roots of a
quadratic polynomial q(t). The number of intersections is the number of
realvalued roots. This is a findintersection query.
For raysphere intersection, it is sufficient to compute the roots of q(t) and keep only those that are nonnegative. The argument is simply that the ray is P+t*D where t >= 0. The top of page 504 describes a method to inexpensively determine if there is an intersection. If the origin P of the ray is inside the sphere, then the ray must intersect the sphere. The point is inside when PC^2 <= R^2. This is a testintersection query. [If you actually need the point(s) of intersection, you still have to compute the roots of q(t). If P is outside the sphere, the ray might or might not intersect the sphere.] The next paragraph starts with "In the case of a line segment, the same approach could be used, but instead checking that the root or roots are bounded between 0 and 1, inclusive". First, the statement is in error about the bound 1. The line segment is P+t*D for 0 <= t <= tmax for some positive value tmax that is not necessarily 1 (it is 1 only if the line segment has length 1 since D = 1). Ignoring this inaccuracy, in my opinion the implication of the statement is that a testintersection query will be formulated. The remaining argument is based on a brief article at this web page. This article, as well as its restatement in the book, strays from the argument, although the web page has a clear thesis statement: If the line containing the line segment does not intersect the sphere, then the line segment does not intersect the sphere . Indeed, this is a testintersection query that attempts to quickly determine there is no intersection. The raysphere testintersection query is an attempt at showing there is an intersection. The essence of the argument given at the web page and in the book is:
Maybe the intent of the article/book was the following:

15 October 2003, page 834.
In the paragraph "The row reduction..." is "The condition for existence
of solutions is...". The equations should be corrected to
0 = a0*b1  a1*b0 = (c0+c1*y)*b1  a1*(d0+d1*y) = (b1*c0a1*d0) + (b1*c1a1*d1)*yLater is "as long as b1*c1b0*d1" which should be "as long as b1*c1a1*d1". Finally, the last equation of the last displayed equation array should be = (c0*f0e0*d0)+((c0*f1e0*d1)+(c1*f0e1*d0))*y + (c1*f1e1*d1)*y^2 
14 October 2003, page 125. The first sentence is: "This makes it completely obvious that the perp dot product reflects not only the angle between two vectors but the direction (that is, the sign) of the angle between them". I think what was intended by this is the following. Let u = (u1,u2) and v = (v1,v2). Then uperp = (u2,u1) and Dot(uperp,v) = u1*v2u2*v1. Write u and v as vectors in 3D, say u' = (u1,u2,0) and v' = (v1,v2,0). Then Cross(u',v') = (0,0,u1*v2u2*v1) = (0,0,Dot(uperp,v)). If you reverse the roles of u and v, then Dot(vperp,u) = Dot(uperp,v) and Cross(v',u') = (0,0,Dot(vperp,u)). The "sign of the angle" is a reference to measuring a counterclockwise (positive) or clockwise (negative) angle between two vectors when you cross the two using the righthand rule. 
14 October 2003, page 135.
The revised proof mentioned in the July 25, 2003 entry is better served with
a couple of additional steps:
T(w) = T(QP) by definition of w = T(Q)  T(P) since T is an affine transformation = (Q+u)  (P+u) applying the definition of T = Q  P applying identities on page 80 and 81 to cancel the u terms = (P+w)  P applying the definition of Q on page 134 = (PP)+w by identity iii on page 81 = zero_vector+w by identity i on page 80 = w zero_vector is the additive identity for vector addition 
14 October 2003, page 135. The clause "that is, the first three rows of the matrix..." should say "that is, the first n rows of the matrix...". 
14 October 2003, page 138. The sentence before the displayed equation for T(v2) ends with "... transformed coordinates of basis vector v1:". This should end with "v2", not "v1". 
14 October 2003, page 141. The last displayed equation has "T(v) = Dot(v,uhat) + ...". This should be "T(v) = Dot(v,uhat) uhat + ...". The derivation of the formula uses identity (iii) on page 91, something that is not mentioned. 
14 October 2003, page 152. The displayed equation that occurs after equation (4.32) ends with "u hat". It should end instead with "d hat perp". 
14 October 2003, page 162. Equation (4.38) starts with "P = ". It should start with "T(P) = ". 
14 October 2003, page 705. The second row of Table 13.1 has the four "Update" entries "i 0 m p". These entries should be "i 0 p m". 
20 September 2003, page 24. The second equation should be "2*x1 + 4 = 12 + 17*x2  5*x2" in order that the second equation in the first displayed set of equations on page 25 be correct. 
20 September 2003, page 51. The third displayed equation has "a(1,1)*x1 + x2*(...)". This should be "a(1,1)*x1 + a(1,2)*(...)". 
20 September 2003, page 55. The second sentence of the second paragraph has "... be a set of basis vectors for a vector space vec(v)" where vec(v) denotes the letter v with an arrow over it. This should not be vec(v), but calligraphy(V) which is the symbol in the first paragraph "Given a Eucliden space calligraphy(V)...". 
20 September 2003, pages 86, 90, 92.
The notation for length of a vector v has been v (double bars).
To be consistent with notation, the following changes should be made.
On page 86, Figure 3.27 has single bars on vector v and should be double bars. On page 90, Equations (3.8), (3.9), and (3.10) have single bars and should be double bars. On page 92, Figure 3.31 has single bars on vector v and should be double bars. 
20 September 2003, page 94.
The construction at the bottom of the page has some inconsistent notation.
A couple of the "reasons" on the righthand side need to be revised.

15 September 2003, page 98. Figure 3.35 has a vector 3*v that should be 3*v_1 (the subscript 1 is missing). 
15 September 2003, pages 128129. On page 128, second paragraph from the bottom, the second affine frame is listed as F_G and has origin O_G. These should be F_B and O_B. Figure 4.6 has the origins swapped. What is typeset as O_A should be O_B, and vice versa. 
07 September 2003, page 95. The equations at the bottom of page 94 use "psi" as the angle between vectors u and v. Figure 3.34 uses "phi". I changed the figure to use "psi". 
04 July 2003, pages 303304.
The equations for the parallel lines at the bottom of the page are missing
an "r" that should multiply the square roots.
a0*x + b0*y + c0 +/ r*sqrt(a0*a0+b0*b0) = 0 a1*x + b1*y + c1 +/ r*sqrt(a1*a1+b1*b1) = 0The equations for x and y on page 304 are correct, but the pseudocode has sign errors. It should be center[0].x = (L1.b*(L0.c+discrm0)L0.b*(L1.c+discrm1))*invDenom; center[0].y = +(L1.a*(L0.c+discrm0)L0.a*(L1.c+discrm1))*invDenom; center[1].x = (L1.b*(L0.c+discrm0)L0.b*(L1.cdiscrm1))*invDenom; center[1].y = +(L1.a*(L0.c+discrm0)L0.a*(L1.cdiscrm1))*invDenom; center[2].x = (L1.b*(L0.cdiscrm0)L0.b*(L1.c+discrm1))*invDenom; center[2].y = +(L1.a*(L0.cdiscrm0)L0.a*(L1.c+discrm1))*invDenom; center[3].x = (L1.b*(L0.cdiscrm0)L0.b*(L1.cdiscrm1))*invDenom; center[3].y = +(L1.a*(L0.cdiscrm0)L0.a*(L1.cdiscrm1))*invDenom; 
02 July 2003, page 304. The pseudocode at the bottom of the page has two occurences of "r" in the body. These should be "radius" to match the input parameter. 
25 June 2003, pages 134135.
The second block of displayed equations on page 134 has vector w in two
places. These should be vector u:
T(origin) = origin + u = u + originThe "proof" at the top of page 135 has problems. A formal proof is T(w) = T(QP) by definition of w = T(Q)  T(P) since T is an affine transformation = (P+w)  P applying the definition of T = (PP)+w by identity iii on page 81 = zero_vector+w by identity i on page 80 = w zero_vector is the additive identity for vector addition 
21 June 2003, page 21. The first displayed equation is a dot product whose second term on the righthand side is a1*b2. This termshould be a2*b2. 
20 June 2003, page 698.
The pseudocode for PointInSubpolygon has a call to GetMiddleIndex that is
missing its last parameter. The call should be
mid = GetMiddleIndex(i0, i1, C.N); 
29 May 2003, pages 535536.
On page 535, the pseudocode for intersection of triangle and plane has
a parameter "p2" in the function signature; it should not be there. On
page 536, the pseudocode has some test for points either all above the
plane or all below the plane. These tests are incorrect. In the style of
the remainder of the pseudocode, the tests are simply stated as
if ( d1 > 0.0 && d2 > 0.0 && d3 > 0.0 ) { // all points above plane return false; } if ( d1 < 0.0 && d2 < 0.0 && d3 < 0.0 ) { // all points below plane return false; } 
28 May 2003, pages 194195.
The displayed equations involving current minimum distance u have the
inequalities in the wrong direction. On page 194 the statement
x(i)  a <= u and x(i+1)  a <= ushould be x(i)  a >= u and x(i+1)  a >= uA similar change for the yequations on page 194: y(i)  b <= u and y(i+1)  b <= ushould be y(i)  b >= u and y(i+1)  b >= uPage 195 has the squared versions that should be x(i)  a^2 >= u^2 and x(i+1)  a^2 >= u^2and y(i)  b^2 >= u^2 and y(i+1)  b^2 >= u^2 
26 May 2003, page 106. The end of the first paragraph of section 3.5.1 says "Of course, since the Q_i are representable in terms of A, we could rewrite R in terms of B". The last part should by "in terms of A". 
25 May 2003, page 95. The first displayed formula for Vol(u,v,w), case Dot(Cross(u,v),w) should have w_parallel parallel to Cross(u,v). However, just as in the May 24th entry below, the use of "parallel" is incorrect terminology. The first case should state "when w_parallel is in the same direction as Cross(u,v)" and the second case should state "when w_parallel is in the opposite direction of Cross(u,v)". 
24 May 2003, pages 9091.
The proof of the distributivity of the dot product over vector addition
mentions two cases involving a vector u_parallel. I suspect what is intended
is that you look at the projections of v, w, and v+w onto the vector u. If
the angle between v+w and u is acute, then either the angle between v and u
is acute or the angle between w and u is acute (possibly both angles are
acute). If angle(v,u) > 0 and angle(w,u) > 0, then
length(v_perp+w_perp) = length(v_perp) + length(w_perp)This is case (a) in the book. If angle(v,u) > 0 and angle (w,u) < 0, then length(v_perp+w_perp) = length(v_perp)  length(w_perp)This is case (b) in the book. Of course you must also handle the case when angle(v+w,u) < 0. Two similar cases occur. I also suspect that u_parallel is meant to be the component of u in the udirection, which is just u itself. The statment u_parallel parallel to v_parallel must mean that angle(u_parallel,v_parallel) = 0. The statment u_parallel not parallel to v_parallel must mean angle (u_parallel,v_parallel) = 180 degrees. The use of the word "parallel" in this context is not good since regardless of the angle you have parallel vectors. 
19 May 2003, page 505.
The displayed equation for the quadratic coefficient c is missing a
subtraction of the squared radius. It should be
c = k*(PxCx)^2 + ell*(PyCy)^2 + m*(PzCz)^2  r^2The pseudocode at the bottom of the page is also missing the r^2, and the scaling factors need to be members of ellipsoid e. It should be c = e.k*(ell.p.x  e.center.x)^2 + e.ell*(ell.p.y  e.center.y)^2 + e.m*(ell.p.z  e.center.z)^2  e.r^2;The pseudocode uses 'l' (lowercase ell) for variable names, making it difficult to distinguish between '1' (one) and 'l' (ell). The standard axisaligned representation for an ellipsoid does not separate out the radius term. The standard form is ((xCx)/Lx)^2 + ((yCy)/Ly)^2 + ((zCz)/Lz)^2 = 1where the center is (Cx,Cy,Cz) and Lx, Ly, and Lz are the lengths measured along the axes from the center to the extreme points of the ellipsoid in those directions. I prefer more verbose names for the member names. The line has an origin and a direction. The ellipsoid has a center and length values. Another pass at the pseudocode is given below. int LineEllipsoidIntersect (Ellipsoid E, Line3 L, float t[2]) { Vector3 diff = L.origin  E.center; float invLx = 1/E.Lx; float invLy = 1/E.Ly; float invLz = 1/E.Lz; float tmp1x = L.direction.x * invLx; float tmp1y = L.direction.y * invLy; float tmp1z = L.direction.z * invLz; float tmp2x = diff.x * invLx; float tmp2y = diff.y * invLy; float tmp2z = diff.z * invLz; float a = tmp1x * tmp1x + tmp1y * tmp1y + tmp1z * tmp1z; float b = 2 * (tmp1x * tmp2x + tmp1y * tmp2y + tmp1z * tmp2z); float c = tmp2x * tmp2x + tmp2y * tmp2y + tmp2z * tmp2z  1; float d = b * b  4 * a * c; if ( d > 0 ) { float invTwoA = 1/(2*a); float rootD = sqrt(d); t[0] = (b + rootD)*invTwoA; t[1] = (b  rootD)*invTwoA; return 2; } else if ( discrm == 0 ) { t[0] = b/(2*a); return 1; } return 0; // no solutions }The return value of the function is the number of solutions. If t[i] is valid (0 <= i < numberOfSolutions), then the corresponding point of intersection is: L.origin + t[i]*L.direction. 
13 May 2003, page 487.
The top of the page makes mention that "u+v <= 1" is necessary for the
intersection, but the pseudocode does not make this check. The block of
code starting with the comparison of v to 0 and 1 and ending with the
assignment of t should be
if ( v < 0.0  v > 1.0 ) return false; if ( u + v > 1.0 ) return false; t = tmp * Dot(e2,q); 
08 Apr 2003, page 410, 411, and 412.
On page 410, the first displayed equation after (10.13) is
v = Q0  Q1 = P0 + sc*d0  P1 + tc*d1should have parentheses v = Q0  Q1 = P0 + sc*d0  (P1 + tc*d1)However, the pair of displayed equations following these are correct. The third line from the bottom defines "b = d0*d1", but should be "b = d0*d1" to be consistent with the pseudocode. The next to last line is missing a parenthesis and should be "f = (P0P1)*(P0P1)". On page 411, with the corrected "b = d0*d1", the formula at the top of the page for tc is wrong. It should be tc = (a*eb*d)/(a*cb*b)The formula is wrong even if you used the incorrect "b = d0*d1". This change propagates to the next displayed equation: L0(sc)L1(tc) = (P0P1)+ ((b*ec*d)d0(a*eb*d)d1)/(a*cb*b)The block of pseudocode for "det < epsilon" has "return d*s+f" which is wrong. It should be "return e*t+f". On page 412, the top line has the correct formula for t despite the errors in the text. However, the return value is wrong and should be return s*(a*sb*t+2*d)+t*(b*s+c*t2*e)+f; 
25 Mar 2003, page 366.
The pseudocode has a line
t = Dot(l.direction, VectorSubtract(q, l.direction));This should be t = Dot(l.direction, VectorSubtract(q, l.origin)); 
15 Feb 2003, page 487.
The pseudocode has a line near the bottom
v = tmp * Dot(d, q);This should be v = tmp * Dot(line.direction, q); 
02 Feb 2003, page 628629. The pseudocode on page 628, the case where "Ray parallel to planes" has a test "P.x < box.xMax". That should be "P.x > box.xMax". The pseudocode on page 629 for swapping t0 and t1 is incorrect. It should be "tmp = t1; t1 = t0; t0 = tmp;" 
02 Feb 2003, page 502.
The quadratic equation in the second line should say "at^2+bt+c=0" to be
consistent with the footnote and the pseudocode. Equation (11.11) should be
t = d*(PC) + sqrt((d*(PC))^2((PC)*(PC)r^2))and is based on d being a unitlength vector (a = 1 in the quadratic equation). Note, though, that the pseudocode does not assume d is unit length. The website does have an implementation of this algorithm, the Intersection 3D page. (The expression in the line after equation (11.11) has an extraneous right parenthesis.) 
09 Jan 2003, page 831. The fourth line from the top of the page is "(1) x1  3*x2 + 2*x3 = 6". It should be "(1) x1  3*x2  2*x3 = 6". The last line of the page lists the solution as (1,3,2). It should be (1,3,2). 
17 Dec 2002, page 67. The distributive rule (iii) should be (a+b)u = au + bu (the text has bv). 
12 Dec 2002, pages 641 and 642. The last line of page 641 has a lowercase subscript 'a' on the unit vector u. That should be an uppercase 'A'. The last line of page 642 also has a lowercase subscript of 'a', but it should be an uppercase 'B'. 
11 Oct 2002, page 77. The last sentence of the first paragraph states "Many values in the range..." but should be "Many values in the domain..." 
Book Corrections Organized by Page Number 
Section 2.2, page 15. The "Matrix" row of Table 2.1 has pseudocode "Matrix3x3 m, m;" that should be "Matrix3x3 m, n;". 
Section 2.3.4, page 21. The first displayed equation is a dot product whose second term on the righthand side is a1*b2. This termshould be a2*b2. 
Section 2.4.1, page 24. The second equation should be "2*x1 + 4 = 12 + 17*x2  5*x2" in order that the second equation in the first displayed set of equations on page 25 be correct. 
Section 2.6.3, page 43. Bullet 'v' in Section 2.6.3 is a repeat of bullet 'ii' and should be removed. 
Section 2.6.5, page 44.
The example combination of vectors in a displayed equation is
2v_{1} + 3v_{3} +1v_{3}but should be 2v_{1} + 3v_{2} +1v_{3}That is, the subscript on the middle term needed to be fixed. 
Section 2.7.1, page 46. The sentence: "We can map real numbers to real numbers: ..." should say "We can map nonnegative real numbers fo nonnegative real numbers: ..." 
Section 2.7.4, page 51. The third displayed equation has "a(1,1)*x1 + x2*(...)". This should be "a(1,1)*x1 + a(1,2)*(...)". 
Section 2.8, page 53. The last sentence should say "... into Equations 2.4 and 2.5, ...". 
Section 2.9.2, page 55. The second sentence of the second paragraph has "... be a set of basis vectors for a vector space vec(v)" where vec(v) denotes the letter v with an arrow over it. This should not be vec(v), but calligraphy(V) which is the symbol in the first paragraph "Given a Eucliden space calligraphy(V)...". 
Section 2.10, page 59.
The third displayed equation should be
D^2 = (b + mx1  y1)^2 + (b + mx2  y2)^2 + (b + mx3  y3)^2The book has yn in the last term of the expression. 
Section 3.1.5, page 67. The distributive rule (iii) should be (a+b)u = au + bu (the text has bv). 
Section 3.2.6, page 77. The last sentence of the first paragraph states "Many values in the range..." but should be "Many values in the domain..." 
Section 3.2.6, page 78.
Figure 3.19 is confusing and has one typographical error. It is an attempt to
describe the linear transformation T(x) = 2*x. In terms of basis vectors u
and v: T(u) = 2*u and T(v) = 2*v. The vector s in the figure should be more
clearly written:
s = T(2*v) = 2*T(v) // because T is linear = 2*(2*v) // by definition of TCall the top vector of the figure r, in which case r = T(3*u) = 3*T(u) // because T is linear = 3*(2*u) // by definition of TThus, the current figure should have "3T(u) = 2(3u)". 
Section 3.2.6, page 79.
Figure 3.20 is confusing and has two typographical errors. It is an attempt to
describe the linear transformation that is defined on the basis vectors u and v:
T(u) = 2*u and T(v) = 1.5*v. Given a vector x = a*u+b*v, the transformation is
T(x) = T(a*u+b*v) = a*T(u)+b*T(v) = a*(2*u)+b*(1.5*v). The vector s in the
figure is incorrect and should be
s = T(2*v) = 2*T(v) // because T is linear = 2*(1.5*v) // by definition of TCall the top vector of the figure r, in which case r = T(3*u) = 3*T(u) // because T is linear = 3*(2*u) // by definition of TThus, the current figure should have "3T(u) = 2(3u)". 
Section 3.3, page 83. The sentence starting with "We've just done two forbidden operations on pointsscaling by a vector and directly adding them." This should say "We've just done two forbidden operations on pointsmultiplying points by scalars and directly adding them." 
Section 3.3, page 84. The first displayed equation has a term "alpha3*(P3P2)". This term should be "alpha3*(P3P1)". 
Section 3.3.1, pages 86, 90, 92.
The notation for length of a vector v has been v (double bars).
To be consistent with notation, the following changes should be made.
On page 86, Figure 3.27 has single bars on vector v and should be double bars. On page 90, Equations (3.8), (3.9), and (3.10) have single bars and should be double bars. On page 92, Figure 3.31 has single bars on vector v and should be double bars. 
Section 3.3.1, pages 9091.
The proof of the distributivity of the dot product over vector addition
mentions two cases involving a vector u_parallel. I suspect what is intended
is that you look at the projections of v, w, and v+w onto the vector u. If
the angle between v+w and u is acute, then either the angle between v and u
is acute or the angle between w and u is acute (possibly both angles are
acute). If angle(v,u) > 0 and angle(w,u) > 0, then
length(v_perp+w_perp) = length(v_perp) + length(w_perp)This is case (a) in the book. If angle(v,u) > 0 and angle (w,u) < 0, then length(v_perp+w_perp) = length(v_perp)  length(w_perp)This is case (b) in the book. Of course you must also handle the case when angle(v+w,u) < 0. Two similar cases occur. I also suspect that u_parallel is meant to be the component of u in the udirection, which is just u itself. The statment u_parallel parallel to v_parallel must mean that angle(u_parallel,v_parallel) = 0. The statment u_parallel not parallel to v_parallel must mean angle (u_parallel,v_parallel) = 180 degrees. The use of the word "parallel" in this context is not good since regardless of the angle you have parallel vectors. 
Section 3.3.1, page 92.
Item "iv. Normal:" should be
v_perp = v  Dot(u,v)*u/Dot(u,u)The book has u as the first vector in the subtraction. 
Section 3.3.2, page 94.
The construction at the bottom of the page has some inconsistent notation.
A couple of the "reasons" on the righthand side need to be revised.

Section 3.3.2, page 95. The equations at the bottom of page 94 use "psi" as the angle between vectors u and v. Figure 3.34 uses "phi". I changed the figure to use "psi". 
Section 3.3.2, page 95. The first displayed formula for Vol(u,v,w), case Dot(Cross(u,v),w) should have w_parallel parallel to Cross(u,v). However, just as in the May 24th entry below, the use of "parallel" is incorrect terminology. The first case should state "when w_parallel is in the same direction as Cross(u,v)" and the second case should state "when w_parallel is in the opposite direction of Cross(u,v)". 
Section 3.3.3, page 98. Figure 3.35 has a vector 3*v that should be 3*v_1 (the subscript 1 is missing). 
Section 3.5.1, page 106. The end of the first paragraph of section 3.5.1 says "Of course, since the Q_i are representable in terms of A, we could rewrite R in terms of B". The last part should by "in terms of A". 
Section 4.4.2, page 120. In the middle of the page, I do not understand why the statement "with, as before, Vector(u) retaining its usual matrix representation" is made. This seems to imply that Vector(u) will be represented as a row matrix as shown at the top of the page. When Cross(Vector(v),Vector(u)) is represented using a skewsymmetric matrix for Vector(v), Vector(u) must be written as a column matrix. Moreover, the skewsymmetric matrix Tilde(v) is the negative of what it should be. Listing the rows from top to bottom, it should be Tilde(v) = ((0,v3,+v2),(+v3,0,v1),(v2,+v1,0)). 
Section 4.4.4, page 125. The first sentence is: "This makes it completely obvious that the perp dot product reflects not only the angle between two vectors but the direction (that is, the sign) of the angle between them". I think what was intended by this is the following. Let u = (u1,u2) and v = (v1,v2). Then uperp = (u2,u1) and Dot(uperp,v) = u1*v2u2*v1. Write u and v as vectors in 3D, say u' = (u1,u2,0) and v' = (v1,v2,0). Then Cross(u',v') = (0,0,u1*v2u2*v1) = (0,0,Dot(uperp,v)). If you reverse the roles of u and v, then Dot(vperp,u) = Dot(uperp,v) and Cross(v',u') = (0,0,Dot(vperp,u)). The "sign of the angle" is a reference to measuring a counterclockwise (positive) or clockwise (negative) angle between two vectors when you cross the two using the righthand rule. 
Section 4.6, pages 128129. On page 128, second paragraph from the bottom, the second affine frame is listed as F_G and has origin O_G. These should be F_B and O_B. Figure 4.6 has the origins swapped. What is typeset as O_A should be O_B, and vice versa. 
Section 4.7.2, pages 134135.
The second block of displayed equations on page 134 has vector w in two
places. These should be vector u:
T(origin) = origin + u = u + originThe "proof" at the top of page 135 has problems. A formal proof is T(w) = T(QP) by definition of w = T(Q)  T(P) since T is an affine transformation = (P+w)  P applying the definition of T = (PP)+w by identity iii on page 81 = zero_vector+w by identity i on page 80 = w zero_vector is the additive identity for vector addition 
Section 4.7.2, page 135.
The revised proof mentioned in the previous entry is better served with
a couple of additional steps:
T(w) = T(QP) by definition of w = T(Q)  T(P) since T is an affine transformation = (Q+u)  (P+u) applying the definition of T = Q  P applying identities on page 80 and 81 to cancel the u terms = (P+w)  P applying the definition of Q on page 134 = (PP)+w by identity iii on page 81 = zero_vector+w by identity i on page 80 = w zero_vector is the additive identity for vector addition 
Section 4.7.2, page 135. The clause "that is, the first three rows of the matrix..." should say "that is, the first n rows of the matrix...". 
Section 4.7.3, page 138. The sentence before the displayed equation for T(v2) ends with "... transformed coordinates of basis vector v1:". This should end with "v2", not "v1". 
Section 4.7.3, page 141. The last displayed equation has "T(v) = Dot(v,uhat) + ...". This should be "T(v) = Dot(v,uhat) uhat + ...". The derivation of the formula uses identity (iii) on page 91, something that is not mentioned. 
Section 4.7.5, page 152. The displayed equation that occurs after equation (4.32) ends with "u hat". It should end instead with "d hat perp". 
Section 4.7.6, pages 157158. For simplicity of notation here, define N = hat(n), N^T is the transpose of N, and V = hat(v). The last row, first column entry of the general shearing matrix should be (tan(theta))*Dot(N,Q)*V. The first line on page 158 makes an appeal to applying the point transformation T to a vector N, which makes no sense mathematically. Even if it did, it is not clear why the author claims T causes N to change at all. Let P = P1 + h*N, where h = Dot(N,PQ). The shearing is P" = P + t*V, where from trigonometry, t = h*tan(theta), so P" = P + (tan(theta))*Dot(N,PQ)*V and P"Q = [I + (tan(theta))*V*N^T](PQ). Add Q to both sides and factor to P" = [I + (tan(theta))*V*N^T]*P  (tan(theta))*V*N^T*Q = [I + (tan(theta))*V*N^T]*P  (tan(theta))*Dot(N,Q)*V. 
Section 4.8.2, page 162. Equation (4.38) starts with "P = ". It should start with "T(P) = ". 
Section 5.6.2, page 187. The displayed equation for the recursion formulas has "B_{i,0}(t)" but should be "B_{i,1}(t)". The displayed equation after that one is the recursion, but should be for "2 <= j <= n" (the text has 1, not 2). 
Section 6.2, pages 194195.
The displayed equations involving current minimum distance u have the
inequalities in the wrong direction. On page 194 the statement
x(i)  a <= u and x(i+1)  a <= ushould be x(i)  a >= u and x(i+1)  a >= uA similar change for the yequations on page 194: y(i)  b <= u and y(i+1)  b <= ushould be y(i)  b >= u and y(i+1)  b >= uPage 195 has the squared versions that should be x(i)  a^2 >= u^2 and x(i+1)  a^2 >= u^2and y(i)  b^2 >= u^2 and y(i+1)  b^2 >= u^2 
Section 6.3.1, pages 208210. The algorithm uses barycentric coordinates to determine closest edge, but this does not work for a triangle with an obtuse angle. Thus, the SquaredDistance pseudocode is incorrect. The online source code does not use this algorithm; it uses the previously discussed (correct) algorithm. 
Section 6.6.4, page 225. The denominator of Equation (6.17) should be dot(perp(d0),d1) rather than dot(perp(d1),d0); this has the effect of negating the righthand side of the equation. 
Section 6.6.6, page 228.
The symbols in Equation (6.23) on page 228 are not typeset correctly. In that equation,
let us say that the first row of the equation is region R0 and the last row of the
equation is R8 (the parameter plane is partitioned into 9 regions). The corrections are
R2: 0 < tilde(t0) < T0 and bar(t1) >= T1 R4: 0 < tilde(t1) < T1 and bar(t0) >= T0 R6: hat(t0) >= T0 and tilde(t1) <= 0 R7: tilde(t0) <= 0 and hat(t1) >= T1 R8: tilde(t0) >= T0 and tilde(t1) >= T1 
Section 6.6.6, page 229. Equation (6.24) has 5 cases (top to bottom). The final inequality in case 2 should be Dot(d0,delta+T0*d0) ≤ 0. The final inequality in case 4 should be Dot(d0,delta+T0*d0T1*d1) ≤ 0. 
Section 8.4, page 290.
The last line of code is missing a left parenthesis. It should be
v0.y = ((r2 * u.y) + rad) / denom; 
Section 8.4, page 291.
The first line of pseudocode with solution[0] has v.x when it should be v0.x:
solution[0].setDir(v0.y, v0.x); 
Section 8.8, pages 303304.
The equations for the parallel lines at the bottom of the page are missing
an "r" that should multiply the square roots.
a0*x + b0*y + c0 +/ r*sqrt(a0*a0+b0*b0) = 0 a1*x + b1*y + c1 +/ r*sqrt(a1*a1+b1*b1) = 0The equations for x and y on page 304 are correct, but the pseudocode has sign errors. It should be center[0].x = (L1.b*(L0.c+discrm0)L0.b*(L1.c+discrm1))*invDenom; center[0].y = +(L1.a*(L0.c+discrm0)L0.a*(L1.c+discrm1))*invDenom; center[1].x = (L1.b*(L0.c+discrm0)L0.b*(L1.cdiscrm1))*invDenom; center[1].y = +(L1.a*(L0.c+discrm0)L0.a*(L1.cdiscrm1))*invDenom; center[2].x = (L1.b*(L0.cdiscrm0)L0.b*(L1.c+discrm1))*invDenom; center[2].y = +(L1.a*(L0.cdiscrm0)L0.a*(L1.c+discrm1))*invDenom; center[3].x = (L1.b*(L0.cdiscrm0)L0.b*(L1.cdiscrm1))*invDenom; center[3].y = +(L1.a*(L0.cdiscrm0)L0.a*(L1.cdiscrm1))*invDenom; 
Section 8.8, page 304. The pseudocode at the bottom of the page has two occurences of "r" in the body. These should be "radius" to match the input parameter. 
Section 8.9, page 305. The fourth sentence of the first paragraph in Section 8.9 has "r0 + 2*r" and "r0  2*r". These should be "r0 + 2*r1" and "r0  2*r1". 
Section 8.14, page 319. Figure 8.32 has a vector d1 for the top line. This should be d0 since the lines are parallel and are using the same direction vector. The first and second displayed equations are identical. I think the first one was intended to be "L1(t) = P1 + t*d0". The third displayed equation should be "L1(t) = P0 + d*hat(d0)^perp + t*hat(d0)". Following this displayed equation is the statement "if L0 is normalized". The term "normalized line" is not a standard one. The intent is that d0 is a normalized vector, meaning it has length 1. 
Section 8.16, page 322. In the first paragraph of Section 8.16, the equation of the circle is listed as "(xx_c)+(yy_c)r=0". It should be "(xx_c)^2+(yy_c)^2r^2=0". 
Section 9.2.1, page 328. The left image of Figure 9.2 has "r = d" and the right image has "r = d". To be consistent with the paragraph on this page "More significantly...", these should be reversed. The origin in the left image has a negative signed distance from the plane (it is "behind" the plane). The origin in the right image has a positive signed distance from the plane (it is "in front of" the plane). 
Section 9.5, page 356.
The first equation for the torus is incorrect. It should be the following
equation. I have factored it in a manner that leads to more efficient
numerical computation:
(x^2 + y^2 + z^2)^2  2*(r0^2 + r1^2)*(x^2 + y^2 + z^2) + 4*r0^2*z^2 + (r0^2  r1^2)^2 = 0My torus class (in Wm3Torus files) has the correct formula listed. 
Section 10.2, page 366.
The pseudocode has a line
t = Dot(l.direction, VectorSubtract(q, l.direction));This should be t = Dot(l.direction, VectorSubtract(q, l.origin)); 
Section 10.3.2, page 377. The displayed equation for "e = e1*(BP)" should not have the leading minus sign. It should be "e = e1*(BP)". 
Section 10.3.2, page 377. The equation at the bottom of the page for bar(t) has the wrong numerator. The equation should be "bar(t) = (b*da*e)/(a*cb*b)". The pseudocode on page 380 does use the correct equation. 
Section 10.2.1, page 367. In the case t0 > 0, the distance from a point Q to a ray P+t*D should be Q(P+t0*D). The book has Q(Q+t0*D), a typographical error. 
Section 10.3, page 380.
The pseudocode at the bottom of the page is missing some conditions that
lead to processing regions 0 and 5. The first part of the pseudocode should
be
if (s+t <= det) { if (s < 0) { if (t < 0) { region 4 } else { region 3 } } else if (t < 0) { region 5 } else { region 0 } }The document from which that was taken Distance Between Point and Triangle in 3D has the correct information. The algorithm is implemented and available at the source code page at this website (see the Distance web page, files WmlDistVec3Tri3.{h,cpp}). 
Section 10.3.3, page 384.
The first displayed equation has a Q that should be a P:
Q' = P + s*e0 + t*e1 
Section 10.3.4, page 385. The equation of the plane is listed "ax+by+cz+d=". It should be "ax+by+cz+d=0". (The zero is missing.) 
Section 10.3.4, page 387. The top line of pseudocode has "Dot(p,n)d", but should be "Dot(p,n)+d". 
Section 10.4.1, page 392. The last paragraph has an inequality "n_i * (P  V_{3i} < 0". It should be "n_i * (P  V_{3i}) < 0". (The right parenthesis is missing.) 
Section 10.5.2, page 403. Equation (10.9) has "gradient q(X)" on the lefthand side. This should be "gradient q(x)". (The function argument is upper case but needs to be lower case.) 
Section 10.7, page 409.
Item 4 for the stopping criteria of Newton's method has
(u_{i+1}  u_i S_u(u_i,v_i)) + (v_{i+1}  v_i S_u(u_i,v_i)) <= epsilon_2This should be (u_{i+1}  u_i) S_u(u_i,v_i) + (v_{i+1}  v_i) S_v(u_i,v_i) <= epsilon_2(Two misplaced parentheses. The subscript on the second partial derivative was u, but needs to be v.) 
Section 10.8.1, page 410, 411, and 412.
On page 410, the first displayed equation after (10.13) is
v = Q0  Q1 = P0 + sc*d0  P1 + tc*d1should have parentheses v = Q0  Q1 = P0 + sc*d0  (P1 + tc*d1)However, the pair of displayed equations following these are correct. The third line from the bottom defines "b = d0*d1", but should be "b = d0*d1" to be consistent with the pseudocode. The next to last line is missing a parenthesis and should be "f = (P0P1)*(P0P1)". On page 411, with the corrected "b = d0*d1", the formula at the top of the page for tc is wrong. It should be tc = (a*eb*d)/(a*cb*b)The formula is wrong even if you used the incorrect "b = d0*d1". This change propagates to the next displayed equation: L0(sc)L1(tc) = (P0P1)+ ((b*ec*d)d0(a*eb*d)d1)/(a*cb*b)The block of pseudocode for "det < epsilon" has "return d*s+f" which is wrong. It should be "return e*t+f". On page 412, the top line has the correct formula for t despite the errors in the text. However, the return value is wrong and should be return s*(a*sb*t+2*d)+t*(b*s+c*t2*e)+f; 
Section 10.8, pages 417426.
These pages cover distance calculations between linear components. The
pseudocode has major problems. Apparently, one of these blocks of code was
written and then cutpastemodified for the other blocks. They all exhibit
the same ills. The difference between closest points needs to be of the
form
v = (base0 + s*direction0)  (base1 + t*direction1)That is, the book parentheses are incorrect. Worse is the handling of the distinction between parallel and nonparallel components. If "det" is ever zero, the pseudocode will do a division by it. Also, the epsilon test is problematic when the segment direction vectors have large length. If D0 and D1 are the line segment lengths, each computed as a difference of the line segment end points, and if A is the angle between the directions, it is the case that det = Cross(D0,D1)^2 = D0^2 * D1^2 * (sin(A))^2If the angle A is very small, but the lengths of the directions are very large, "det" can be quite large, causing you to classify the segments are nonparallel when they are nearly parallel. My own implementation uses a line segment represention of C+s*D, where C is the midpoint of the segment, D is a unitlength direction, and s <= e where e is half the length of the segment. This leads to a test det = (sin(A))^2so that the epsilon test is a relative error measurement, not an absolute error measurement. Section 10.8.3 of the book describes how I think of distance measurements for linear components, and I much prefer this approach based on actually writing code and seeing how it performs on real data. 
Section 10.8.2, page 419.
The pseudocode for c should be
c = Dot(ray.direction, ray.direction); 
Section 10.9.1, page 437.
The displayed matrix equation is incorrect. Setting the gradient to zero
leads to the system of equations
a00*u + a01*v + a02*t + b0 = 0 a01*u + a11*v + a12*t + b1 = 0 a02*u + a12*v + a22*t + b2 = 0The pseudocode on that page essentially solves this system using Cramer's rule (computes cofactors and determinant). 
Section 10.8.3, page 428. The paragraph about the case of region 2 has "the gradient at the corner (1,1) cannot point inside the unit square". This should say "outside" rather than "inside". 
Section 10.9.4, page 454.
The pseudocode after the "Closest point on the box" comment should be
qPrime.x = box.extents.x; qPrime.y = line.origin.y; qPrime.z = line.origin.z; 
Section 11.1, page 487.
The pseudocode has a line near the bottom
v = tmp * Dot(d, q);This should be v = tmp * Dot(line.direction, q); 
Section 11.1, page 487.
The top of the page makes mention that "u+v <= 1" is necessary for the
intersection, but the pseudocode does not make this check. The block of
code starting with the comparison of v to 0 and 1 and ending with the
assignment of t should be
if ( v < 0.0  v > 1.0 ) return false; if ( u + v > 1.0 ) return false; t = tmp * Dot(e2,q); 
Section 11.1.3, page 490.
The discussion in this section is about linear components (lines, rays, or
line segments) and polygons. The text sometimes talks about lines, sometimes
about rays. The pseudocode has a test "if (t < 0 ) { return false; }" which
indicates it is about raypolygon intersection, but is incorrectly named
LinePolygonIntersection. Also, there is a line
p = line.origin + t * r.d;which should be p = line.origin + t * line.direction;another indication that the author had a ray in mind. Near the end of the pseudocode is intersection = line.origin + t * line.direction;However, the intersection p was calculated earlier, so just use intersection = p; 
Section 11.1.4, page 493.
The pseudocode has the declaration
Point3D intersectionThis is not needed. 
Section 11.3, page 502.
The quadratic equation in the second line should say "at^2+bt+c=0" to be
consistent with the footnote and the pseudocode. Equation (11.11) should be
t = d*(PC) + sqrt((d*(PC))^2((PC)*(PC)r^2))and is based on d being a unitlength vector (a = 1 in the quadratic equation). Note, though, that the pseudocode does not assume d is unit length. The website does have an implementation of this algorithm, the Intersection page. (The expression in the line after equation (11.11) has an extraneous right parenthesis.) 
Section 11.3.2, page 504.
Section 11.3.2 starts out with a construction to find the points of
intersection of a line and a sphere. The line is parameterized by
P+t*D where D is a unitlength vector and t is a real number. The sphere
has center C and radius R. A point X on the sphere satisfies XC^2 = R^2.
The intersection calculation is reduced to computing the roots of a
quadratic polynomial q(t). The number of intersections is the number of
realvalued roots. This is a findintersection query.
For raysphere intersection, it is sufficient to compute the roots of q(t) and keep only those that are nonnegative. The argument is simply that the ray is P+t*D where t >= 0. The top of page 504 describes a method to inexpensively determine if there is an intersection. If the origin P of the ray is inside the sphere, then the ray must intersect the sphere. The point is inside when PC^2 <= R^2. This is a testintersection query. [If you actually need the point(s) of intersection, you still have to compute the roots of q(t). If P is outside the sphere, the ray might or might not intersect the sphere.] The next paragraph starts with "In the case of a line segment, the same approach could be used, but instead checking that the root or roots are bounded between 0 and 1, inclusive". First, the statement is in error about the bound 1. The line segment is P+t*D for 0 <= t <= tmax for some positive value tmax that is not necessarily 1 (it is 1 only if the line segment has length 1 since D = 1). Ignoring this inaccuracy, in my opinion the implication of the statement is that a testintersection query will be formulated. The remaining argument is based on a brief article at this web page . This article, as well as its restatement in the book, strays from the argument, although the web page has a clear thesis statement: If the line containing the line segment does not intersect the sphere, then the line segment does not intersect the sphere . Indeed, this is a testintersection query that attempts to quickly determine there is no intersection. The raysphere testintersection query is an attempt at showing there is an intersection. The essence of the argument given at the web page and in the book is:
Maybe the intent of the article/book was the following:

Section 11.3.3, page 505.
The displayed equation for the quadratic coefficient c is missing a
subtraction of the squared radius. It should be
c = k*(PxCx)^2 + ell*(PyCy)^2 + m*(PzCz)^2  r^2The pseudocode at the bottom of the page is also missing the r^2, and the scaling factors need to be members of ellipsoid e. It should be c = e.k*(ell.p.x  e.center.x)^2 + e.ell*(ell.p.y  e.center.y)^2 + e.m*(ell.p.z  e.center.z)^2  e.r^2;The pseudocode uses 'l' (lowercase ell) for variable names, making it difficult to distinguish between '1' (one) and 'l' (ell). The standard axisaligned representation for an ellipsoid does not separate out the radius term. The standard form is ((xCx)/Lx)^2 + ((yCy)/Ly)^2 + ((zCz)/Lz)^2 = 1where the center is (Cx,Cy,Cz) and Lx, Ly, and Lz are the lengths measured along the axes from the center to the extreme points of the ellipsoid in those directions. I prefer more verbose names for the member names. The line has an origin and a direction. The ellipsoid has a center and length values. Another pass at the pseudocode is given below. int LineEllipsoidIntersect (Ellipsoid E, Line3 L, float t[2]) { Vector3 diff = L.origin  E.center; float invLx = 1/E.Lx; float invLy = 1/E.Ly; float invLz = 1/E.Lz; float tmp1x = L.direction.x * invLx; float tmp1y = L.direction.y * invLy; float tmp1z = L.direction.z * invLz; float tmp2x = diff.x * invLx; float tmp2y = diff.y * invLy; float tmp2z = diff.z * invLz; float a = tmp1x * tmp1x + tmp1y * tmp1y + tmp1z * tmp1z; float b = 2 * (tmp1x * tmp2x + tmp1y * tmp2y + tmp1z * tmp2z); float c = tmp2x * tmp2x + tmp2y * tmp2y + tmp2z * tmp2z  1; float d = b * b  4 * a * c; if ( d > 0 ) { float invTwoA = 1/(2*a); float rootD = sqrt(d); t[0] = (b + rootD)*invTwoA; t[1] = (b  rootD)*invTwoA; return 2; } else if ( discrm == 0 ) { t[0] = b/(2*a); return 1; } return 0; // no solutions }The return value of the function is the number of solutions. If t[i] is valid (0 <= i < numberOfSolutions), then the corresponding point of intersection is: L.origin + t[i]*L.direction. 
Section 11.3.4, page 510.
The solution for the quadratic formula is missing a minus sign on the "b"
term. It should be
t = (b +/ sqrt(b^24*a*c))/(2*a)The pseudocode has a = (tline.direction.x + tline.direction.y)^2;This should be a = tline.direction.x^2 + tline.direction.y^2; 
Section 11.3.4, page 511.
Lines 4 and 5 are missing right parentheses. The pseudocode should be
t[0] = (b + sqrt(discrm))/(2*a); t[1] = (b  sqrt(discrm))/(2*a); 
Section 11.3.5, page 518. The page has 7 references to array elements "[3]". These should be "[2]". 
Section 11.5.1, page 531.
The pseudocode for b at the bottom of the page should be
b = (s1 * n1n2dot  s2 * n1normsqr) / (n1n2dot^2  n1normsqr * n2normsqr);That is, the first occurrence of n2normsqr in the book should have been n1normsqr. 
Section 11.5.2, page 534.
The equations for Cramer's rule have errors. The inverse determinant should
be
invDet = 1/(a0*BC  b0*AC + c0*AB)The Z value should be Z = (b0*AD + a0*DB  d0*AB)*invDetI do not have a copy of Bowyer and Woodwark, the reference from which the equations are supposedly taken, so I do not know if these errors occur in that reference or were introduced in transcribing. 
Section 11.5.3, pages 535536.
On page 535, the pseudocode for intersection of triangle and plane has
a parameter "p2" in the function signature; it should not be there. On
page 536, the pseudocode has some test for points either all above the
plane or all below the plane. These tests are incorrect. In the style of
the remainder of the pseudocode, the tests are simply stated as
if ( d1 > 0.0 && d2 > 0.0 && d3 > 0.0 ) { // all points above plane return false; } if ( d1 < 0.0 && d2 < 0.0 && d3 < 0.0 ) { // all points below plane return false; } 
Section 11.7.1, page 548. In matrix M, entry 6 of row 2 says "2(u_{x} v_{y} + v_{x} v_{y})". This should be "2(u_{x} v_{y} + v_{x} u_{y})". 
Section 11.7.2, page 548. Equation (11.20) has "r" when it should be "r^2" (r squared). 
Section 11.7.2, page 550. To be consistent with the other pseudocode, "if (b < r)" should be written as "if (b < sphere.radius)". 
Section 11.7.2, page 550.
The last line of pseudocode on this page should be
float radius = sqrt(sphere.radius^2  b^2);The pseudocode also rejects the "just touching" case as "not an intersection"; that is, "b < sphere.radius" should be "b ≤ sphere.radius". In fact, since the pseudocode allows a small amount of error ("radius < epsilon"), it might even be good to support an error term in the first test: "if (b ≤ sphere.radius + epsilon)". 
Section 11.7.2, page 550.
The pseudocode computes b as an absolute value, but the absolute value should
only be used in the comparison to the sphere radius on the line after
computing b. The pseudocode on this page has had other corrections, so the
revision is listed next. I do not like the book's formulationthe revision
is more like how I write the actual code. This is the mathematically correct
formulation, without concern for "epsilon" arithmetic to deal with numerical
roundoff errors. The plane normal is assumed to be a unitlength vector.
The Circle3 data structure is for a circle in three dimensions.
bool PlaneSphereIntersection (Plane plane, Sphere sphere, Circle3& circle) { float b = dotProd(plane.normal, sphere.center  plane.pointOnPlane); float absB = fabs(b); circle.center = sphere.center  b * plane.normal; circle.normal = plane.normal; if (absB <= sphere.radius) { // The sphere intersects the plane in a circle. The circle is // degenerate when fabs(b) is equal to sphere.radius, in which // case the circle radius is zero. circle.radius = sqrt(sphere.radius * sphere.radius  b * b); return true; } return false; }The circle parameter is an output of the function. It is valid only when the function's return value is true. The plane's normal vector is copied to the circle's normal so that you can actually reconstruct the circle. If the circle center is K, the normal is N, and the radius is R, then choose vectors U and V so that {U,V,N} is an orthonormal set (vectors are unit length and mutually perpendicular). A circle parameterization is X(t) = K + R*(cos(t)*U + sin(t)*V), 0 <= t < 2*piA numerically robust way to compute U and V is if (fabs(N.x) >= fabs(N.y)) { // N.x or N.z is the largest magnitude component, swap them. invLength = 1/sqrt(N.x * N.x + N.z * N.z); U = Vector(N.z * invLength, 0, N.x * invLength); V = Vector(N.y * U.z, N.z * U.x  N.x * U.z, N.y * U.x); // = Cross(N,U) } else { // N.y or N.z is the largest magnitude component, swap them. invLength = 1/sqrt(N.y * N.y + N.z * N.z); U = Vector(0, N.z * invLength, N.y * invLength); V = Vector(N.y * U.z  N.z * U.y, N.x * U.z, N.x * U.y); // = Cross(N,U) }If you want to allow for numerical roundoff error, say, by requiring the intersection to be a single point when b is within epsilon of the radius, use bool PlaneSphereIntersection (Plane plane, Sphere sphere, Circle3& circle) { float b = dotProd(plane.normal, sphere.center  plane.pointOnPlane); float absB = fabs(b); circle.center = sphere.center  b * plane.normal; circle.normal = plane.normal; if (absB <= sphere.radius  epsilon) { // The sphere intersects the plane in a circle. The circle is // degenerate when fabs(b) is equal to sphere.radius, in which // case the circle radius is zero. circle.radius = sqrt(sphere.radius * sphere.radius  b * b); return true; } else if (absB <= sphere.radius + epsilon) { // The sphere is _nearly_ tangent to the plane. Choose the // intersection to be the point of tangency. circle.radius = 0; return true; } return false; }The magnitude of the nonnegative number epsilon is up to the programmer. Observe that if you choose epsilon to be zero, you effectively have the pseudocode for the mathematically correct case. 
Section 11.7.3, page 562.
The pseudocode has a few issues. The first line of code that computes d
should be a signed distance calculation. The first comparison of
cos(theta) should have the abs about the cosine term, not
about the epsilon term. If an epsilon comparison is made for when
cos(theta) is nearly zero, to be consistent an epsilon comparison
should also be made for when cos(theta) is nearly one (the circle
case). Here is my own rewrite of the pseudocode without epsilon tests.
The following data structures are assumed.
struct Plane { Point origin; Vector normal; }; struct Cylinder { Point origin; Vector direction; float radius; }; struct Line { Point origin; Vector direction; }; struct Circle { Point center; Vector normal; float radius; }; struct Ellipse { Point center; Vector normal, major, minor; float majorLength, minorLength; }; enum Type { EMPTY_SET, ONE_LINE, TWO_LINES, CIRCLE, ELLIPSE }; struct Intersection { Type type; Line line; // valid when type = ONE_LINE Line line1, line2; // valid when type = TWO_LINES Circle circle; // valid when type = CIRCLE Ellipse ellipse; // valid when type = ELLIPSE };The plane normal, cylinder direction, line direction, circle normal, ellipse normal, ellipse major direction, and ellipse minor direction are all required to be unitlength vectors. The output is intr. The returned bool value is true when there is an intersection, in which case intr is valid (of the type specified in its Type member). The returned value is false when there is no intersection, in which case the Type member has value EMPTY_SET. bool PlaneCylinderIntersection (Plane plane, Cylinder cylinder, Intersection& intr) { float d = Dot(plane.normal, cylinder.origin  plane.origin); Point projectedCylinderOrigin = cylinder.origin  d * plane.normal; float cosTheta = Dot(cylinder.direction, plane.normal); if (abs(cosTheta) > 0) // cylinder axis intersects plane in a unique point { if (abs(cosTheta) < 1) // intersection is an ellipse { intr.type = ELLIPSE; intr.ellipse.normal = plane.normal; intr.ellipse.center = projectedCylinderOrigin  (d/cosTheta) * cylinder.direction; intr.ellipse.major = cylinder.direction  cosTheta * plane.normal; intr.ellipse.minor = Cross(plane.normal, intr.ellipse.major); intr.ellipse.majorLength = cylinder.radius / abs(cosTheta); intr.ellipse.minorLength = cylinder.radius; Normalize(intr.ellipse.major); Normalize(intr.ellipse.minor); return true; } else // abs(cosTheta) == 1, intersection is a circle { intr.type = CIRCLE; intr.circle.normal = plane.normal; intr.circle.center = projectedCylinderOrigin; intr.circle.radius = cylinder.radius; return true; } } else // abs(cosTheta) == 0, cylinder is parallel to plane { if (abs(d) < cylinder.radius) // cylinder intersects plane in two lines { Vector offset = Cross(cylinder.direction, plane.normal); float e = sqrt(cylinder.radius * cylinder.radius  d * d); intr.type = TWO_LINES; intr.line1.origin = projectedCylinderOrigin  e * offset; intr.line1.direction = cylinder.direction; intr.line2.origin = projectedCylinderOrigin + e * offset; intr.line2.direction = cylinder.direction; return true; } else if (abs(d) == cylinder.radius) // cylinder intersects plane in one line { intr.type = ONE_LINE; intr.line.origin = projectedCylinderOrigin; intr.line.direction = cylinder.direction; return true; } else // abs(d) > cylinder.radius, cylinder does not intersect plane { intr.type = EMPTY_SET; return false; } } } 
Section 11.7.4, pages 564569.
The first paragraph of Intersection Detection mentions a base point B
and refers you to Figure 11.14, which shows the base point labeled as C.
Later in Intersection Detection, there is a section Finite Cone
that refers to the base point using name B. This point is as shown in Figure
11.14 and is not the same as the apex A. On page 567, a vector hat(w) is defined. The author used the hat notation to denote unitlength vectors. However, the vector Cross(hat(d),Cross(hat(n),hat(d))) is not necessarily unit length. For example, if hat(n) is nearly parallel to hat(d), Cross(hat(n),hat(d)) is nearly the zero vector, in which case hat(w) is nearly the zero vector. The correct equation is hat(w) = Cross(hat(d),Cross(hat(n),hat(d)))/Length(Cross(hat(n),hat(d)))On pages 567, 568, and 569, the author has occurrences of arrow(w), but these should be hat(w). The variable a is the distance between points B and Ia, so the last equation on page 567 should be a = Length(Ia  B)That is, the subtraction of h is incorrect. Just to be clear, the variable b is the distance between points B and Ic, and the variable c is the distance between points Ia and Ic. On pages 568 and 569, any occurrence of Length(IaB)h should be replaced by Length(IaB). Observe that with the correct hat(w) equation, it must be that Dot(hat(n),hat(w)) = Length(Cross(hat(n),hat(d)), so it is not necessary to explicitly compute hat(w) in the source code. The equation for b^{2} is b^{2} = a^{2} [1/Length(Cross(hat(n),hat(d)))^{2}  1]However, notice that this equation has numerical problems when hat(n) and hat(d) are nearly parallel. The length of the cross product is nearly zero, so computing the reciprocal length is a problem. Instead, you want to code up the test as v = Cross(hat(n),hat(d)); lensqr = Dot(v,v); a*a*(1  lensqr) ≤ lensqr*r*r;Here is an alternative for the intersection test of plane and finite cone. Let the plane have unitlength normal N and contain point P. The plane equation is Dot(N,XP) = 0. The idea is to compute the interval of projection of the finite cone onto the line containing P with direction N. This inteval has endpoints P + t0*N and P + t1*N, where t0 < t1. If the interval [t0,t1] contains 0, then the cone and plane intersect. The projection interval endpoints are determined by one of the following.
X = C + h*D + h*tan(A)*(cos(B)*U + sin(B)*V)for B in [0,2*pi). The vectors U and V are unit length and {U,V,D} is an orthonormal set (vectors are unit length and mutually perpendicular). The projection of X onto the normal line is P + t(B)*N, where t(B) = Dot(N,XP) = Dot(N,CP) + h*Dot(N,D) + h*tan(A)*(cos(B)*Dot(N,U) + sin(B)*Dot(N,V))I have indicated that t is a function of B (different circle points project to different normal line points). The projection of C onto the normal line is Dot(N,CP). The minimum value of t(B) is t0 and the maximum value of t(B) is t1. Define d = Dot(N,CP) and r = sqrt(Dot(N,U)^2 + Dot(N,V)^2) = sqrt(1  Dot(N,D)^2)t0 is the minimum, and t0 is the maximum, of the set {d, d + h*Dot(N,D)  h*tan(A)*r, d + h*Dot(N,D) + h*tan(A)*r}The pseudocode is bool ConeIntersectsPlane (Cone cone, Plane plane) { Vector CmP = cone.C  plane.P; float d = Dot(plane.N, CmP); float NdotD = Dot(plane.N, cone.D); float r = sqrt(1  NdotD*NdotD); float tmp0 = d + cone.h*NdotD; float tmp1 = cone.h*cone.tanAngle*r; float vm = tmp0  tmp1; float vp = tmp0 + tmp1; float t0 = min(d, min(vm, vp)); float t1 = max(d, max(vm, vp)); return t0 <= 0 && 0 <= t1; }This algorithm does not require you to explicitly determine which of the bookmentioned cases 1, 2, 3a, or 3b applies. 
Section 11.7.4, page 571. The end of the last sentence of the first paragraph of the subsection Nondegenerate PlaneCone Intersections should be: cos(theta) = Dot(hat(a),hat(n)). 
Section 11.7.4, page 575.
The style of the pseudocode differs from other pages in this section. In the
first block starting with ``float d = '', remove the 'float' types (to be
consistent). In the section about the parabola, consistent notation would be
parabola.axis = Normalize(cone.axis  cosTheta * plane.normal); parabola.perpendicularAxis = Cross(parabola.axis, plane.normal); parabola.vertex = cone.vertex + d* plane.normal + e*parabola.axis; parabola.focalLength = d/2 * cosTheta;Also on this page, the first paragraph has a condition on the dot product of the plane normal hat(n) and the cone axis hat(a) in order that the two vectors be parallel. The condition should be: Dot(hat(n),hat(a) = 1 (for the true mathematical classification) or Dot(hat(n),hat(a)) >= 1epsilon (if you want to worry about numerical roundoff errors). 
Section 11.7.4, page 580.
The pseudocode has a line that computes the minor axis. It should be
minorAxis = Cross(plane.normal, majorAxis);because nowhere is there defined ``ellipse.u''. 
Section 11.11.1, page 612.
The last line of pseudocode is
if (value < min) min = value; else max = value;This should be if (value < min) min = value; else if (value > max) max = value; 
Section 11.12, page 628629. The pseudocode on page 628, the case where "Ray parallel to planes" has a test "P.x < box.xMax". That should be "P.x > box.xMax". The pseudocode on page 629 for swapping t0 and t1 is incorrect. It should be "tmp = t1; t1 = t0; t0 = tmp;" 
Section 11.12, pages 628629.
The pseudocode for rayAABB intersection has a loop that illustrates clipping
only against the planes x=xmin and x=xmax. This makes it unclear how the
the full set of tests fit in with the rest of the loop logic. The following
pseudocode is the replacement. The Point and Vector data structures are
assumed to have operator[] functions to access the coordinates. Index 0
corresponds to x (P[0] and P.x are equivalent), index 1 corresponds to y,
and index 2 corresponds to z. It is my opinion that this function should
also return the interval of intersection (first point of intersection
and last point of intersection). This change shows up in the last input to
the function. It is guaranteed that tIntersect[0] <= tIntersect[1]. If the
ray just grazes the AABB, we will return the same value in both parameters.
boolean RayIntersectAABB (Point P, Vector d, AABB box, float tIntersect[2]) { // The first part of this code tests for lineAABB intersection. tNear = INFINITY; tFar = INFINITY; for (i = 0; i < 3; i++) { // Check if the line is parallel to the min[i] and max[i] planes. if (abs(d[i]) < epsilon) { // The line is parallel to the min[i] and max[i] planes. if (P[i] < box.min[i]  P[i] > box.max[i]) { // The line is outside the slab defined by min[i] and max[i]. return false; } // else: The line is inside the slab defined by min[i] and max[i]. // Just continue testing with the next slab. continue; } // The line is not parallel to the slab planes, so compute the // parameters of intersection. The line segment of intersection is // P+t*d, where t0 <= t <= t1. t0 = (box.min[i]  P[i])/d[i]; t1 = (box.max[i]  P[i])/d[i]; if (t0 > t1) { tmp = t1; t1 = t0; t0 = tmp; } // Compare with current values. The current slab may have increased // tNear and/or decreased tFar. if (t0 > tNear) { tNear = t0; } if (t1 < tFar) { tFar = t1; } // Check if the line misses the AABB entirely. if (tNear > tFar) { return false; } } // At this place in the code, the line intersects the AABB in the line // segment P+t*d, tNear <= t <= tFar. If all you care about is lineAABB // intersection, then use the following three lines of code and omit the // remaining pseudocode. // // tIntersect[0] = tNear; // tIntersect[1] = tFar; // return true; // The ray is (P+t*d, t >= 0). We need to compute the intersection of // the interval tNear <= t <= tFar with the interval t >= 0. if (tFar < 0) { return false; } if (tNear < 0) { tNear = 0; } tIntersect[0] = tNear; tIntersect[1] = tFar; return true; } 
Section 11.12.3, page 633. The pseudocode for the rayOBB intersection is absurd. The first line of code is missing two right parentheses after "box.axis[i]". The inequality "r + box.halfLength[i] > 0" should be "r + box.halfLength[i] < 0". The rest of the pseudocode is worthless. The value s can very well be zero, leading to a division by zero when computing t0 and t1. If you find worthless pseudocode in the book, first let me know and, second, look at the code at my website. I wrote the real code (not the pseudocode), so if you find errors in the real code, please let me know. 
Section 11.12.4, page 635.
The return statement in the pseudocode is incorrect. It appears to be
designed for culling an AABB against the plane rather than indicating if the
AABB is intersected by the plane. The return statement also indicates that
the plane and AABB are considered not to be intersecting if the box
is just touching the plane at a single vertex or single edge or single face.
If the plane intersects the box at any point, I suggest replacing the return
statement by
return Dot(plane.normal,dMin)+plane.d <= 0 && Dot(plane.normal,dMax)+plane.d >= 0;I much prefer the "test intersection" that I use for planeOBB, which specializes easily when you have an AABB. No branching, just a few arithmetic calculations. Let N = (nx,ny,nz) be the plane normal and d be the plane constant. The plane is Dot(N,X)+d = 0. Let the AABB minima be x0, y0, and z0. Let the maxima be x1, y1, and z1. The test for intersection is nx*(x1+x0)+ny*(y1+y0)+nz*(z1+z0)+2*d <= nx*(x1x0)+ny*(y1y0)+nx*(z1z0) 
Section 11.12.7, pages 641 and 642. The last line of page 641 has a lowercase subscript 'a' on the unit vector u. That should be an uppercase 'A'. The last line of page 642 also has a lowercase subscript of 'a', but it should be an uppercase 'B'. 
Section 11.12.8, page 645.
The upper bounds on the displayed inequalities are wrong. They should be
P_{min,x} <= Q_x <= P_{max,x} P_{min,y} <= Q_y <= P_{max,y} P_{min,z} <= Q_z <= P_{max,z}The pseudocode on that same page is inconsistent about its reference to the sphere center. I changed "sphere.c" to "sphere.center" in two places. 
Section 11.12.10, page 660.
The equation for c0 is missing a power 2. It should be
c0 = (Dot(P,P)  (R^2+r^2))^2  4*R^2*(r^2  P_z^2) 
Section 13.1.4, page 686. The pseudocode has two sets named CoPos and CoNeg. These should be CoSame and CoDiff, respectively. 
Section 13.3.2, page 698.
The pseudocode for PointInSubpolygon has a call to GetMiddleIndex that is
missing its last parameter. The call should be
mid = GetMiddleIndex(i0, i1, C.N); 
Section 13.3.3, page 705. The second row of Table 13.1 has the four "Update" entries "i 0 m p". These entries should be "i 0 p m". 
Section 13.8.1, page 758759.
The input N is number of points, but the first mesh insert call returns a
new triangle also named N. These are, of course, different things. For
clarity, use Nj as the return value of the two "mesh.Insert" calls, one on
page 758 and one on page 759. Use Nj as the argument to the two
"stack.Push" calls, one on page 758 and one on page 759.
In the bottom block of pseudocode on page 759, the figure reference should be 13.46, not 13.45. After that line the pseudocode should be compute i0, i1, i2, i3 wwith T.v(i1) = A.v(i1) and T.v(i2) = A.v(i2);In the "if" clause, two lines must be added to the end of the clause mesh.Remove(T); mesh.Remove(A);After the two new triangles N0 and N1 are inserted into the mesh, the old triangles T and A must be removed. 
Section 13.9.1, pages 767768. The visibility graph does not contain an arc (i,i), where i is an index of a polygon vertex. The text should also state that adjacent polygon vertices are not considered to be visible to each other. Therefore, the visibility graph does not contain arcs (i,(i+1) modulo n), where n is the number of polygon vertices and 0 <= i <= n1. On page 768, the text states that the adjacency matrix for a convex polygon has 3 zero diagonals, the rest of the entries are ones. However, the entries at (row,column) of (0,n1) and (n1,0) are also zero since these correspond to adjacent polygon vertices. 
Section 13.9.1, page 769. The sentence that starts with "If < V_{i1}, V_{i+2} > is a diagonal" should have V_{i+1}, not V_{i+2}. 
Section 13.9.1, page 771.
In the SegmentInCode listing, the case for the "vertex is convex"
has the inequalities reversed. The line of code should be
return Kross(diff, edgeR) < 0 and Kross(diff, edgeL) > 0 
Section 13.9.3, page 777. The first sentence of the last paragraph has a subjectverb mismatch. It should be "An edge s_{k} has ...". In the same paragraph, the end of the second line has "s_{1}.y" but should be "s_{1}.y_{0}". 
Appendix A.1, page 831. The fourth line from the top of the page is "(1) x1  3*x2 + 2*x3 = 6". It should be "(1) x1  3*x2  2*x3 = 6". The last line of the page lists the solution as (1,3,2). It should be (1,3,2). 
Appendix A.1, page 831. The "09 January 2003" entry mentions that the first displayed equation "(1) x13*x2+2*x3=6" should be "(1) x13*x22x3=6". The same correction must be made in the next occurrence of this equation near the bottom of the page (after "which gives us a new system"). 
Appendix A.1, page 832.
The pseudocode at the top of the page has a loop "(for k = 1 to n)". This
should be "for (k = j+1 to n)". In that loop after the multiplier "m" is
computed, you could set a_{k,j} = 0 to make it clear you are eliminating
the column entries. If an implementation allows for the rowreduced matrix
to be used in further calculations, the assignment of zero should occur in
the code. But if the only goal is to solve the system, an implementation need
not assign zero since a_{k,j} is never used again in the code.
The first line of pseudocode for back substitution has "x_n = c_n/a_{n,h}". This should be "x_n = c_n/a_{n,n}". 
Appendix A.1, page 832.
In the elimination step, the line of code
a_{k,i} = a_{k,i} + m * a_{k,j}should be a_{k,i} = a_{k,i} + m * a_{j,i} 
Appendix A.2, page 834.
In the paragraph "The row reduction..." is "The condition for existence
of solutions is...". The equations should be corrected to
0 = a0*b1  a1*b0 = (c0+c1*y)*b1  a1*(d0+d1*y) = (b1*c0a1*d0) + (b1*c1a1*d1)*yLater is "as long as b1*c1b0*d1" which should be "as long as b1*c1a1*d1". Finally, the last equation of the last displayed equation array should be = (c0*f0e0*d0)+((c0*f1e0*d1)+(c1*f0e1*d0))*y + (c1*f1e1*d1)*y^2 
Appendix A.2, page 838. The last set of equations with the d[i][j] should have alpha[2][2] replaced by alpha[2][0] and beta[2][2] replaced by beta[2][0]. 
Appendix A.2, pages 844 and 846. The second line from the bottom on page 844 should have "e2 = 98". The fifth line from the bottom of page 846 should have "(1.11596,1.14945,1.63258)". 
Appendix A.4.3, page 861.
A bunch of errors are in "Quaternion to Matrix" and "Matrix to Quaternion".
This material was from "3D Game Engine Design", where the same errors occur,
but I had posted the correction for that book in February 2001. If R = [rij],
then equation (A.8) should have
r00 = 12y^22z^2 r01 = 2xy2wz r02 = 2xz+2wy r10 = 2xy+2wz r11 = 12x^22z^2 r12 = 2yz2wx r20 = 2xz2wy r21 = 2yz+2wx r22 = 12x^22y^2The last sentence of the first paragraph of "Matrix to Quaternion" has y = (r20r02)/(4w)but it should be y = (r02r20)/(4w)In the second paragraph of that section. When r00 is the maximum, the book has w = (r12r21)/(4x)but it should be w = (r21r12)/(4x)When r11 is the maximum, the book has w = (r20r02)/(4y)but it should be w = (r02r20)/(4y)When r22 is the maximum, the book has w = (r01r10)/(4z)but it should be w = (r10r01)/(4z) 
Appendix A.5.1, page 872.
The first displayed equation is
max{a_0, 1+a_1, ..., 1+a_{n1}} = 1 + max{a_0, ..., a_{n1}}but the last equaltity should be an inequality max{a_0, 1+a_1, ..., 1+a_{n1}} ≤ 1 + max{a_0, ..., a_{n1}} 
Appendix A.7.7, page 888.
The last displayed equation is
E(hat(c)) = (sum_{i=0}^n hat(c) * vector(v))^2 = hat(c)^T M hat(c)which is the square of a sum, but it should be a sum of squares, E(hat(c)) = sum_{i=0}^n (hat(c) * vector(v))^2 = hat(c)^T M hat(c) 
Appendix A.7.8, page 889.
Line 5 of the first paragraph of section A.7.8 has
M = sum_{i=0}^2 vector(v) vector(v)^TThe upper limit should be n, M = sum_{i=0}^n vector(v) vector(v)^T 
Appendix B.3.3, page 941.
All the product formulas are incorrect. They should be
sin(a)*sin(b) = (cos(ab)  cos(a+b))/2 sin(a)*cos(b) = (sin(a+b) + sin(ab))/2 cos(a)*cos(b) = (cos(ab) + cos(a+b))/2The last formula in the book appears to be an attempt to cover all the cases, but the first three equations are all the cases. 