src/rt/o_mesh.c

#ifndef lint
static const char RCSid[] = "$Id: o_mesh.c,v 2.6 2003/06/20 00:25:50 greg Exp $";
#endif
/*
 *  Routines for computing ray intersections with meshes.
 *
 *  Intersection with a triangle mesh is based on Segura and Feito's
 *  WSCG 2001 paper, "Algorithms to Test Ray-Triangle Intersection,
 *  Comparative Study."  This method avoids additional storage
 *  requirements, floating divides, and allows some savings by
 *  caching ray-edge comparisons that are otherwise repeated locally
 *  in typical mesh geometries.  (This is our own optimization.)
 *
 *  The code herein is quite similar to that in o_instance.c, the
 *  chief differences being the custom triangle intersection routines
 *  and the fact that an "OBJECT" in the mesh octree is not an index
 *  into the Radiance OBJREC list, but a mesh triangle index.  We still
 *  utilize the standard octree traversal code by setting the hitf
 *  function pointer in the RAY struct to our custom mesh_hit() call.
 */

#include  "copyright.h"

#include  "ray.h"

#include  "mesh.h"

#include  "tmesh.h"


#define  EDGE_CACHE_SIZ         251     /* length of mesh edge cache */

#define  curmi                  (edge_cache.mi)
#define  curmsh                 (curmi->msh)


/* Cache of signed volumes for this ray and this mesh */
struct EdgeCache {
        OBJREC          *o;     /* mesh object */
        MESHINST        *mi;    /* current mesh instance */
        struct EdgeSide {
                int32   v1i, v2i;       /* vertex indices (lowest first) */
                short   signum;         /* signed volume */
        }               cache[EDGE_CACHE_SIZ];
}       edge_cache;


static void
prep_edge_cache(o)              /* get instance and clear edge cache */
OBJREC  *o;
{
                                        /* get mesh instance */
        edge_cache.mi = getmeshinst(edge_cache.o = o, IO_ALL);
                                        /* clear edge cache */
        bzero((void *)edge_cache.cache, sizeof(edge_cache.cache));
}


static int
signed_volume(r, v1, v2)        /* get signed volume for ray and edge */
register RAY    *r;
int32           v1, v2;
{
        int                             reversed = 0;
        register struct EdgeSide        *ecp;
        
        if (v1 > v2) {
                int32   t = v2; v2 = v1; v1 = t;
                reversed = 1;
        }
        ecp = &edge_cache.cache[((v2<<11 ^ v1) & 0x7fffffff) % EDGE_CACHE_SIZ];
        if (ecp->v1i != v1 || ecp->v2i != v2) {
                MESHVERT        tv1, tv2;       /* compute signed volume */
                FVECT           v2d;
                double          vol;
                if (!getmeshvert(&tv1, edge_cache.mi->msh, v1, MT_V) ||
                            !getmeshvert(&tv2, edge_cache.mi->msh, v2, MT_V))
                        objerror(edge_cache.o, INTERNAL,
                                "missing mesh vertex in signed_volume");
                VSUB(v2d, tv2.v, r->rorg);
                vol = (tv1.v[0] - r->rorg[0]) *
                                (v2d[1]*r->rdir[2] - v2d[2]*r->rdir[1]);
                vol += (tv1.v[1] - r->rorg[1]) *
                                (v2d[2]*r->rdir[0] - v2d[0]*r->rdir[2]);
                vol += (tv1.v[2] - r->rorg[2]) *
                                (v2d[0]*r->rdir[1] - v2d[1]*r->rdir[0]);
                if (vol > .0)
                        ecp->signum = 1;
                else if (vol < .0)
                        ecp->signum = -1;
                else
                        ecp->signum = 0;
                ecp->v1i = v1;
                ecp->v2i = v2;
        }
        return(reversed ? -ecp->signum : ecp->signum);
}


