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

Comparing ray/src/rt/srcsamp.c (file contents):
Revision 2.17 by greg, Sat Oct 23 18:13:54 2010 UTC vs.
Revision 2.18 by greg, Wed Dec 28 18:39:36 2011 UTC

# Line 16 | Line 16 | static const char      RCSid[] = "$Id$";
16   #include  "random.h"
17  
18  
19 < static int  cyl_partit(), flt_partit();
19 > static int
20 > srcskip(                        /* pre-emptive test for out-of-range glow */
21 >        SRCREC  *sp,
22 >        FVECT  orig
23 > )
24 > {
25 >        if (sp->sflags & SSKIP)
26 >                return(1);
27  
28 +        if ((sp->sflags & (SPROX|SDISTANT)) != SPROX)
29 +                return(0);
30  
31 +        return(dist2(orig, sp->sloc) >
32 +                        (sp->sl.prox + sp->srad)*(sp->sl.prox + sp->srad));
33 + }
34 +
35   double
36 < nextssamp(r, si)                /* compute sample for source, rtn. distance */
37 < register RAY  *r;               /* origin is read, direction is set */
38 < register SRCINDEX  *si;         /* source index (modified to current) */
36 > nextssamp(                      /* compute sample for source, rtn. distance */
37 >        RAY  *r,                /* origin is read, direction is set */
38 >        SRCINDEX  *si           /* source index (modified to current) */\
39 > )
40   {
41          int  cent[3], size[3], parr[2];
42          SRCREC  *srcp;
43          FVECT  vpos;
44          double  d;
45 <        register int  i;
45 >        int  i;
46   nextsample:
47          while (++si->sp >= si->np) {    /* get next sample */
48                  if (++si->sn >= nsources)
49                          return(0.0);    /* no more */
50 <                if (source[si->sn].sflags & SSKIP)
50 >                if (srcskip(source+si->sn, r->rorg))
51                          si->np = 0;
52                  else if (srcsizerat <= FTINY)
53                          nopart(si, r);
# Line 131 | Line 145 | nextsample:
145  
146  
147   int
148 < skipparts(ct, sz, pp, pt)               /* skip to requested partition */
149 < int  ct[3], sz[3];              /* center and size of partition (returned) */
150 < register int  pp[2];            /* current index, number to skip (modified) */
151 < unsigned char  *pt;             /* partition array */
148 > skipparts(                      /* skip to requested partition */
149 >        int  ct[3],
150 >        int  sz[3],             /* center and size of partition (returned) */
151 >        int  pp[2],             /* current index, number to skip (modified) */
152 >        unsigned char  *pt      /* partition array */
153 > )
154   {
155 <        register int  p;
155 >        int  p;
156                                          /* check this partition */
157          p = spart(pt, pp[0]);
158          pp[0]++;
# Line 164 | Line 180 | unsigned char  *pt;            /* partition array */
180  
181  
182   void
183 < nopart(si, r)                   /* single source partition */
184 < register SRCINDEX  *si;
185 < RAY  *r;
183 > nopart(                         /* single source partition */
184 >        SRCINDEX  *si,
185 >        RAY  *r
186 > )
187   {
188          clrpart(si->spt);
189          setpart(si->spt, 0, S0);
# Line 174 | Line 191 | RAY  *r;
191   }
192  
193  
177 void
178 cylpart(si, r)                  /* partition a cylinder */
179 SRCINDEX  *si;
180 register RAY  *r;
181 {
182        double  dist2, safedist2, dist2cent, rad2;
183        FVECT  v;
184        register SRCREC  *sp;
185        int  pi;
186                                        /* first check point location */
187        clrpart(si->spt);
188        sp = source + si->sn;
189        rad2 = 1.365 * DOT(sp->ss[SV],sp->ss[SV]);
190        v[0] = r->rorg[0] - sp->sloc[0];
191        v[1] = r->rorg[1] - sp->sloc[1];
192        v[2] = r->rorg[2] - sp->sloc[2];
193        dist2 = DOT(v,sp->ss[SU]);
194        safedist2 = DOT(sp->ss[SU],sp->ss[SU]);
195        dist2 *= dist2 / safedist2;
196        dist2cent = DOT(v,v);
197        dist2 = dist2cent - dist2;
198        if (dist2 <= rad2) {            /* point inside extended cylinder */
199                si->np = 0;
200                return;
201        }
202        safedist2 *= 4.*r->rweight*r->rweight/(srcsizerat*srcsizerat);
203        if (dist2 <= 4.*rad2 ||         /* point too close to subdivide */
204                        dist2cent >= safedist2) {       /* or too far */
205                setpart(si->spt, 0, S0);
206                si->np = 1;
207                return;
208        }
209        pi = 0;
210        si->np = cyl_partit(r->rorg, si->spt, &pi, MAXSPART,
211                        sp->sloc, sp->ss[SU], safedist2);
212 }
213
214
194   static int
195 < cyl_partit(ro, pt, pi, mp, cent, axis, d2)      /* slice a cylinder */
196 < FVECT  ro;
197 < unsigned char  *pt;
198 < register int  *pi;
199 < int  mp;
200 < FVECT  cent, axis;
201 < double  d2;
195 > cyl_partit(                             /* slice a cylinder */
196 >        FVECT  ro,
197 >        unsigned char  *pt,
198 >        int  *pi,
199 >        int  mp,
200 >        FVECT  cent,
201 >        FVECT  axis,
202 >        double  d2
203 > )
204   {
205          FVECT  newct, newax;
206          int  npl, npu;
# Line 252 | Line 233 | double  d2;
233  
234  
235   void
236 < flatpart(si, r)                         /* partition a flat source */
237 < register SRCINDEX  *si;
238 < register RAY  *r;
236 > cylpart(                        /* partition a cylinder */
237 >        SRCINDEX  *si,
238 >        RAY  *r
239 > )
240   {
241 <        register RREAL  *vp;
241 >        double  dist2, safedist2, dist2cent, rad2;
242          FVECT  v;
243 <        double  du2, dv2;
243 >        SRCREC  *sp;
244          int  pi;
245 <
245 >                                        /* first check point location */
246          clrpart(si->spt);
247 <        vp = source[si->sn].sloc;
248 <        v[0] = r->rorg[0] - vp[0];
249 <        v[1] = r->rorg[1] - vp[1];
250 <        v[2] = r->rorg[2] - vp[2];
251 <        vp = source[si->sn].snorm;
252 <        if (DOT(v,vp) <= 0.) {          /* behind source */
247 >        sp = source + si->sn;
248 >        rad2 = 1.365 * DOT(sp->ss[SV],sp->ss[SV]);
249 >        v[0] = r->rorg[0] - sp->sloc[0];
250 >        v[1] = r->rorg[1] - sp->sloc[1];
251 >        v[2] = r->rorg[2] - sp->sloc[2];
252 >        dist2 = DOT(v,sp->ss[SU]);
253 >        safedist2 = DOT(sp->ss[SU],sp->ss[SU]);
254 >        dist2 *= dist2 / safedist2;
255 >        dist2cent = DOT(v,v);
256 >        dist2 = dist2cent - dist2;
257 >        if (dist2 <= rad2) {            /* point inside extended cylinder */
258                  si->np = 0;
259                  return;
260          }
261 <        dv2 = 2.*r->rweight/srcsizerat;
262 <        dv2 *= dv2;
263 <        vp = source[si->sn].ss[SU];
264 <        du2 = dv2 * DOT(vp,vp);
265 <        vp = source[si->sn].ss[SV];
266 <        dv2 *= DOT(vp,vp);
261 >        safedist2 *= 4.*r->rweight*r->rweight/(srcsizerat*srcsizerat);
262 >        if (dist2 <= 4.*rad2 ||         /* point too close to subdivide */
263 >                        dist2cent >= safedist2) {       /* or too far */
264 >                setpart(si->spt, 0, S0);
265 >                si->np = 1;
266 >                return;
267 >        }
268          pi = 0;
269 <        si->np = flt_partit(r->rorg, si->spt, &pi, MAXSPART,
270 <                source[si->sn].sloc,
283 <                source[si->sn].ss[SU], source[si->sn].ss[SV], du2, dv2);
269 >        si->np = cyl_partit(r->rorg, si->spt, &pi, MAXSPART,
270 >                        sp->sloc, sp->ss[SU], safedist2);
271   }
272  
273  
274   static int
275 < flt_partit(ro, pt, pi, mp, cent, u, v, du2, dv2)        /* partition flatty */
276 < FVECT  ro;
277 < unsigned char  *pt;
278 < register int  *pi;
279 < int  mp;
280 < FVECT  cent, u, v;
281 < double  du2, dv2;
275 > flt_partit(                             /* partition flatty */
276 >        FVECT  ro,
277 >        unsigned char  *pt,
278 >        int  *pi,
279 >        int  mp,
280 >        FVECT  cent,
281 >        FVECT  u,
282 >        FVECT  v,
283 >        double  du2,
284 >        double  dv2
285 > )
286   {
287          double  d2;
288          FVECT  newct, newax;
# Line 335 | Line 326 | double  du2, dv2;
326   }
327  
328  
329 + void
330 + flatpart(                               /* partition a flat source */
331 +        SRCINDEX  *si,
332 +        RAY  *r
333 + )
334 + {
335 +        RREAL  *vp;
336 +        FVECT  v;
337 +        double  du2, dv2;
338 +        int  pi;
339 +
340 +        clrpart(si->spt);
341 +        vp = source[si->sn].sloc;
342 +        v[0] = r->rorg[0] - vp[0];
343 +        v[1] = r->rorg[1] - vp[1];
344 +        v[2] = r->rorg[2] - vp[2];
345 +        vp = source[si->sn].snorm;
346 +        if (DOT(v,vp) <= 0.) {          /* behind source */
347 +                si->np = 0;
348 +                return;
349 +        }
350 +        dv2 = 2.*r->rweight/srcsizerat;
351 +        dv2 *= dv2;
352 +        vp = source[si->sn].ss[SU];
353 +        du2 = dv2 * DOT(vp,vp);
354 +        vp = source[si->sn].ss[SV];
355 +        dv2 *= DOT(vp,vp);
356 +        pi = 0;
357 +        si->np = flt_partit(r->rorg, si->spt, &pi, MAXSPART,
358 +                source[si->sn].sloc,
359 +                source[si->sn].ss[SU], source[si->sn].ss[SV], du2, dv2);
360 + }
361 +
362 +
363   double
364 < scylform(sn, dir)               /* compute cosine for cylinder's projection */
365 < int  sn;
366 < register FVECT  dir;            /* assume normalized */
364 > scylform(                       /* compute cosine for cylinder's projection */
365 >        int  sn,
366 >        FVECT  dir              /* assume normalized */
367 > )
368   {
369 <        register RREAL  *dv;
369 >        RREAL  *dv;
370          double  d;
371  
372          dv = source[sn].ss[SU];

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines