ViewVC Help

View File | Revision Log | Show Annotations | Download File | Root Listing

root/Development/ray/src/cv/bsdfrep.c

Comparing ray/src/cv/bsdfrep.c (file contents):
Revision 2.17 by greg, Thu Oct 24 16:11:38 2013 UTC vs.
Revision 2.25 by greg, Thu Aug 21 10:33:48 2014 UTC

#	Line 14 | Line 14 | static const char RCSid[] = "$Id$";
14		#include "rtio.h"
15		#include "resolu.h"
16		#include "bsdfrep.h"
17	+	/* name and manufacturer if known */
18	+	char                    bsdf_name[256];
19	+	char                    bsdf_manuf[256];
20		/* active grid resolution */
21		int                     grid_res = GRIDRES;
22
#	Line 31 | Line 34 | unsigned long          bsdf_hist[HISTLEN];
34
35		/* BSDF value for boundary regions */
36		double                  bsdf_min = 0;
37	+	float                   bsdf_spec_peak = 0;
38	+	int                     bsdf_spec_crad = 0;
39
40		/* processed incident DSF measurements */
41		RBFNODE                 *dsf_list = NULL;
#	Line 195 | Line 200 | rotate_rbf(RBFNODE *rbf, const FVECT invec)
200		int                     pos[2];
201		int                     n;
202
203	<	for (n = ((-.01 > phi) | (phi > .01))*rbf->nrbf; n-- > 0; ) {
203	>	for (n = (cos(phi) < 1.-FTINY)*rbf->nrbf; n-- > 0; ) {
204		ovec_from_pos(outvec, rbf->rbfa[n].gx, rbf->rbfa[n].gy);
205		spinvector(outvec, outvec, vnorm, phi);
206		pos_from_vec(pos, outvec);
#	Line 254 | Line 259 | rbf_volume(const RBFVAL *rbfp)
259		return(integ);
260		}
261
262	<	/* Evaluate RBF for DSF at the given normalized outgoing direction */
262	>	/* Evaluate BSDF at the given normalized outgoing direction */
263		double
264		eval_rbfrep(const RBFNODE *rp, const FVECT outvec)
265		{
266		const double    rfact2 = (38./M_PI/M_PI)*(grid_res*grid_res);
262	–	double          minval = bsdf_min*output_orient*outvec[2];
267		int             pos[2];
268		double          res = 0;
269		const RBFVAL    *rbfp;
#	Line 271 | Line 275 | eval_rbfrep(const RBFNODE *rp, const FVECT outvec)
275		return(.0);
276		/* use minimum if no information avail. */
277		if (rp == NULL)
278	<	return(minval);
278	>	return(bsdf_min);
279		/* optimization for fast lobe culling */
280		pos_from_vec(pos, outvec);
281		/* sum radial basis function */
#	Line 286 | Line 290 | eval_rbfrep(const RBFNODE *rp, const FVECT outvec)
290		ovec_from_pos(odir, rbfp->gx, rbfp->gy);
291		res += rbfp->peak * exp((DOT(odir,outvec) - 1.) / rad2);
292		}
293	<	if (res < minval)       /* never return less than minval */
294	<	return(minval);
293	>	res /= output_orient*outvec[2];
294	>	if (res < bsdf_min)     /* never return less than bsdf_min */
295	>	return(bsdf_min);
296		return(res);
297		}
298
#	Line 301 | Line 306 | insert_dsf(RBFNODE *newrbf)
306		for (rbf = dsf_list; rbf != NULL; rbf = rbf->next)
307		if (DOT(rbf->invec, newrbf->invec) >= 1.-FTINY) {
308		fprintf(stderr,
309	<	"%s: Duplicate incident measurement (ignored)\n",
310	<	progname);
309	>	"%s: Duplicate incident measurement ignored at (%.1f,%.1f)\n",
310	>	progname, get_theta180(newrbf->invec),
311	>	get_phi360(newrbf->invec));
312		free(newrbf);
313		return(-1);
314		}
#	Line 390 | Line 396 | get_triangles(RBFNODE *rbfv[2], const MIGRATION *mig)
396		return((rbfv[0] != NULL) + (rbfv[1] != NULL));
397		}
398
399	+	/* Return single-lobe specular RBF for the given incident direction */
400	+	RBFNODE *
401	+	def_rbf_spec(const FVECT invec)
402	+	{
403	+	RBFNODE         *rbf;
404	+	FVECT           ovec;
405	+	int             pos[2];
406	+
407	+	if (input_orient > 0 ^ invec[2] > 0)    /* wrong side? */
408	+	return(NULL);
409	+	if ((bsdf_spec_peak <= bsdf_min) | (bsdf_spec_crad <= 0))
410	+	return(NULL);                   /* nothing set */
411	+	rbf = (RBFNODE *)malloc(sizeof(RBFNODE));
412	+	if (rbf == NULL)
413	+	return(NULL);
414	+	ovec[0] = -invec[0];
415	+	ovec[1] = -invec[1];
416	+	ovec[2] = invec[2]*(2*(input_orient==output_orient) - 1);
417	+	pos_from_vec(pos, ovec);
418	+	rbf->ord = 0;
419	+	rbf->next = NULL;
420	+	rbf->ejl = NULL;
421	+	VCOPY(rbf->invec, invec);
422	+	rbf->nrbf = 1;
423	+	rbf->rbfa[0].peak = bsdf_spec_peak * output_orient*ovec[2];
424	+	rbf->rbfa[0].crad = bsdf_spec_crad;
425	+	rbf->rbfa[0].gx = pos[0];
426	+	rbf->rbfa[0].gy = pos[1];
427	+	rbf->vtotal = rbf_volume(rbf->rbfa);
428	+	return(rbf);
429	+	}
430	+
431	+	/* Advect and allocate new RBF along edge (internal call) */
432	+	RBFNODE *
433	+	e_advect_rbf(const MIGRATION *mig, const FVECT invec, int lobe_lim)
434	+	{
435	+	double          cthresh = FTINY;
436	+	RBFNODE         *rbf;
437	+	int             n, i, j;
438	+	double          t, full_dist;
439	+	/* get relative position */
440	+	t = Acos(DOT(invec, mig->rbfv[0]->invec));
441	+	if (t < M_PI/grid_res) {                /* near first DSF */
442	+	n = sizeof(RBFNODE) + sizeof(RBFVAL)*(mig->rbfv[0]->nrbf-1);
443	+	rbf = (RBFNODE *)malloc(n);
444	+	if (rbf == NULL)
445	+	goto memerr;
446	+	memcpy(rbf, mig->rbfv[0], n);   /* just duplicate */
447	+	rbf->next = NULL; rbf->ejl = NULL;
448	+	return(rbf);
449	+	}
450	+	full_dist = acos(DOT(mig->rbfv[0]->invec, mig->rbfv[1]->invec));
451	+	if (t > full_dist-M_PI/grid_res) {      /* near second DSF */
452	+	n = sizeof(RBFNODE) + sizeof(RBFVAL)*(mig->rbfv[1]->nrbf-1);
453	+	rbf = (RBFNODE *)malloc(n);
454	+	if (rbf == NULL)
455	+	goto memerr;
456	+	memcpy(rbf, mig->rbfv[1], n);   /* just duplicate */
457	+	rbf->next = NULL; rbf->ejl = NULL;
458	+	return(rbf);
459	+	}
460	+	t /= full_dist;
461	+	tryagain:
462	+	n = 0;                                  /* count migrating particles */
463	+	for (i = 0; i < mtx_nrows(mig); i++)
464	+	for (j = 0; j < mtx_ncols(mig); j++)
465	+	n += (mtx_coef(mig,i,j) > cthresh);
466	+	/* are we over our limit? */
467	+	if ((lobe_lim > 0) & (n > lobe_lim)) {
468	+	cthresh = cthresh*2. + 10.*FTINY;
469	+	goto tryagain;
470	+	}
471	+	#ifdef DEBUG
472	+	fprintf(stderr, "Input RBFs have %d, %d nodes -> output has %d\n",
473	+	mig->rbfv[0]->nrbf, mig->rbfv[1]->nrbf, n);
474	+	#endif
475	+	rbf = (RBFNODE *)malloc(sizeof(RBFNODE) + sizeof(RBFVAL)*(n-1));
476	+	if (rbf == NULL)
477	+	goto memerr;
478	+	rbf->next = NULL; rbf->ejl = NULL;
479	+	VCOPY(rbf->invec, invec);
480	+	rbf->nrbf = n;
481	+	rbf->vtotal = 1.-t + t*mig->rbfv[1]->vtotal/mig->rbfv[0]->vtotal;
482	+	n = 0;                                  /* advect RBF lobes */
483	+	for (i = 0; i < mtx_nrows(mig); i++) {
484	+	const RBFVAL        *rbf0i = &mig->rbfv[0]->rbfa[i];
485	+	const float         peak0 = rbf0i->peak;
486	+	const double        rad0 = R2ANG(rbf0i->crad);
487	+	FVECT               v0;
488	+	float               mv;
489	+	ovec_from_pos(v0, rbf0i->gx, rbf0i->gy);
490	+	for (j = 0; j < mtx_ncols(mig); j++)
491	+	if ((mv = mtx_coef(mig,i,j)) > cthresh) {
492	+	const RBFVAL    *rbf1j = &mig->rbfv[1]->rbfa[j];
493	+	double          rad2;
494	+	FVECT           v;
495	+	int             pos[2];
496	+	rad2 = R2ANG(rbf1j->crad);
497	+	rad2 = rad0*rad0*(1.-t) + rad2*rad2*t;
498	+	rbf->rbfa[n].peak = peak0 * mv * rbf->vtotal *
499	+	rad0*rad0/rad2;
500	+	rbf->rbfa[n].crad = ANG2R(sqrt(rad2));
501	+	ovec_from_pos(v, rbf1j->gx, rbf1j->gy);
502	+	geodesic(v, v0, v, t, GEOD_REL);
503	+	pos_from_vec(pos, v);
504	+	rbf->rbfa[n].gx = pos[0];
505	+	rbf->rbfa[n].gy = pos[1];
506	+	++n;
507	+	}
508	+	}
509	+	rbf->vtotal *= mig->rbfv[0]->vtotal;    /* turn ratio into actual */
510	+	return(rbf);
511	+	memerr:
512	+	fprintf(stderr, "%s: Out of memory in e_advect_rbf()\n", progname);
513	+	exit(1);
514	+	return(NULL);   /* pro forma return */
515	+	}
516	+
517		/* Clear our BSDF representation and free memory */
518		void
519		clear_bsdf_rep(void)
#	Line 404 | Line 528 | clear_bsdf_rep(void)
528		dsf_list = rbf->next;
529		free(rbf);
530		}
531	+	bsdf_name[0] = '\0';
532	+	bsdf_manuf[0] = '\0';
533		inp_coverage = 0;
534		single_plane_incident = -1;
535		input_orient = output_orient = 0;
536		grid_res = GRIDRES;
537	+	bsdf_min = 0;
538	+	bsdf_spec_peak = 0;
539	+	bsdf_spec_crad = 0;
540		}
541
542		/* Write our BSDF mesh interpolant out to the given binary stream */
#	Line 418 | Line 547 | save_bsdf_rep(FILE *ofp)
547		MIGRATION       *mig;
548		int             i, n;
549		/* finish header */
550	+	if (bsdf_name[0])
551	+	fprintf(ofp, "NAME=%s\n", bsdf_name);
552	+	if (bsdf_manuf[0])
553	+	fprintf(ofp, "MANUFACT=%s\n", bsdf_manuf);
554		fprintf(ofp, "SYMMETRY=%d\n", !single_plane_incident * inp_coverage);
555		fprintf(ofp, "IO_SIDES= %d %d\n", input_orient, output_orient);
556		fprintf(ofp, "GRIDRES=%d\n", grid_res);
557		fprintf(ofp, "BSDFMIN=%g\n", bsdf_min);
558	+	if ((bsdf_spec_peak > bsdf_min) & (bsdf_spec_crad > 0))
559	+	fprintf(ofp, "BSDFSPEC= %f %f\n", bsdf_spec_peak,
560	+	R2ANG(bsdf_spec_crad));
561		fputformat(BSDFREP_FMT, ofp);
562		fputc('\n', ofp);
563		/* write each DSF */
#	Line 474 | Line 610 | headline(char *s, void *p)
610		{
611		char    fmt[32];
612
613	+	if (!strncmp(s, "NAME=", 5)) {
614	+	strcpy(bsdf_name, s+5);
615	+	bsdf_name[strlen(bsdf_name)-1] = '\0';
616	+	}
617	+	if (!strncmp(s, "MANUFACT=", 9)) {
618	+	strcpy(bsdf_manuf, s+9);
619	+	bsdf_manuf[strlen(bsdf_manuf)-1] = '\0';
620	+	}
621		if (!strncmp(s, "SYMMETRY=", 9)) {
622		inp_coverage = atoi(s+9);
623		single_plane_incident = !inp_coverage;
#	Line 491 | Line 635 | headline(char *s, void *p)
635		sscanf(s+8, "%lf", &bsdf_min);
636		return(0);
637		}
638	+	if (!strncmp(s, "BSDFSPEC=", 9)) {
639	+	float   bsdf_spec_rad = 0;
640	+	sscanf(s+9, "%f %f", &bsdf_spec_peak, &bsdf_spec_rad);
641	+	bsdf_spec_crad = ANG2R(bsdf_spec_rad);
642	+	return(0);
643	+	}
644		if (formatval(fmt, s) && strcmp(fmt, BSDFREP_FMT))
645		return(-1);
646		return(0);
#	Line 507 | Line 657 | load_bsdf_rep(FILE *ifp)
657		clear_bsdf_rep();
658		if (ifp == NULL)
659		return(0);
660	<	if (getheader(ifp, headline, NULL) < 0 || single_plane_incident < 0 |
661	<	!input_orient | !output_orient) {
660	>	if (getheader(ifp, headline, NULL) < 0 || (single_plane_incident < 0) |
661	>	!input_orient | !output_orient |
662	>	(grid_res < 16) | (grid_res > 256)) {
663		fprintf(stderr, "%s: missing/bad format for BSDF interpolant\n",
664		progname);
665		return(0);
666		}
667	<	rbfh.next = NULL;               /* read each DSF */
517	<	rbfh.ejl = NULL;
667	>	memset(&rbfh, 0, sizeof(rbfh)); /* read each DSF */
668		while ((rbfh.ord = getint(4, ifp)) >= 0) {
669		RBFNODE         *newrbf;
670
#	Line 531 | Line 681 | load_bsdf_rep(FILE *ifp)
681		sizeof(RBFVAL)*(rbfh.nrbf-1));
682		if (newrbf == NULL)
683		goto memerr;
684	<	memcpy(newrbf, &rbfh, sizeof(RBFNODE)-sizeof(RBFVAL));
684	>	*newrbf = rbfh;
685		for (i = 0; i < rbfh.nrbf; i++) {
686		newrbf->rbfa[i].peak = getflt(ifp);
687		newrbf->rbfa[i].crad = getint(2, ifp) & 0xffff;

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines (old)
>	Changed lines (new)