Back to index

radiance  4R0+20100331
Functions | Variables
srcsamp.c File Reference
#include "copyright.h"
#include "ray.h"
#include "source.h"
#include "random.h"

Go to the source code of this file.

Functions

static int cyl_partit ()
static int flt_partit ()
double nextssamp (RAY *r, SRCINDEX *si)
int skipparts (ct, sz, pp, unsigned char *pt)
void nopart (SRCINDEX *si, RAY *r)
void cylpart (SRCINDEX *si, RAY *r)
static int cyl_partit (FVECT ro, unsigned char *pt, int *pi, int mp, FVECT cent, FVECT axis, double d2)
void flatpart (SRCINDEX *si, RAY *r)
static int flt_partit (FVECT ro, unsigned char *pt, int *pi, int mp, FVECT cent, FVECT u, FVECT v, double du2, double dv2)
double scylform (int sn, FVECT dir)

Variables

static const char RCSid [] = "$Id: srcsamp.c,v 2.16 2009/06/06 02:11:44 greg Exp $"

Function Documentation

static int cyl_partit ( ) [static]

Here is the caller graph for this function:

static int cyl_partit ( FVECT  ro,
unsigned char *  pt,
int *  pi,
int  mp,
FVECT  cent,
FVECT  axis,
double  d2 
) [static]

Definition at line 218 of file srcsamp.c.

{
       FVECT  newct, newax;
       int  npl, npu;

       if (mp < 2 || dist2(ro, cent) >= d2) {    /* hit limit? */
              setpart(pt, *pi, S0);
              (*pi)++;
              return(1);
       }
                                   /* subdivide */
       setpart(pt, *pi, SU);
       (*pi)++;
       newax[0] = .5*axis[0];
       newax[1] = .5*axis[1];
       newax[2] = .5*axis[2];
       d2 *= 0.25;
                                   /* lower half */
       newct[0] = cent[0] - newax[0];
       newct[1] = cent[1] - newax[1];
       newct[2] = cent[2] - newax[2];
       npl = cyl_partit(ro, pt, pi, mp/2, newct, newax, d2);
                                   /* upper half */
       newct[0] = cent[0] + newax[0];
       newct[1] = cent[1] + newax[1];
       newct[2] = cent[2] + newax[2];
       npu = cyl_partit(ro, pt, pi, mp/2, newct, newax, d2);
                                   /* return total */
       return(npl + npu);
}

Here is the call graph for this function:

void cylpart ( SRCINDEX si,
RAY r 
)

Definition at line 180 of file srcsamp.c.

{
       double  dist2, safedist2, dist2cent, rad2;
       FVECT  v;
       register SRCREC  *sp;
       int  pi;
                                   /* first check point location */
       clrpart(si->spt);
       sp = source + si->sn;
       rad2 = 1.365 * DOT(sp->ss[SV],sp->ss[SV]);
       v[0] = r->rorg[0] - sp->sloc[0];
       v[1] = r->rorg[1] - sp->sloc[1];
       v[2] = r->rorg[2] - sp->sloc[2];
       dist2 = DOT(v,sp->ss[SU]);
       safedist2 = DOT(sp->ss[SU],sp->ss[SU]);
       dist2 *= dist2 / safedist2;
       dist2cent = DOT(v,v);
       dist2 = dist2cent - dist2;
       if (dist2 <= rad2) {        /* point inside extended cylinder */
              si->np = 0;
              return;
       }
       safedist2 *= 4.*r->rweight*r->rweight/(srcsizerat*srcsizerat);
       if (dist2 <= 4.*rad2 ||            /* point too close to subdivide */
                     dist2cent >= safedist2) {   /* or too far */
              setpart(si->spt, 0, S0);
              si->np = 1;
              return;
       }
       pi = 0;
       si->np = cyl_partit(r->rorg, si->spt, &pi, MAXSPART,
                     sp->sloc, sp->ss[SU], safedist2);
}

Here is the call graph for this function:

Here is the caller graph for this function:

void flatpart ( SRCINDEX si,
RAY r 
)

Definition at line 257 of file srcsamp.c.

{
       register RREAL  *vp;
       FVECT  v;
       double  du2, dv2;
       int  pi;

       clrpart(si->spt);
       vp = source[si->sn].sloc;
       v[0] = r->rorg[0] - vp[0];
       v[1] = r->rorg[1] - vp[1];
       v[2] = r->rorg[2] - vp[2];
       vp = source[si->sn].snorm;
       if (DOT(v,vp) <= 0.) {             /* behind source */
              si->np = 0;
              return;
       }
       dv2 = 2.*r->rweight/srcsizerat;
       dv2 *= dv2;
       vp = source[si->sn].ss[SU];
       du2 = dv2 * DOT(vp,vp);
       vp = source[si->sn].ss[SV];
       dv2 *= DOT(vp,vp);
       pi = 0;
       si->np = flt_partit(r->rorg, si->spt, &pi, MAXSPART,
              source[si->sn].sloc,
              source[si->sn].ss[SU], source[si->sn].ss[SV], du2, dv2);
}

Here is the call graph for this function:

Here is the caller graph for this function:

static int flt_partit ( )

Here is the caller graph for this function:

static int flt_partit ( FVECT  ro,
unsigned char *  pt,
int *  pi,
int  mp,
FVECT  cent,
FVECT  u,
FVECT  v,
double  du2,
double  dv2 
) [static]

Definition at line 290 of file srcsamp.c.

{
       double  d2;
       FVECT  newct, newax;
       int  npl, npu;

       if (mp < 2 || ((d2 = dist2(ro, cent)) >= du2
                     && d2 >= dv2)) {     /* hit limit? */
              setpart(pt, *pi, S0);
              (*pi)++;
              return(1);
       }
       if (du2 > dv2) {                   /* subdivide in U */
              setpart(pt, *pi, SU);
              (*pi)++;
              newax[0] = .5*u[0];
              newax[1] = .5*u[1];
              newax[2] = .5*u[2];
              u = newax;
              du2 *= 0.25;
       } else {                           /* subdivide in V */
              setpart(pt, *pi, SV);
              (*pi)++;
              newax[0] = .5*v[0];
              newax[1] = .5*v[1];
              newax[2] = .5*v[2];
              v = newax;
              dv2 *= 0.25;
       }
                                   /* lower half */
       newct[0] = cent[0] - newax[0];
       newct[1] = cent[1] - newax[1];
       newct[2] = cent[2] - newax[2];
       npl = flt_partit(ro, pt, pi, mp/2, newct, u, v, du2, dv2);
                                   /* upper half */
       newct[0] = cent[0] + newax[0];
       newct[1] = cent[1] + newax[1];
       newct[2] = cent[2] + newax[2];
       npu = flt_partit(ro, pt, pi, mp/2, newct, u, v, du2, dv2);
                            /* return total */
       return(npl + npu);
}

Here is the call graph for this function:

double nextssamp ( RAY r,
SRCINDEX si 
)

Definition at line 23 of file srcsamp.c.

