ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/virtuals.c
(Generate patch)

Comparing ray/src/rt/virtuals.c (file contents):
Revision 1.3 by greg, Thu Jun 20 13:43:38 1991 UTC vs.
Revision 1.8 by greg, Tue Jun 25 08:52:04 1991 UTC

# Line 11 | Line 11 | static char SCCSid[] = "$SunId$ LBL";
11  
12   #include  "ray.h"
13  
14 + #include  "octree.h"
15 +
16   #include  "otypes.h"
17  
18   #include  "source.h"
19  
20 < #include  "cone.h"
20 > #include  "random.h"
21  
20 #include  "face.h"
22  
23 + double  getdisk();
24  
23 double  intercircle();
24 SRCREC  *makevsrc();
25
25   static OBJECT  *vobject;                /* virtual source objects */
26   static int  nvobjects = 0;              /* number of virtual source objects */
27  
# Line 55 | Line 54 | markvirtuals()                 /* find and mark virtual sources */
54          }
55          if (nvobjects == 0)
56                  return;
57 + #ifdef DEBUG
58 +        fprintf(stderr, "found %d virtual source objects\n", nvobjects);
59 + #endif
60                                          /* append virtual sources */
61          for (i = nsources; i-- > 0; )
62 <                if (!(source[i].sflags & SSKIP))
61 <                        addvirtuals(&source[i], directrelay);
62 >                addvirtuals(i, directrelay);
63                                          /* done with our object list */
64          free((char *)vobject);
65          nvobjects = 0;
66   }
67  
68  
69 < addvirtuals(sr, nr)             /* add virtual sources associated with sr */
70 < SRCREC  *sr;
69 > addvirtuals(sn, nr)             /* add virtuals associated with source */
70 > int  sn;
71   int  nr;
72   {
73          register int  i;
74                                  /* check relay limit first */
75          if (nr <= 0)
76                  return;
77 +        if (source[sn].sflags & SSKIP)
78 +                return;
79                                  /* check each virtual object for projection */
80          for (i = 0; i < nvobjects; i++)
81                                          /* vproject() calls us recursively */
82 <                vproject(objptr(i), sr, nr-1);
82 >                vproject(objptr(vobject[i]), sn, nr-1);
83   }
84  
85  
86 < vproject(o, s, n)               /* create projected source(s) if they exist */
86 > vproject(o, sn, n)              /* create projected source(s) if they exist */
87   OBJREC  *o;
88 < SRCREC  *s;
88 > int  sn;
89   int  n;
90   {
91          register int  i;
92          register VSMATERIAL  *vsmat;
93          MAT4  proj;
94 <        SRCREC  *ns;
94 >        int  ns;
95 >
96 >        if (o == source[sn].so) /* objects cannot project themselves */
97 >                return;
98                                  /* get virtual source material */
99          vsmat = sfun[objptr(o->omod)->otype].mf;
100                                  /* project virtual sources */
101          for (i = 0; i < vsmat->nproj; i++)
102 <                if ((*vsmat->vproj)(proj, o, s, i))
103 <                        if ((ns = makevsrc(o, s, proj)) != NULL)
102 >                if ((*vsmat->vproj)(proj, o, &source[sn], i))
103 >                        if ((ns = makevsrc(o, sn, proj)) >= 0) {
104 > #ifdef DEBUG
105 >                                virtverb(ns, stderr);
106 > #endif
107                                  addvirtuals(ns, n);
108 +                        }
109   }
110  
111  
112 < SRCREC *
113 < makevsrc(op, sp, pm)            /* make virtual source if reasonable */
112 > int
113 > makevsrc(op, sn, pm)            /* make virtual source if reasonable */
114   OBJREC  *op;
115 < register SRCREC  *sp;
115 > register int  sn;
116   MAT4  pm;
117   {
118 <        register SRCREC  *newsrc;
109 <        FVECT  nsloc, ocent, nsnorm;
110 <        int  nsflags;
118 >        FVECT  nsloc, nsnorm, ocent;
119          double  maxrad2;
120 +        int  nsflags;
121          double  d1;
122          SPOT  theirspot, ourspot;
123          register int  i;
124  
125 <        nsflags = (sp->sflags|(SVIRTUAL|SFOLLOW)) & ~SSPOT;
125 >        nsflags = source[sn].sflags | (SVIRTUAL|SSPOT|SFOLLOW);
126                                          /* get object center and max. radius */
127 <        if (sfun[op->otype].of->getdisk != NULL) {
128 <                maxrad2 = (*sfun[op->otype].of->getdisk)(ocent, op);
129 <                if (maxrad2 <= FTINY)                   /* too small? */
121 <                        return(NULL);
122 <                nsflags |= SSPOT;
123 <        }
127 >        maxrad2 = getdisk(ocent, op, sn);
128 >        if (maxrad2 <= FTINY)                   /* too small? */
129 >                return(-1);
130                                          /* get location and spot */
131 <        if (sp->sflags & SDISTANT) {            /* distant source */
132 <                if (sp->sflags & SPROX)
133 <                        return(NULL);           /* should never get here! */
134 <                multv3(nsloc, sp->sloc, pm);
135 <                if (nsflags & SSPOT) {
136 <                        VCOPY(ourspot.aim, ocent);
137 <                        ourspot.siz = PI*maxrad2;
138 <                        ourspot.flen = 0.;
131 >        if (source[sn].sflags & SDISTANT) {             /* distant source */
132 >                if (source[sn].sflags & SPROX)
133 >                        return(-1);             /* should never get here! */
134 >                multv3(nsloc, source[sn].sloc, pm);
135 >                VCOPY(ourspot.aim, ocent);
136 >                ourspot.siz = PI*maxrad2;
137 >                ourspot.flen = 0.;
138 >                if (source[sn].sflags & SSPOT) {
139 >                        copystruct(&theirspot, source[sn].sl.s);
140 >                        multp3(theirspot.aim, source[sn].sl.s->aim, pm);
141 >                        if (!commonbeam(&ourspot, &theirspot, nsloc))
142 >                                return(-1);             /* no overlap */
143                  }
134                if (sp->sflags & SSPOT) {
135                        copystruct(&theirspot, sp->sl.s);
136                        multp3(theirspot.aim, sp->sl.s->aim, pm);
137                        if (nsflags & SSPOT &&
138                                !commonbeam(&ourspot, &theirspot, nsloc))
139                                return(NULL);           /* no overlap */
140                }
144          } else {                                /* local source */
145 <                multp3(nsloc, sp->sloc, pm);
146 <                if (nsflags & SSPOT) {
147 <                        for (i = 0; i < 3; i++)
148 <                                ourspot.aim[i] = ocent[i] - nsloc[i];
149 <                        if ((d1 = normalize(ourspot.aim)) == 0.)
150 <                                return(NULL);           /* at source!! */
151 <                        if (sp->sflags & SPROX && d1 > sp->sl.prox)
152 <                                return(NULL);           /* too far away */
153 <                        ourspot.siz = 2.*PI*(1. - d1/sqrt(d1*d1+maxrad2));
154 <                        ourspot.flen = 0.;
155 <                } else if (sp->sflags & SPROX) {
156 <                        FVECT  norm;
157 <                        double  offs;
158 <                                                /* use distance from plane */
159 <                        offs = (*sfun[op->otype].of->getpleq)(norm, op);
157 <                        d1 = DOT(norm, nsloc) - offs;
158 <                        if (d1 > sp->sl.prox || d1 < -sp->sl.prox)
159 <                                return(NULL);           /* too far away */
145 >                multp3(nsloc, source[sn].sloc, pm);
146 >                for (i = 0; i < 3; i++)
147 >                        ourspot.aim[i] = ocent[i] - nsloc[i];
148 >                if ((d1 = normalize(ourspot.aim)) == 0.)
149 >                        return(-1);             /* at source!! */
150 >                if (source[sn].sflags & SPROX && d1 > source[sn].sl.prox)
151 >                        return(-1);             /* too far away */
152 >                ourspot.siz = 2.*PI*(1. - d1/sqrt(d1*d1+maxrad2));
153 >                ourspot.flen = 0.;
154 >                if (source[sn].sflags & SSPOT) {
155 >                        copystruct(&theirspot, source[sn].sl.s);
156 >                        multv3(theirspot.aim, source[sn].sl.s->aim, pm);
157 >                        if (!commonspot(&ourspot, &theirspot, nsloc))
158 >                                return(-1);     /* no overlap */
159 >                        ourspot.flen = theirspot.flen;
160                  }
161 <                if (sp->sflags & SSPOT) {
162 <                        copystruct(&theirspot, sp->sl.s);
163 <                        multv3(theirspot.aim, sp->sl.s->aim, pm);
164 <                        if (nsflags & SSPOT) {
165 <                                if (!commonspot(&ourspot, &theirspot, nsloc))
166 <                                        return(NULL);   /* no overlap */
167 <                                ourspot.flen = theirspot.flen;
168 <                        }
161 >                if (source[sn].sflags & SFLAT) {        /* behind source? */
162 >                        multv3(nsnorm, source[sn].snorm, pm);
163 >                        if (checkspot(&ourspot, nsnorm) < 0)
164 >                                return(-1);
165                  }
170                if (sp->sflags & SFLAT) {       /* check for behind source */
171                        multv3(nsnorm, sp->snorm, pm);
172                        if (nsflags & SSPOT && checkspot(&ourspot, nsnorm) < 0)
173                                return(NULL);
174                }
166          }
167 <                                        /* everything is OK, make source */
168 <        if ((newsrc = newsource()) == NULL)
167 >                                        /* pretest visibility */
168 >        nsflags = vstestvis(nsflags, op, ocent, maxrad2, sn);
169 >        if (nsflags & SSKIP)
170 >                return(-1);     /* obstructed */
171 >                                        /* it all checks out, so make it */
172 >        if ((i = newsource()) < 0)
173                  goto memerr;
174 <        newsrc->sflags = nsflags;
175 <        VCOPY(newsrc->sloc, nsloc);
174 >        source[i].sflags = nsflags;
175 >        VCOPY(source[i].sloc, nsloc);
176          if (nsflags & SFLAT)
177 <                VCOPY(newsrc->snorm, nsnorm);
178 <        newsrc->ss = sp->ss; newsrc->ss2 = sp->ss2;
179 <        if ((nsflags | sp->sflags) & SSPOT) {
180 <                if ((newsrc->sl.s = (SPOT *)malloc(sizeof(SPOT))) == NULL)
181 <                        goto memerr;
187 <                if (nsflags & SSPOT)
188 <                        copystruct(newsrc->sl.s, &ourspot);
189 <                else
190 <                        copystruct(newsrc->sl.s, &theirspot);
191 <                newsrc->sflags |= SSPOT;
192 <        }
177 >                VCOPY(source[i].snorm, nsnorm);
178 >        source[i].ss = source[sn].ss; source[i].ss2 = source[sn].ss2;
179 >        if ((source[i].sl.s = (SPOT *)malloc(sizeof(SPOT))) == NULL)
180 >                goto memerr;
181 >        copystruct(source[i].sl.s, &ourspot);
182          if (nsflags & SPROX)
183 <                newsrc->sl.prox = sp->sl.prox;
184 <        newsrc->sa.svnext = sp - source;
185 <        return(newsrc);
183 >                source[i].sl.prox = source[sn].sl.prox;
184 >        source[i].sa.svnext = sn;
185 >        source[i].so = op;
186 >        return(i);
187   memerr:
188          error(SYSTEM, "out of memory in makevsrc");
189   }
190  
191  
192 < commonspot(sp1, sp2, org)       /* set sp1 to intersection of sp1 and sp2 */
193 < register SPOT  *sp1, *sp2;
194 < FVECT  org;
192 > double
193 > getdisk(oc, op, sn)             /* get visible object disk */
194 > FVECT  oc;
195 > OBJREC  *op;
196 > register int  sn;
197   {
198 <        FVECT  cent;
199 <        double  rad2, cos1, cos2;
200 <
201 <        cos1 = 1. - sp1->siz/(2.*PI);
202 <        cos2 = 1. - sp2->siz/(2.*PI);
203 <        if (sp2->siz >= 2.*PI-FTINY)            /* BIG, just check overlap */
204 <                return(DOT(sp1->aim,sp2->aim) >= cos1*cos2 -
205 <                                        sqrt((1.-cos1*cos1)*(1.-cos2*cos2)));
206 <                                /* compute and check disks */
207 <        rad2 = intercircle(cent, sp1->aim, sp2->aim,
208 <                        1./(cos1*cos1) - 1.,  1./(cos2*cos2) - 1.);
209 <        if (rad2 <= FTINY || normalize(cent) == 0.)
210 <                return(0);
211 <        VCOPY(sp1->aim, cent);
212 <        sp1->siz = 2.*PI*(1. - 1./sqrt(1.+rad2));
213 <        return(1);
198 >        double  rad2, roffs, offs, d, rd, rdoto;
199 >        FVECT  rnrm, nrm;
200 >                                /* first, use object getdisk function */
201 >        rad2 = (*sfun[op->otype].of->getdisk)(oc, op);
202 >        if (!(source[sn].sflags & SVIRTUAL))
203 >                return(rad2);           /* all done for normal source */
204 >                                /* check for correct side of relay surface */
205 >        roffs = (*sfun[source[sn].so->otype].of->getpleq)(rnrm, source[sn].so);
206 >        rd = DOT(rnrm, source[sn].sloc);        /* source projection */
207 >        if (!(source[sn].sflags & SDISTANT))
208 >                rd -= roffs;
209 >        d = DOT(rnrm, oc) - roffs;      /* disk distance to relay plane */
210 >        if ((d > 0.) ^ (rd > 0.))
211 >                return(rad2);           /* OK if opposite sides */
212 >        if (d*d >= rad2)
213 >                return(.0);             /* no relay is possible */
214 >                                /* we need a closer look */
215 >        offs = (*sfun[op->otype].of->getpleq)(nrm, op);
216 >        rdoto = DOT(rnrm, nrm);
217 >        if (d*d >= rad2*(1.-rdoto*rdoto))
218 >                return(0.);             /* disk entirely on projection side */
219 >                                /* should shrink disk but I'm lazy */
220 >        return(rad2);
221   }
222  
223  
224 < commonbeam(sp1, sp2, dir)       /* set sp1 to intersection of sp1 and sp2 */
225 < register SPOT  *sp1, *sp2;
226 < FVECT  dir;
224 > int
225 > vstestvis(f, o, oc, or2, sn)            /* pretest source visibility */
226 > int  f;                 /* virtual source flags */
227 > OBJREC  *o;             /* relay object */
228 > FVECT  oc;              /* relay object center */
229 > double  or2;            /* relay object radius squared */
230 > register int  sn;       /* target source number */
231   {
232 <        FVECT  cent, c1, c2;
233 <        double  rad2, d;
234 <        register int  i;
235 <                                        /* move centers to common plane */
236 <        d = DOT(sp1->aim, dir);
237 <        for (i = 0; i < 3; i++)
238 <                c1[i] = sp1->aim[i] - d*dir[i];
239 <        d = DOT(sp2->aim, dir);
240 <        for (i = 0; i < 3; i++)
241 <                c2[i] = sp2->aim[i] - d*dir[i];
242 <                                        /* compute overlap */
243 <        rad2 = intercircle(cent, c1, c2, sp1->siz/PI, sp2->siz/PI);
244 <        if (rad2 <= FTINY)
245 <                return(0);
246 <        VCOPY(sp1->aim, cent);
247 <        sp1->siz = PI*rad2;
248 <        return(1);
232 >        RAY  sr;
233 >        FVECT  onorm;
234 >        FVECT  offsdir;
235 >        double  or, d;
236 >        int  infront;
237 >        int  ssn;
238 >        int  nok, nhit;
239 >        register int  i, n;
240 >                                /* return if pretesting disabled */
241 >        if (vspretest <= 0)
242 >                return(f);
243 >                                /* get surface normal */
244 >        (*sfun[o->otype].of->getpleq)(onorm, o);
245 >                                /* set number of rays to sample */
246 >        if (source[sn].sflags & SDISTANT) {
247 >                n = (2./3.*PI*PI)*or2/(thescene.cusize*thescene.cusize)*
248 >                                vspretest + .5;
249 >                infront = DOT(onorm, source[sn].sloc) > 0.;
250 >        } else {
251 >                n = or2/dist2(oc,source[sn].sloc)*vspretest + .5;
252 >                for (i = 0; i < 3; i++)
253 >                        offsdir[i] = source[sn].sloc[i] - oc[i];
254 >                infront = DOT(onorm, offsdir) > 0.;
255 >        }
256 >        if (n < 1) n = 1;
257 >                                /* sample */
258 >        or = sqrt(or2);
259 >        ssn = 7*n;
260 >        nhit = nok = 0;
261 >        while (n-- > 0) {
262 >                                        /* get sample point */
263 >                do {
264 >                        if (--ssn < 0)
265 >                                return(f);      /* too small a target! */
266 >                        for (i = 0; i < 3; i++)
267 >                                offsdir[i] = or*(1. -
268 >                                                2.*urand(931*i+5827+ssn));
269 >                        for (i = 0; i < 3; i++)
270 >                                sr.rorg[i] = oc[i] + offsdir[i];
271 >                        d = DOT(offsdir,onorm);
272 >                        if (infront)
273 >                                for (i = 0; i < 3; i++) {
274 >                                        sr.rorg[i] -= (d-.0001)*onorm[i];
275 >                                        sr.rdir[i] = -onorm[i];
276 >                                }
277 >                        else
278 >                                for (i = 0; i < 3; i++) {
279 >                                        sr.rorg[i] -= (d+.0001)*onorm[i];
280 >                                        sr.rdir[i] = onorm[i];
281 >                                }
282 >                        rayorigin(&sr, NULL, PRIMARY, 1.0);
283 >                } while (!(*ofun[o->otype].funp)(o, &sr));
284 >                                        /* check against source */
285 >                samplendx++;
286 >                if (srcray(&sr, NULL, sn) == 0.0)
287 >                        continue;
288 >                sr.revf = srcvalue;
289 >                rayvalue(&sr);
290 >                if (bright(sr.rcol) <= FTINY)
291 >                        continue;
292 >                nok++;
293 >                                        /* check against obstructions */
294 >                srcray(&sr, NULL, sn);
295 >                rayvalue(&sr);
296 >                if (bright(sr.rcol) <= FTINY)
297 >                        continue;
298 >                nhit++;
299 >        }
300 >                                /* interpret results */
301 >        if (nhit == 0)
302 >                return(f | SSKIP);      /* 0% hit rate:  totally occluded */
303 >        if (nhit == nok)
304 >                return(f & ~SFOLLOW);   /* 100% hit rate:  no occlusion */
305 >        return(f);              /* no comment */
306   }
307 +        
308  
309 <
310 < checkspot(sp, nrm)              /* check spotlight for behind source */
311 < register SPOT  *sp;
312 < FVECT  nrm;
309 > #ifdef DEBUG
310 > virtverb(sn, fp)        /* print verbose description of virtual source */
311 > register int  sn;
312 > FILE  *fp;
313   {
253        double  d, d1;
254
255        d = DOT(sp->aim, nrm);
256        if (d > FTINY)                  /* center in front? */
257                return(0);
258                                        /* else check horizon */
259        d1 = 1. - sp->siz/(2.*PI);
260        return(1.-FTINY-d*d > d1*d1);
261 }
262
263
264 double
265 intercircle(cc, c1, c2, r1s, r2s)       /* intersect two circles */
266 FVECT  cc;                      /* midpoint (return value) */
267 FVECT  c1, c2;                  /* circle centers */
268 double  r1s, r2s;               /* radii squared */
269 {
270        double  a2, d2, l;
271        FVECT  disp;
314          register int  i;
315  
316 <        for (i = 0; i < 3; i++)
317 <                disp[i] = c2[i] - c1[i];
318 <        d2 = DOT(disp,disp);
319 <                                        /* circle within overlap? */
320 <        if (r1s < r2s) {
321 <                if (r2s >= r1s + d2) {
322 <                        VCOPY(cc, c1);
323 <                        return(r1s);
324 <                }
325 <        } else {
326 <                if (r1s >= r2s + d2) {
327 <                        VCOPY(cc, c2);
328 <                        return(r2s);
329 <                }
330 <        }
331 <        a2 = .25*(2.*(r1s+r2s) - d2 - (r2s-r1s)*(r2s-r1s)/d2);
332 <                                        /* no overlap? */
291 <        if (a2 <= 0.)
292 <                return(0.);
293 <                                        /* overlap, compute center */
294 <        l = sqrt((r1s - a2)/d2);
295 <        for (i = 0; i < 3; i++)
296 <                cc[i] = c1[i] + l*disp[i];
297 <        return(a2);
316 >        fprintf(fp, "%s virtual source %d in %s %s\n",
317 >                        source[sn].sflags & SDISTANT ? "distant" : "local",
318 >                        sn, ofun[source[sn].so->otype].funame,
319 >                        source[sn].so->oname);
320 >        fprintf(fp, "\tat (%f,%f,%f)\n",
321 >                source[sn].sloc[0], source[sn].sloc[1], source[sn].sloc[2]);
322 >        fprintf(fp, "\tlinked to source %d (%s)\n",
323 >                source[sn].sa.svnext, source[source[sn].sa.svnext].so->oname);
324 >        if (source[sn].sflags & SFOLLOW)
325 >                fprintf(fp, "\talways followed\n");
326 >        else
327 >                fprintf(fp, "\tnever followed\n");
328 >        if (!(source[sn].sflags & SSPOT))
329 >                return;
330 >        fprintf(fp, "\twith spot aim (%f,%f,%f) and size %f\n",
331 >                        source[sn].sl.s->aim[0], source[sn].sl.s->aim[1],
332 >                        source[sn].sl.s->aim[2], source[sn].sl.s->siz);
333   }
334 + #endif

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines