File: ObjectCastImplB.mesa
Created January 2, 1983 12:17 am
Copyright © 1984 by Xerox Corporation. All rights reserved.
Last edited by Bier on February 18, 1987 6:29:23 pm PST
Contents: Code for casting rays at cones
DIRECTORY
Algebra3d, CastRays, ObjectCast, CSG, SV3d, SVRayTypes, SVVector3d;
ObjectCastImplB:
CEDAR PROGRAM
IMPORTS Algebra3d, CastRays, CSG
EXPORTS ObjectCast =
BEGIN
Classification: TYPE = SVRayTypes.Classification;
Point3d: TYPE = SV3d.Point3d;
Primitive: TYPE = SVRayTypes.Primitive;
Ray: TYPE = CSG.Ray;
SurfaceArray: TYPE = SVRayTypes.SurfaceArray;
Vector3d: TYPE = SV3d.Vector3d;
HitRec: TYPE = REF HitRecObj;
HitRecObj:
TYPE =
RECORD [
t: REAL,
normal: Vector3d];
HitArray: TYPE = REF HitArrayObj;
HitArrayObj: TYPE = ARRAY[1..2] OF HitRec;
globalPlaneHit: HitRec;
globalConeHitArray: HitArray;
globalConeHitCount: NAT;
almostZero: REAL = 1.0e-12;
d1d1, d2d2, d3d3: REAL;
p1p1, p2p2, p3p3: REAL;
p1d1, p2d2, p3d3: REAL;
xd2, zd2: REAL;
a, b, c: REAL;
x0, y0, z0: REAL;
normal: Vector3d;
PositiveRoots:
PROCEDURE [a, b, c:
REAL]
RETURNS [rootArray:
ARRAY[1..2]
OF
REAL, rootCount:
NAT] = {
Use only the positive roots.
i: NAT ← 1;
[rootArray, rootCount] ← Algebra3d.QuadraticFormula[a, b, c];
IF rootCount = 0 THEN RETURN;
WHILE i <= rootCount
DO
IF rootArray[i] <= 0
THEN {
FOR j:
NAT
IN [i..rootCount-1]
DO
rootArray[j] ← rootArray[j+1];
ENDLOOP;
rootCount ← rootCount - 1;
}
ELSE {i ← i + 1};
ENDLOOP;
};
ConeCast:
PUBLIC
PROC [localRay: Ray, surfaces: SurfaceArray, prim: Primitive]
RETURNS [class: Classification] = {
1 ray-disk test and 1 ray-cone test.
The disk is at y=0. radius = 1. Height = 1;
As usual, ray equation is:x = p[1] + t*d[1]
y = p[2] + t*d[2]
z = p[3] + t*d[3]
Disk test is thus p[2] + t*d[2] = 0; (t = -p[2]/d[2]) where d[2] # 0;
Then x = p[1] + (-p[2])*d[1]/d[2]
and z = p[3] + (-p[2])*d[3]/d[2]
x*d[2] = p[1]*d[2] -p[2]*d[1]; and z*d[2] = p[3]*d[2] -p[2]*d[3]
we require x^2 + z^2 <= 1;
ie (x*d[2])^2 + (z*d[2])^2 <= d[2]^2;
d: Vector3d;
p: Point3d;
hitDisk: BOOL ← FALSE;
planeT: REAL;
coneT1: REAL;
realRoots: ARRAY[1..2] OF REAL;
rootCount: NAT;
[p, d] ← CSG.GetLocalRay[localRay];
class ← CastRays.GetClassFromPool[];
IF Abs[d[2]] > almostZero
THEN {
xd2 ← p[1]*d[2] -p[2]*d[1];
zd2 ← p[3]*d[2] -p[2]*d[3];
IF xd2*xd2 + zd2*zd2 <= d[2]*d[2]
THEN {
planeT ← -p[2]/d[2];
IF planeT > 0
THEN {
hitDisk ← TRUE;
globalPlaneHit^ ← [planeT, [0,-1.0, 0] ];
};
};
};
Disk test complete. Do cone test.
I claim that the parameter t at a ray-cone intersection must satisfy the quadratic equation:
at^2 + bt + c = 0; where
a = d1d1 + d3d3 - d2d2.
b = 2*(p1d1 + p3d3 + d[2] - p2d2).
c = p1p1 + p3p3 - 1 + 2*p[2] + p2p2.
d1d1 ← d[1]*d[1];
d2d2 ← d[2]*d[2];
d3d3 ← d[3]*d[3];
p1p1 ← p[1]*p[1];
p2p2 ← p[2]*p[2];
p3p3 ← p[3]*p[3];
p1d1 ← p[1]*d[1];
p2d2 ← p[2]*d[2];
p3d3 ← p[3]*d[3];
a ← d1d1 + d3d3 - d2d2;
b ← 2*(p1d1 + p3d3 + d[2] - p2d2);
c ← p1p1 + p3p3 - 1 + 2*p[2] - p2p2;
Solve this quadratic.
[realRoots, rootCount] ← PositiveRoots[a, b, c];
Find how many of the cone hits are above the disk and below the apex.
globalConeHitCount ← 0;
FOR i:
NAT
IN [1..rootCount]
DO
y0 ← p[2] + realRoots[i]*d[2];
IF y0 > 0.0
AND y0 <= 1
THEN {
Calculate the surface normal. It is just [x0, 1-y0, z0].
x0 ← p[1] + realRoots[i]*d[1];
z0 ← p[3] + realRoots[i]*d[3];
normal ← [x0, 1-y0, z0];
globalConeHitCount ← globalConeHitCount+1;
globalConeHitArray[globalConeHitCount]^ ← [realRoots[i], normal];
};
ENDLOOP;
SELECT globalConeHitCount
FROM
0 => {
-- ray missed cone completely or scraped an edge. Count as miss
class.count ← 0;
class.classifs[1] ← FALSE;
};
1 => {
We have a hit on the shroud. There should be a hit on the plane
as well. If not, assume that this is only a glancing hit and count as a miss.
coneT1 ← globalConeHitArray[1].t;
IF hitDisk
THEN {
class.count ← 2;
IF coneT1 < planeT
THEN {
class.params[1] ← coneT1;
class.params[2] ← planeT;
class.normals[1] ← globalConeHitArray[1].normal;
class.normals[2] ← globalPlaneHit.normal;
class.surfaces[1] ← surfaces[1];
class.surfaces[2] ← surfaces[2];}
ELSE {
class.params[2] ← coneT1;
class.params[1] ← planeT;
class.normals[2] ← globalConeHitArray[1].normal;
class.normals[1] ← globalPlaneHit.normal;
class.surfaces[2] ← surfaces[1];
class.surfaces[1] ← surfaces[2];
};
class.classifs[1] ← FALSE; class.classifs[2] ← TRUE; class.classifs[3] ← FALSE;
class.primitives[1] ← class.primitives[2] ← prim;
}
ELSE {
-- count as a miss
class.count ← 0;
class.classifs[1] ← FALSE;
};
}; -- end of one shroud hit.
2 => {
Two shroud hits. There should be no disk hits unless ray goes through an edge.
If so, ignore the disk hit. Hence don't even check for a disk hit.
class.count ← 2;
IF globalConeHitArray[1].t < globalConeHitArray[2].t
THEN {
class.params[1] ← globalConeHitArray[1].t;
class.params[2] ← globalConeHitArray[2].t;
class.normals[1] ← globalConeHitArray[1].normal;
class.normals[2] ← globalConeHitArray[2].normal}
ELSE {
class.params[2] ← globalConeHitArray[1].t;
class.params[1] ← globalConeHitArray[2].t;
class.normals[2] ← globalConeHitArray[1].normal;
class.normals[1] ← globalConeHitArray[2].normal};
class.surfaces[1] ← surfaces[1]; class.surfaces[2] ← surfaces[1];
class.classifs[1] ← FALSE; class.classifs[2] ← TRUE; class.classifs[3] ← FALSE;
class.primitives[1] ← class.primitives[2] ← prim;
};
3, 4 => {
I don't know how this happens. Count as a miss for now
class.count ← 0;
class.classifs[1] ← FALSE;
};
ENDCASE => ERROR;
}; -- end of ConeCast
Abs:
PROC [r:
REAL]
RETURNS [
REAL] = {
RETURN [IF r>=0 THEN r ELSE -r];
};
Init:
PROC = {
globalPlaneHit ← NEW[HitRecObj];
globalConeHitArray ← NEW[HitArrayObj];
FOR i:
NAT
IN[1..2]
DO
globalConeHitArray[i] ← NEW[HitRecObj];
ENDLOOP;
};
Init[];
END.