{
       int  cent[3], size[3], parr[2];
       SRCREC  *srcp;
       FVECT  vpos;
       double  d;
       register int  i;
nextsample:
       while (++si->sp >= si->np) {       /* get next sample */
              if (++si->sn >= nsources)
                     return(0.0);  /* no more */
              if (source[si->sn].sflags & SSKIP)
                     si->np = 0;
              else if (srcsizerat <= FTINY)
                     nopart(si, r);
              else {
                     for (i = si->sn; source[i].sflags & SVIRTUAL;
                                   i = source[i].sa.sv.sn)
                            ;             /* partition source */
                     (*sfun[source[i].so->otype].of->partit)(si, r);
              }
              si->sp = -1;
       }
                                   /* get partition */
       cent[0] = cent[1] = cent[2] = 0;
       size[0] = size[1] = size[2] = MAXSPART;
       parr[0] = 0; parr[1] = si->sp;
       if (!skipparts(cent, size, parr, si->spt))
              error(CONSISTENCY, "bad source partition in nextssamp");
                                   /* compute sample */
       srcp = source + si->sn;
       if (dstrsrc > FTINY) {                    /* jitter sample */
              dimlist[ndims] = si->sn + 8831;
              dimlist[ndims+1] = si->sp + 3109;
              d = urand(ilhash(dimlist,ndims+2)+samplendx);
              if (srcp->sflags & SFLAT) {
                     multisamp(vpos, 2, d);
                     vpos[SW] = 0.5;
              } else
                     multisamp(vpos, 3, d);
              for (i = 0; i < 3; i++)
                     vpos[i] = dstrsrc * (1. - 2.*vpos[i]) *
                                   (double)size[i]/MAXSPART;
       } else
              vpos[0] = vpos[1] = vpos[2] = 0.0;

       for (i = 0; i < 3; i++)
              vpos[i] += (double)cent[i]/MAXSPART;
                                   /* avoid circular aiming failures */
       if ((srcp->sflags & SCIR) && (si->np > 1 || dstrsrc > 0.7)) {
              FVECT  trim;
              if (srcp->sflags & (SFLAT|SDISTANT)) {
                     d = 1.12837917;             /* correct setflatss() */
                     trim[SU] = d*sqrt(1.0 - 0.5*vpos[SV]*vpos[SV]);
                     trim[SV] = d*sqrt(1.0 - 0.5*vpos[SU]*vpos[SU]);
                     trim[SW] = 0.0;
              } else {
                     trim[SW] = trim[SU] = vpos[SU]*vpos[SU];
                     d = vpos[SV]*vpos[SV];
                     if (d > trim[SW]) trim[SW] = d;
                     trim[SU] += d;
                     d = vpos[SW]*vpos[SW];
                     if (d > trim[SW]) trim[SW] = d;
                     trim[SU] += d;
                     if (trim[SU] > FTINY*FTINY) {
                            d = 1.0/0.7236;      /* correct sphsetsrc() */
                            trim[SW] = trim[SV] = trim[SU] =
                                          d*sqrt(trim[SW]/trim[SU]);
                     } else
                            trim[SW] = trim[SV] = trim[SU] = 0.0;
              }
              for (i = 0; i < 3; i++)
                     vpos[i] *= trim[i];
       }
                                   /* compute direction */
       for (i = 0; i < 3; i++)
              r->rdir[i] = srcp->sloc[i] +
                            vpos[SU]*srcp->ss[SU][i] +
                            vpos[SV]*srcp->ss[SV][i] +
                            vpos[SW]*srcp->ss[SW][i];

       if (!(srcp->sflags & SDISTANT))
              for (i = 0; i < 3; i++)
                     r->rdir[i] -= r->rorg[i];
                                   /* compute distance */
       if ((d = normalize(r->rdir)) == 0.0)
              goto nextsample;            /* at source! */

                                   /* compute sample size */
       if (srcp->sflags & SFLAT) {
              si->dom = sflatform(si->sn, r->rdir);
              si->dom *= size[SU]*size[SV]/(MAXSPART*(double)MAXSPART);
       } else if (srcp->sflags & SCYL) {
              si->dom = scylform(si->sn, r->rdir);
              si->dom *= size[SU]/(double)MAXSPART;
       } else {
              si->dom = size[SU]*size[SV]*(double)size[SW] /
                            (MAXSPART*MAXSPART*(double)MAXSPART) ;
       }
       if (srcp->sflags & SDISTANT) {
              si->dom *= srcp->ss2;
              return(FHUGE);
       }
       if (si->dom <= 1e-4)
              goto nextsample;            /* behind source? */
       si->dom *= srcp->ss2/(d*d);
       return(d);           /* sample OK, return distance */
}

Here is the call graph for this function:

Here is the caller graph for this function:

void nopart ( SRCINDEX si,
RAY r 
)

Definition at line 169 of file srcsamp.c.

{
       clrpart(si->spt);
       setpart(si->spt, 0, S0);
       si->np = 1;
}

Here is the caller graph for this function:

double scylform ( int  sn,
FVECT  dir 
)

Definition at line 341 of file srcsamp.c.

{
       register RREAL  *dv;
       double  d;

       dv = source[sn].ss[SU];
       d = DOT(dir, dv);
       d *= d / DOT(dv,dv);
       return(sqrt(1. - d));
}

Here is the caller graph for this function:

int skipparts ( ct  ,
sz  ,
pp  ,
unsigned char *  pt 
)

Definition at line 136 of file srcsamp.c.

{
       register int  p;
                                   /* check this partition */
       p = spart(pt, pp[0]);
       pp[0]++;
       if (p == S0) {                     /* leaf partition */
              if (pp[1]) {
                     pp[1]--;
                     return(0);    /* not there yet */
              } else
                     return(1);    /* we've arrived */
       }
                            /* else check lower */
       sz[p] >>= 1;
       ct[p] -= sz[p];
       if (skipparts(ct, sz, pp, pt))
              return(1);    /* return hit */
                            /* else check upper */
       ct[p] += sz[p] << 1;
       if (skipparts(ct, sz, pp, pt))
              return(1);    /* return hit */
                            /* else return to starting position */
       ct[p] -= sz[p];
       sz[p] <<= 1;
       return(0);           /* return miss */
}

Here is the call graph for this function:


Variable Documentation

const char RCSid[] = "$Id: srcsamp.c,v 2.16 2009/06/06 02:11:44 greg Exp $" [static]

Definition at line 2 of file srcsamp.c.