static void
mesh_hit(oset, r)               /* intersect ray with mesh triangle(s) */
OBJECT  *oset;
RAY     *r;
{
        int32           tvi[3];
        int             sv1, sv2, sv3;
        MESHVERT        tv[3];
        OBJECT          tmod;
        FVECT           va, vb, nrm;
        double          d;
        int             i;
                                        /* check each triangle */
        for (i = oset[0]; i > 0; i--) {
                if (!getmeshtrivid(tvi, &tmod, curmsh, oset[i]))
                        objerror(edge_cache.o, INTERNAL,
                                "missing triangle vertices in mesh_hit");
                sv1 = signed_volume(r, tvi[0], tvi[1]);
                sv2 = signed_volume(r, tvi[1], tvi[2]);
                sv3 = signed_volume(r, tvi[2], tvi[0]);
                if (sv1 != sv2 || sv2 != sv3)   /* compare volume signs */
                        if (sv1 && sv2 && sv3)
                                continue;
                                                /* compute intersection */
                getmeshvert(&tv[0], curmsh, tvi[0], MT_V);
                getmeshvert(&tv[1], curmsh, tvi[1], MT_V);
                getmeshvert(&tv[2], curmsh, tvi[2], MT_V);
                VSUB(va, tv[0].v, tv[2].v);
                VSUB(vb, tv[1].v, tv[0].v);
                VCROSS(nrm, va, vb);
                d = DOT(r->rdir, nrm);
                if (d == 0.0)
                        continue;               /* ray is tangent */
                VSUB(va, tv[0].v, r->rorg);
                d = DOT(va, nrm) / d;
                if (d <= FTINY || d >= r->rot)
                        continue;               /* not good enough */
                r->robj = oset[i];              /* else record hit */
                r->ro = edge_cache.o;
                r->rot = d;
                VSUM(r->rop, r->rorg, r->rdir, d);
                VCOPY(r->ron, nrm);
                /* normalize(r->ron) called & r->rod set in o_mesh() */
        }
}


