--- ray/src/util/rfluxmtx.c 2020/03/02 22:00:05 2.51 +++ ray/src/util/rfluxmtx.c 2024/03/19 10:48:05 2.57 @@ -1,5 +1,5 @@ #ifndef lint -static const char RCSid[] = "$Id: rfluxmtx.c,v 2.51 2020/03/02 22:00:05 greg Exp $"; +static const char RCSid[] = "$Id: rfluxmtx.c,v 2.57 2024/03/19 10:48:05 greg Exp $"; #endif /* * Calculate flux transfer matrix or matrices using rcontrib @@ -41,6 +41,7 @@ char *shirchiufn = "disk2square.cal"; char *kfullfn = "klems_full.cal"; char *khalffn = "klems_half.cal"; char *kquarterfn = "klems_quarter.cal"; +char *ciefn = "cieskyscan.cal"; /* string indicating parameters */ const char PARAMSTART[] = "@rfluxmtx"; @@ -499,6 +500,14 @@ finish_receiver(void) calfn = kquarterfn; kquarterfn = NULL; binf = "kqbin"; nbins = "Nkqbins"; + } else if (!strcasecmp(curparams.hemis, "cie")) { + calfn = ciefn; ciefn = NULL; + sprintf(sbuf, "rNx=%g,rNy=%g,rNz=%g,Ux=%g,Uy=%g,Uz=%g,RHS=%c1", + curparams.nrm[0], curparams.nrm[1], curparams.nrm[2], + curparams.vup[0], curparams.vup[1], curparams.vup[2], + curparams.sign); + binv = "cbin"; + nbins = "Ncbins"; } else { fprintf(stderr, "%s: unrecognized hemisphere sampling: h=%s\n", progname, curparams.hemis); @@ -720,9 +729,9 @@ sample_origin(PARAMS *p, FVECT orig, const FVECT rdir, /* compute projected areas */ for (i = 0, sp = p->slist; sp != NULL; i++, sp = sp->next) { projsa[i] = -DOT(sp->snrm, rdir) * sp->area; - tarea += projsa[i] *= (double)(projsa[i] > FTINY); + tarea += projsa[i] *= (double)(projsa[i] > 0); } - if (tarea < 0) { /* wrong side of sender? */ + if (tarea < FTINY*FTINY) { /* wrong side of sender? */ fputs(progname, stderr); fputs(": internal - sample behind all sender elements!\n", stderr); @@ -740,8 +749,7 @@ sample_uniform(PARAMS *p, int b, FILE *fp) { int n = sampcnt; double samp3[3]; - double duvw[3]; - FVECT orig_dir[2]; + FVECT duvw, orig_dir[2]; int i; if (fp == NULL) /* just requesting number of bins? */ @@ -749,7 +757,7 @@ sample_uniform(PARAMS *p, int b, FILE *fp) while (n--) { /* stratified hemisphere sampling */ SDmultiSamp(samp3, 3, (n+frandom())/sampcnt); - SDsquare2disk(duvw, samp3[1], samp3[2]); + square2disk(duvw, samp3[1], samp3[2]); duvw[2] = -sqrt(1. - duvw[0]*duvw[0] - duvw[1]*duvw[1]); for (i = 3; i--; ) orig_dir[1][i] = duvw[0]*p->udir[i] + @@ -769,8 +777,7 @@ sample_shirchiu(PARAMS *p, int b, FILE *fp) { int n = sampcnt; double samp3[3]; - double duvw[3]; - FVECT orig_dir[2]; + FVECT duvw, orig_dir[2]; int i; if (fp == NULL) /* just requesting number of bins? */ @@ -778,7 +785,7 @@ sample_shirchiu(PARAMS *p, int b, FILE *fp) while (n--) { /* stratified sampling */ SDmultiSamp(samp3, 3, (n+frandom())/sampcnt); - SDsquare2disk(duvw, (b/p->hsiz + samp3[1])/curparams.hsiz, + square2disk(duvw, (b/p->hsiz + samp3[1])/curparams.hsiz, (b%p->hsiz + samp3[2])/curparams.hsiz); duvw[2] = sqrt(1. - duvw[0]*duvw[0] - duvw[1]*duvw[1]); for (i = 3; i--; ) @@ -997,7 +1004,7 @@ finish_polygon(SURF *p) VCOPY(e1, e2); } p->area = normalize(p->snrm)*0.5; - return(p->area > FTINY); + return(p->area > FTINY*FTINY); } /* Add a surface to our current parameters */ @@ -1060,7 +1067,7 @@ add_surface(int st, const char *oname, FILE *fp) snew->area *= PI*snew->area; break; } - if ((snew->area <= FTINY) & (verbose >= 0)) { + if ((snew->area <= FTINY*FTINY) & (verbose >= 0)) { fprintf(stderr, "%s: warning - zero area for surface '%s'\n", progname, oname); free(snew); @@ -1285,12 +1292,22 @@ main(int argc, char *argv[]) yrs = argv[++a]; na = 0; continue; - case 'c': /* number of samples */ - sampcnt = atoi(argv[++a]); - if (sampcnt <= 0) - goto userr; - na = 0; /* we re-add this later */ - continue; + case 'c': /* spectral sampling or count */ + switch (argv[a][2]) { + case 's': + na = 2; + break; + case 'w': + na = 3; + break; + case '\0': + sampcnt = atoi(argv[++a]); + if (sampcnt <= 0) + goto userr; + na = 0; /* we re-add this later */ + continue; + } + break; case 'I': /* only for pass-through mode */ case 'i': iropt = argv[a]; @@ -1308,6 +1325,7 @@ main(int argc, char *argv[]) case 'n': /* options with 1 argument */ case 's': case 'o': + case 't': na = 2; break; case 'b': /* special case */