Back to index

radiance  4R0+20100331
srcobstr.c
Go to the documentation of this file.
00001 #ifndef lint
00002 static const char RCSid[] = "$Id: srcobstr.c,v 2.15 2009/12/12 00:03:42 greg Exp $";
00003 #endif
00004 /*
00005  * Source occlusion caching routines
00006  */
00007 
00008 #include  "ray.h"
00009 
00010 #include  "otypes.h"
00011 
00012 #include  "otspecial.h"
00013 
00014 #include  "rtotypes.h"
00015 
00016 #include  "source.h"
00017 
00018 #define ABS(x)  ((x)>0 ? (x) : -(x))
00019 
00020 
00021 #if  SHADCACHE                     /* preemptive shadow checking */
00022 
00023 
00024 OBJECT *      antimodlist = NULL;  /* set of clipped materials */
00025 
00026 
00027 static int                         /* cast source ray to first blocker */
00028 castshadow(int sn, FVECT rorg, FVECT rdir)
00029 {
00030        RAY     rt;
00031        
00032        VCOPY(rt.rorg, rorg);
00033        VCOPY(rt.rdir, rdir);
00034        rt.rmax = 0;
00035        rayorigin(&rt, PRIMARY, NULL, NULL);
00036                                    /* check for intersection */
00037        while (localhit(&rt, &thescene)) {
00038               RAY    rt1 = rt;     /* pretend we were aimed at source */
00039               rt1.crtype |= rt1.rtype = SHADOW;
00040               rt1.rdir[0] = -rt.rdir[0];
00041               rt1.rdir[1] = -rt.rdir[1];
00042               rt1.rdir[2] = -rt.rdir[2];
00043               rt1.rod = -rt.rod;
00044               VSUB(rt1.rorg, rt.rop, rt.rdir);
00045               rt1.rot = 1.;
00046               rt1.rsrc = sn;
00047                                    /* record blocker */
00048               if (srcblocker(&rt1))
00049                      return(1);
00050                                    /* move past failed blocker */
00051               VSUM(rt.rorg, rt.rop, rt.rdir, FTINY);
00052               rayclear(&rt);              /* & try again... */
00053        }
00054        return(0);                  /* found no blockers */
00055 }
00056 
00057 
00058 void                               /* initialize occlusion cache */
00059 initobscache(int sn)
00060 {
00061        register SRCREC *srcp = &source[sn];
00062        int           cachelen;
00063        FVECT         rorg, rdir;
00064        RREAL         d;
00065        int           i, j, k;
00066        int           ax, ax1, ax2;
00067 
00068        if (srcp->sflags & (SSKIP|SPROX|SSPOT|SVIRTUAL))
00069               return;                     /* don't cache these */
00070        if (srcp->sflags & SDISTANT)
00071               cachelen = 4*SHADCACHE*SHADCACHE;
00072        else if (srcp->sflags & SFLAT)
00073               cachelen = SHADCACHE*SHADCACHE*3 + (SHADCACHE&1)*SHADCACHE*4;
00074        else /* spherical distribution */
00075               cachelen = SHADCACHE*SHADCACHE*6;
00076                                    /* allocate cache */
00077        srcp->obscache = (OBSCACHE *)malloc(sizeof(OBSCACHE) +
00078                                           sizeof(OBJECT)*(cachelen-1));
00079        if (srcp->obscache == NULL)
00080               error(SYSTEM, "out of memory in initobscache()");
00081                                    /* set parameters */
00082        if (srcp->sflags & SDISTANT) {
00083               RREAL   amax = 0;
00084               for (ax1 = 3; ax1--; )
00085                      if (ABS(srcp->sloc[ax1]) > amax) {
00086                             amax = ABS(srcp->sloc[ax1]);
00087                             ax = ax1;
00088                      }
00089               srcp->obscache->p.d.ax = ax;
00090               ax1 = (ax+1)%3;
00091               ax2 = (ax+2)%3;
00092               VCOPY(srcp->obscache->p.d.o, thescene.cuorg);
00093               if (srcp->sloc[ax] > 0)
00094                      srcp->obscache->p.d.o[ax] += thescene.cusize;
00095               if (srcp->sloc[ax1] < 0)
00096                      srcp->obscache->p.d.o[ax1] += thescene.cusize *
00097                                    srcp->sloc[ax1] / amax;
00098               if (srcp->sloc[ax2] < 0)
00099                      srcp->obscache->p.d.o[ax2] += thescene.cusize *
00100                                    srcp->sloc[ax2] / amax;
00101               srcp->obscache->p.d.e1 = 1. / (thescene.cusize*(1. +
00102                             fabs(srcp->sloc[ax1])/amax));
00103               srcp->obscache->p.d.e2 = 1. / (thescene.cusize*(1. +
00104                             fabs(srcp->sloc[ax2])/amax));
00105        } else if (srcp->sflags & SFLAT) {
00106               VCOPY(srcp->obscache->p.f.u, srcp->ss[SU]);
00107               normalize(srcp->obscache->p.f.u);
00108               fcross(srcp->obscache->p.f.v,
00109                             srcp->snorm, srcp->obscache->p.f.u);
00110        }
00111                                    /* clear cache */
00112        for (i = cachelen; i--; )
00113               srcp->obscache->obs[i] = OVOID;
00114                                    /* cast shadow rays */
00115        if (srcp->sflags & SDISTANT) {
00116               for (k = 3; k--; )
00117                      rdir[k] = -srcp->sloc[k];
00118               for (i = 2*SHADCACHE; i--; )
00119                      for (j = 2*SHADCACHE; j--; ) {
00120                             VCOPY(rorg, srcp->obscache->p.d.o);
00121                             rorg[ax1] += (i+.5) /
00122                                    (2*SHADCACHE*srcp->obscache->p.d.e1);
00123                             rorg[ax2] += (j+.5) /
00124                                    (2*SHADCACHE*srcp->obscache->p.d.e2);
00125                             castshadow(sn, rorg, rdir);
00126                      }
00127        } else if (srcp->sflags & SFLAT) {
00128               d = 0.01*srcp->srad;
00129               VSUM(rorg, srcp->sloc, srcp->snorm, d);
00130               for (i = SHADCACHE; i--; )
00131                      for (j = SHADCACHE; j--; ) {
00132                             d = 2./SHADCACHE*(i+.5) - 1.;
00133                             VSUM(rdir, srcp->snorm,
00134                                           srcp->obscache->p.f.u, d);
00135                             d = 2./SHADCACHE*(j+.5) - 1.;
00136                             VSUM(rdir, rdir, srcp->obscache->p.f.v, d);
00137                             normalize(rdir);
00138                             castshadow(sn, rorg, rdir);
00139                      }
00140               for (k = 2; k--; )
00141                   for (i = SHADCACHE; i--; )
00142                      for (j = SHADCACHE>>1; j--; ) {
00143                             d = 2./SHADCACHE*(i+.5) - 1.;
00144                             if (k)
00145                                    VSUM(rdir, srcp->obscache->p.f.u,
00146                                           srcp->obscache->p.f.v, d);
00147                             else
00148                                    VSUM(rdir, srcp->obscache->p.f.v,
00149                                           srcp->obscache->p.f.u, d);
00150                             d = 1. - 2./SHADCACHE*(j+.5);
00151                             VSUM(rdir, rdir, srcp->snorm, d);
00152                             normalize(rdir);
00153                             castshadow(sn, rorg, rdir);
00154                             d = 2.*DOT(rdir, srcp->snorm);
00155                             rdir[0] = d*srcp->snorm[0] - rdir[0];
00156                             rdir[1] = d*srcp->snorm[1] - rdir[1];
00157                             rdir[2] = d*srcp->snorm[2] - rdir[2];
00158                             castshadow(sn, rorg, rdir);
00159                      }
00160        } else /* spherical distribution */
00161            for (k = 6; k--; ) {
00162               ax = k%3;
00163               ax1 = (k+1)%3;
00164               ax2 = (k+2)%3;
00165               for (i = SHADCACHE; i--; )
00166                      for (j = SHADCACHE; j--; ) {
00167                             rdir[0]=rdir[1]=rdir[2] = 0.;
00168                             rdir[ax] = k<3 ? 1. : -1.;
00169                             rdir[ax1] = 2./SHADCACHE*(i+.5) - 1.;
00170                             rdir[ax2] = 2./SHADCACHE*(j+.5) - 1.;
00171                             normalize(rdir);
00172                             d = 1.05*srcp->srad;
00173                             VSUM(rorg, srcp->sloc, rdir, d);
00174                             castshadow(sn, rorg, rdir);
00175                      }
00176            }
00177 }
00178 
00179 
00180 static OBJECT *                    /* return occluder cache entry */
00181 srcobstructp(register RAY *r)
00182 {
00183        static RNUMBER       lastrno = ~0;
00184        static OBJECT   noobs;
00185        static OBJECT *lastobjp;
00186        SRCREC        *srcp;
00187        int           ondx;
00188 
00189        noobs = OVOID;
00190        if (r->rno == lastrno)
00191               return lastobjp;     /* just recall last pointer */
00192        DCHECK(r->rsrc < 0, CONSISTENCY,
00193                      "srcobstructp() called with unaimed ray");
00194        lastrno = r->rno;
00195        lastobjp = &noobs;
00196        srcp = &source[r->rsrc];
00197        if (srcp->sflags & (SSKIP|SPROX|SSPOT|SVIRTUAL))
00198               return(&noobs);             /* don't cache these */
00199        if (srcp->obscache == NULL)     /* initialize cache */
00200               initobscache(r->rsrc);
00201                                    /* compute cache index */
00202        if (srcp->sflags & SDISTANT) {
00203               int     ax, ax1, ax2;
00204               double  t;
00205               ax = srcp->obscache->p.d.ax;
00206               if ((ax1 = ax+1) >= 3) ax1 -= 3;
00207               if ((ax2 = ax+2) >= 3) ax2 -= 3;
00208               t = (srcp->obscache->p.d.o[ax] - r->rorg[ax]) / srcp->sloc[ax];
00209               if (t <= FTINY)
00210                      return(&noobs); /* could happen if ray is outside */
00211               ondx = 2*SHADCACHE*(int)(2*SHADCACHE*srcp->obscache->p.d.e1 *
00212                             (r->rorg[ax1] + t*srcp->sloc[ax1] -
00213                                    srcp->obscache->p.d.o[ax1]));
00214               ondx += (int)(2*SHADCACHE*srcp->obscache->p.d.e2 *
00215                             (r->rorg[ax2] + t*srcp->sloc[ax2] -
00216                                    srcp->obscache->p.d.o[ax2]));
00217               if ((ondx < 0) | (ondx >= 4*SHADCACHE*SHADCACHE))
00218                      return(&noobs); /* could happen if ray is outside */
00219        } else if (srcp->sflags & SFLAT) {
00220               FVECT   sd;
00221               RREAL   sd0m, sd1m;
00222               sd[0] = -DOT(r->rdir, srcp->obscache->p.f.u);
00223               sd[1] = -DOT(r->rdir, srcp->obscache->p.f.v);
00224               sd[2] = -DOT(r->rdir, srcp->snorm);
00225               if (sd[2] < 0)
00226                      return(&noobs); /* shouldn't happen */
00227               sd0m = ABS(sd[0]);
00228               sd1m = ABS(sd[1]);
00229               if (sd[2] >= sd0m && sd[2] >= sd1m) {
00230                      ondx = SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
00231                                    (1. + sd[0]/sd[2]));
00232                      ondx += (int)(SHADCACHE*(.5-FTINY) *
00233                                    (1. + sd[1]/sd[2]));
00234               } else if (sd0m >= sd1m) {
00235                      ondx = SHADCACHE*SHADCACHE;
00236                      if (sd[0] < 0)
00237                             ondx += ((SHADCACHE+1)>>1)*SHADCACHE;
00238                      ondx += SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
00239                                    (1. - sd[2]/sd0m));
00240                      ondx += (int)(SHADCACHE*(.5-FTINY) *
00241                                    (1. + sd[1]/sd0m));
00242               } else /* sd1m > sd0m */ {
00243                      ondx = SHADCACHE*SHADCACHE +
00244                                    ((SHADCACHE+1)>>1)*SHADCACHE*2;
00245                      if (sd[1] < 0)
00246                             ondx += ((SHADCACHE+1)>>1)*SHADCACHE;
00247                      ondx += SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
00248                                    (1. - sd[2]/sd1m));
00249                      ondx += (int)(SHADCACHE*(.5-FTINY) *
00250                                    (1. + sd[0]/sd1m));
00251               }
00252               DCHECK((ondx < 0) | (ondx >= SHADCACHE*SHADCACHE*3 +
00253                             (SHADCACHE&1)*SHADCACHE*4), CONSISTENCY,
00254                             "flat source cache index out of bounds");
00255        } else /* spherical distribution */ {
00256               int     ax, ax1, ax2;
00257               RREAL   amax = 0;
00258               for (ax1 = 3; ax1--; )
00259                      if (ABS(r->rdir[ax1]) > amax) {
00260                             amax = ABS(r->rdir[ax1]);
00261                             ax = ax1;
00262                      }
00263               if ((ax1 = ax+1) >= 3) ax1 -= 3;
00264               if ((ax2 = ax+2) >= 3) ax2 -= 3;
00265               ondx = 2*SHADCACHE*SHADCACHE * ax;
00266               if (r->rdir[ax] < 0)
00267                      ondx += SHADCACHE*SHADCACHE;
00268               ondx += SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
00269                                    (1. + r->rdir[ax1]/amax));
00270               ondx += (int)(SHADCACHE*(.5-FTINY) *
00271                             (1. + r->rdir[ax2]/amax));
00272               DCHECK((ondx < 0) | (ondx >= SHADCACHE*SHADCACHE*6), CONSISTENCY,
00273                             "radial source cache index out of bounds");
00274        }
00275                                    /* return cache pointer */
00276        return(lastobjp = &srcp->obscache->obs[ondx]);
00277 }
00278 
00279 
00280 void                        /* free obstruction cache */
00281 freeobscache(SRCREC *srcp)
00282 {
00283        if (srcp->obscache == NULL)
00284               return;
00285        free((void *)srcp->obscache);
00286        srcp->obscache = NULL;
00287 }
00288 
00289        
00290 int                         /* record a source blocker */
00291 srcblocker(register RAY *r)
00292 {
00293        OBJREC  *m;
00294 
00295        if (r->robj == OVOID || objptr(r->robj) != r->ro ||
00296                      isvolume(r->ro->otype))
00297               return(0);           /* don't record complex blockers */
00298        if (r->rsrc < 0 || source[r->rsrc].so == r->ro)
00299               return(0);           /* just a mistake, that's all */
00300        if (antimodlist != NULL && inset(antimodlist, r->ro->omod))
00301               return(0);           /* could be clipped */
00302        m = findmaterial(r->ro);
00303        if (m == NULL)
00304               return(0);           /* no material?! */
00305        if (!isopaque(m->otype))
00306               return(0);           /* material not a reliable blocker */
00307        *srcobstructp(r) = r->robj;     /* else record obstructor */
00308        return(1);
00309 }
00310 
00311 
00312 int                         /* check ray against cached blocker */
00313 srcblocked(RAY *r)
00314 {
00315        OBJECT  obs = *srcobstructp(r);
00316        OBJREC  *op;
00317 
00318        if (obs == OVOID)
00319               return(0);
00320        op = objptr(obs);           /* check blocker intersection */
00321        if (!(*ofun[op->otype].funp)(op, r))
00322               return(0);
00323        if (source[r->rsrc].sflags & SDISTANT)
00324               return(1);
00325        op = source[r->rsrc].so;    /* check source intersection */
00326        if (!(*ofun[op->otype].funp)(op, r))
00327               return(1);
00328        rayclear(r);
00329        return(0);                  /* source in front */
00330 }
00331 
00332 
00333 void                        /* record potentially clipped materials */
00334 markclip(OBJREC *m)
00335 {
00336        OBJECT  *set2add, *oldset;
00337 
00338        m_clip(m, NULL);            /* initialize modifier list */
00339        if ((set2add = (OBJECT *)m->os) == NULL || !set2add[0])
00340               return;
00341 
00342        if (antimodlist == NULL) {  /* start of list */
00343               antimodlist = setsave(set2add);
00344               return;
00345        }
00346                                    /* else add to previous list */
00347        oldset = antimodlist;
00348        antimodlist = (OBJECT *)malloc((oldset[0]+set2add[0]+1)*sizeof(OBJECT));
00349        if (antimodlist == NULL)
00350               error(SYSTEM, "out of memory in markclip");
00351        setunion(antimodlist, oldset, set2add);
00352        free((void *)oldset);
00353 }
00354 
00355 
00356 #else  /* SHADCACHE */
00357 
00358 
00359 void                        /* no-op also avoids linker warning */
00360 markclip(OBJREC *m)
00361 {
00362        (void)m;
00363 }
00364 
00365 
00366 #endif  /* SHADCACHE */