int
o_mesh(o, r)                    /* compute ray intersection with a mesh */
OBJREC          *o;
register RAY    *r;
{
        RAY             rcont;
        int             flags;
        MESHVERT        tv[3];
        OBJECT          tmod;
        RREAL           wt[3];
        int             i;
                                        /* get the mesh instance */
        prep_edge_cache(o);
                                        /* copy and transform ray */
        copystruct(&rcont, r);
        multp3(rcont.rorg, r->rorg, curmi->x.b.xfm);
        multv3(rcont.rdir, r->rdir, curmi->x.b.xfm);
        for (i = 0; i < 3; i++)
                rcont.rdir[i] /= curmi->x.b.sca;
        rcont.rmax *= curmi->x.b.sca;
                                        /* clear and trace ray */
        rayclear(&rcont);
        rcont.hitf = mesh_hit;
        if (!localhit(&rcont, &curmi->msh->mcube))
                return(0);                      /* missed */
        if (rcont.rot * curmi->x.f.sca >= r->rot)
                return(0);                      /* not close enough */
                                        /* transform ray back */
        r->rot = rcont.rot * curmi->x.f.sca;
        multp3(r->rop, rcont.rop, curmi->x.f.xfm);
        multv3(r->ron, rcont.ron, curmi->x.f.xfm);
        normalize(r->ron);
        r->rod = -DOT(r->rdir, r->ron);
                                        /* get triangle */
        flags = getmeshtri(tv, &tmod, curmsh, rcont.robj, MT_ALL);
        if (!(flags & MT_V))
                objerror(o, INTERNAL, "missing mesh vertices in o_mesh");
        r->robj = objndx(o);            /* set object and material */
        if (o->omod == OVOID && tmod != OVOID) {
                r->ro = getmeshpseudo(curmsh, tmod);
                r->rox = &curmi->x;
        } else
                r->ro = o;
                                        /* compute barycentric weights */
        if (flags & (MT_N|MT_UV))
                if (get_baryc(wt, rcont.rop, tv[0].v, tv[1].v, tv[2].v) < 0) {
                        objerror(o, WARNING, "bad triangle in o_mesh");
                        flags &= ~(MT_N|MT_UV);
                }
        if (flags & MT_N) {             /* interpolate normal */
                for (i = 0; i < 3; i++)
                        rcont.pert[i] = wt[0]*tv[0].n[i] +
                                        wt[1]*tv[1].n[i] +
                                        wt[2]*tv[2].n[i];
                multv3(r->pert, rcont.pert, curmi->x.f.xfm);
                if (normalize(r->pert) != 0.0)
                        for (i = 0; i < 3; i++)
                                r->pert[i] -= r->ron[i];
        } else
                r->pert[0] = r->pert[1] = r->pert[2] = .0;

        if (flags & MT_UV)              /* interpolate uv coordinates */
                for (i = 0; i < 2; i++)
                        r->uv[i] = wt[0]*tv[0].uv[i] +
                                        wt[1]*tv[1].uv[i] +
                                        wt[2]*tv[2].uv[i];
        else
                r->uv[0] = r->uv[1] = .0;

                                        /* return hit */
        return(1);
}
Revision:	2.7
Committed:	Thu Jun 26 00:58:10 2003 UTC (20 years, 10 months ago) by schorsch
Content type:	text/plain
Branch:	MAIN
Changes since 2.6:	+2 -2 lines
Log Message:	Abstracted process and path handling for Windows. Renamed FLOAT to RREAL because of conflict on Windows. Added conditional compiles for some signal handlers.
#	User	Rev	Content
1	greg	2.1	#ifndef lint
2	schorsch	2.7	static const char RCSid[] = "$Id: o_mesh.c,v 2.6 2003/06/20 00:25:50 greg Exp $";
3	greg	2.1	#endif
4			/*
5			* Routines for computing ray intersections with meshes.
6			*
7			* Intersection with a triangle mesh is based on Segura and Feito's
8			* WSCG 2001 paper, "Algorithms to Test Ray-Triangle Intersection,
9			* Comparative Study." This method avoids additional storage
10			* requirements, floating divides, and allows some savings by
11			* caching ray-edge comparisons that are otherwise repeated locally
12			* in typical mesh geometries. (This is our own optimization.)
13			*
14			* The code herein is quite similar to that in o_instance.c, the
15			* chief differences being the custom triangle intersection routines
16			* and the fact that an "OBJECT" in the mesh octree is not an index
17			* into the Radiance OBJREC list, but a mesh triangle index. We still
18			* utilize the standard octree traversal code by setting the hitf
19			* function pointer in the RAY struct to our custom mesh_hit() call.
20			*/
21
22			#include "copyright.h"
23
24			#include "ray.h"
25
26			#include "mesh.h"
27
28			#include "tmesh.h"
29
30
31	greg	2.2	#define EDGE_CACHE_SIZ 251 /* length of mesh edge cache */
32	greg	2.1
33			#define curmi (edge_cache.mi)
34			#define curmsh (curmi->msh)
35
36
37			/* Cache of signed volumes for this ray and this mesh */
38			struct EdgeCache {
39			OBJREC o; / mesh object */
40			MESHINST mi; / current mesh instance */
41			struct EdgeSide {
42	greg	2.6	int32 v1i, v2i; /* vertex indices (lowest first) */
43	greg	2.1	short signum; /* signed volume */
44			} cache[EDGE_CACHE_SIZ];
45			} edge_cache;
46
47
48			static void
49			prep_edge_cache(o) /* get instance and clear edge cache */
50			OBJREC *o;
51			{
52			/* get mesh instance */
53			edge_cache.mi = getmeshinst(edge_cache.o = o, IO_ALL);
54			/* clear edge cache */
55			bzero((void *)edge_cache.cache, sizeof(edge_cache.cache));
56			}
57
58
59			static int
60			signed_volume(r, v1, v2) /* get signed volume for ray and edge */
61			register RAY *r;
62	greg	2.6	int32 v1, v2;
63	greg	2.1	{
64			int reversed = 0;
65			register struct EdgeSide *ecp;
66
67			if (v1 > v2) {
68	greg	2.6	int32 t = v2; v2 = v1; v1 = t;
69	greg	2.1	reversed = 1;
70			}
71	greg	2.2	ecp = &edge_cache.cache[((v2<<11 ^ v1) & 0x7fffffff) % EDGE_CACHE_SIZ];
72	greg	2.1	if (ecp->v1i != v1 \|\| ecp->v2i != v2) {
73			MESHVERT tv1, tv2; /* compute signed volume */
74	greg	2.4	FVECT v2d;
75	greg	2.1	double vol;
76			if (!getmeshvert(&tv1, edge_cache.mi->msh, v1, MT_V) \|\|
77	greg	2.2	!getmeshvert(&tv2, edge_cache.mi->msh, v2, MT_V))
78	greg	2.1	objerror(edge_cache.o, INTERNAL,
79	greg	2.2	"missing mesh vertex in signed_volume");
80	greg	2.4	VSUB(v2d, tv2.v, r->rorg);
81	greg	2.1	vol = (tv1.v[0] - r->rorg[0]) *
82	greg	2.4	(v2d[1]r->rdir[2] - v2d[2]r->rdir[1]);
83	greg	2.1	vol += (tv1.v[1] - r->rorg[1]) *
84	greg	2.4	(v2d[2]r->rdir[0] - v2d[0]r->rdir[2]);
85	greg	2.1	vol += (tv1.v[2] - r->rorg[2]) *
86	greg	2.4	(v2d[0]r->rdir[1] - v2d[1]r->rdir[0]);
87	greg	2.1	if (vol > .0)
88			ecp->signum = 1;
89			else if (vol < .0)
90			ecp->signum = -1;
91			else
92			ecp->signum = 0;
93			ecp->v1i = v1;
94			ecp->v2i = v2;
95			}
96			return(reversed ? -ecp->signum : ecp->signum);
97			}
98
99
100			static void
101			mesh_hit(oset, r) /* intersect ray with mesh triangle(s) */
102			OBJECT *oset;
103			RAY *r;
104			{
105	greg	2.6	int32 tvi[3];
106	greg	2.1	int sv1, sv2, sv3;
107			MESHVERT tv[3];
108	greg	2.5	OBJECT tmod;
109	greg	2.1	FVECT va, vb, nrm;
110			double d;
111			int i;
112			/* check each triangle */
113			for (i = oset[0]; i > 0; i--) {
114	greg	2.5	if (!getmeshtrivid(tvi, &tmod, curmsh, oset[i]))
115	greg	2.1	objerror(edge_cache.o, INTERNAL,
116	greg	2.2	"missing triangle vertices in mesh_hit");
117	greg	2.1	sv1 = signed_volume(r, tvi[0], tvi[1]);
118			sv2 = signed_volume(r, tvi[1], tvi[2]);
119			sv3 = signed_volume(r, tvi[2], tvi[0]);
120			if (sv1 != sv2 \|\| sv2 != sv3) /* compare volume signs */
121			if (sv1 && sv2 && sv3)
122			continue;
123			/* compute intersection */
124			getmeshvert(&tv[0], curmsh, tvi[0], MT_V);
125			getmeshvert(&tv[1], curmsh, tvi[1], MT_V);
126			getmeshvert(&tv[2], curmsh, tvi[2], MT_V);
127			VSUB(va, tv[0].v, tv[2].v);
128			VSUB(vb, tv[1].v, tv[0].v);
129			VCROSS(nrm, va, vb);
130			d = DOT(r->rdir, nrm);
131			if (d == 0.0)
132			continue; /* ray is tangent */
133			VSUB(va, tv[0].v, r->rorg);
134			d = DOT(va, nrm) / d;
135			if (d <= FTINY \|\| d >= r->rot)
136			continue; /* not good enough */
137			r->robj = oset[i]; /* else record hit */
138			r->ro = edge_cache.o;
139			r->rot = d;
140			VSUM(r->rop, r->rorg, r->rdir, d);
141			VCOPY(r->ron, nrm);
142	greg	2.2	/* normalize(r->ron) called & r->rod set in o_mesh() */
143	greg	2.1	}
144			}
145
146
147			int
148			o_mesh(o, r) /* compute ray intersection with a mesh */
149			OBJREC *o;
150			register RAY *r;
151			{
152			RAY rcont;
153			int flags;
154			MESHVERT tv[3];
155	greg	2.5	OBJECT tmod;
156	schorsch	2.7	RREAL wt[3];
157	greg	2.1	int i;
158			/* get the mesh instance */
159			prep_edge_cache(o);
160			/* copy and transform ray */
161			copystruct(&rcont, r);
162			multp3(rcont.rorg, r->rorg, curmi->x.b.xfm);
163			multv3(rcont.rdir, r->rdir, curmi->x.b.xfm);
164			for (i = 0; i < 3; i++)
165			rcont.rdir[i] /= curmi->x.b.sca;
166			rcont.rmax *= curmi->x.b.sca;
167			/* clear and trace ray */
168			rayclear(&rcont);
169			rcont.hitf = mesh_hit;
170			if (!localhit(&rcont, &curmi->msh->mcube))
171			return(0); /* missed */
172			if (rcont.rot * curmi->x.f.sca >= r->rot)
173			return(0); /* not close enough */
174			/* transform ray back */
175			r->rot = rcont.rot * curmi->x.f.sca;
176			multp3(r->rop, rcont.rop, curmi->x.f.xfm);
177			multv3(r->ron, rcont.ron, curmi->x.f.xfm);
178			normalize(r->ron);
179			r->rod = -DOT(r->rdir, r->ron);
180	greg	2.5	/* get triangle */
181			flags = getmeshtri(tv, &tmod, curmsh, rcont.robj, MT_ALL);
182	greg	2.1	if (!(flags & MT_V))
183			objerror(o, INTERNAL, "missing mesh vertices in o_mesh");
184	greg	2.5	r->robj = objndx(o); /* set object and material */
185			if (o->omod == OVOID && tmod != OVOID) {
186			r->ro = getmeshpseudo(curmsh, tmod);
187			r->rox = &curmi->x;
188			} else
189			r->ro = o;
190			/* compute barycentric weights */
191	greg	2.1	if (flags & (MT_N\|MT_UV))
192			if (get_baryc(wt, rcont.rop, tv[0].v, tv[1].v, tv[2].v) < 0) {
193			objerror(o, WARNING, "bad triangle in o_mesh");
194			flags &= ~(MT_N\|MT_UV);
195			}
196			if (flags & MT_N) { /* interpolate normal */
197			for (i = 0; i < 3; i++)
198			rcont.pert[i] = wt[0]*tv[0].n[i] +
199			wt[1]*tv[1].n[i] +
200			wt[2]*tv[2].n[i];
201			multv3(r->pert, rcont.pert, curmi->x.f.xfm);
202			if (normalize(r->pert) != 0.0)
203			for (i = 0; i < 3; i++)
204			r->pert[i] -= r->ron[i];
205	greg	2.3	} else
206			r->pert[0] = r->pert[1] = r->pert[2] = .0;
207
208	greg	2.1	if (flags & MT_UV) /* interpolate uv coordinates */
209			for (i = 0; i < 2; i++)
210			r->uv[i] = wt[0]*tv[0].uv[i] +
211			wt[1]*tv[1].uv[i] +
212			wt[2]*tv[2].uv[i];
213	greg	2.3	else
214			r->uv[0] = r->uv[1] = .0;
215	greg	2.1
216			/* return hit */
217			return(1);
218			}