[Bf-blender-cvs] SVN commit: /data/svn/bf-blender [27242] branches/render25/source/blender: Render Branch: move hemi-polygon form factor code form approximate AO

Brecht Van Lommel brecht at blender.org
Wed Mar 3 11:35:50 CET 2010


Revision: 27242
          http://projects.blender.org/plugins/scmsvn/viewcvs.php?view=rev&root=bf-blender&revision=27242
Author:   blendix
Date:     2010-03-03 11:35:50 +0100 (Wed, 03 Mar 2010)

Log Message:
-----------
Render Branch: move hemi-polygon form factor code form approximate AO
into math library.

Modified Paths:
--------------
    branches/render25/source/blender/blenlib/BLI_math_geom.h
    branches/render25/source/blender/blenlib/intern/math_geom.c
    branches/render25/source/blender/render/intern/include/lamp.h
    branches/render25/source/blender/render/intern/source/diskocclusion.c

Modified: branches/render25/source/blender/blenlib/BLI_math_geom.h
===================================================================
--- branches/render25/source/blender/blenlib/BLI_math_geom.h	2010-03-03 10:21:04 UTC (rev 27241)
+++ branches/render25/source/blender/blenlib/BLI_math_geom.h	2010-03-03 10:35:50 UTC (rev 27242)
@@ -193,6 +193,11 @@
 MINLINE float eval_shv3(float r[9], float v[3]);
 MINLINE void vec_fac_to_sh(float r[9], float v[3], float fac);
 
+/********************************* Form Factor *******************************/
+
+float form_factor_hemi_poly(float p[3], float n[3],
+	float v1[3], float v2[3], float v3[3], float v4[3]);
+
 #ifdef __cplusplus
 }
 #endif

Modified: branches/render25/source/blender/blenlib/intern/math_geom.c
===================================================================
--- branches/render25/source/blender/blenlib/intern/math_geom.c	2010-03-03 10:21:04 UTC (rev 27241)
+++ branches/render25/source/blender/blenlib/intern/math_geom.c	2010-03-03 10:35:50 UTC (rev 27242)
@@ -2051,3 +2051,424 @@
 	}
 }
 
+/******************************* Form Factor *********************************/
+
+static void vec_add_dir(float r[3], float v1[3], float v2[3], float fac)
+{
+	r[0]= v1[0] + fac*(v2[0] - v1[0]);
+	r[1]= v1[1] + fac*(v2[1] - v1[1]);
+	r[2]= v1[2] + fac*(v2[2] - v1[2]);
+}
+
+static int ff_visible_quad(float p[3], float n[3], float v0[3], float v1[3], float v2[3], float q0[3], float q1[3], float q2[3], float q3[3])
+{
+	static const float epsilon = 1e-6f;
+	float c, sd[3];
+	
+	c= dot_v3v3(n, p);
+
+	/* signed distances from the vertices to the plane. */
+	sd[0]= dot_v3v3(n, v0) - c;
+	sd[1]= dot_v3v3(n, v1) - c;
+	sd[2]= dot_v3v3(n, v2) - c;
+
+	if(fabsf(sd[0]) < epsilon) sd[0] = 0.0f;
+	if(fabsf(sd[1]) < epsilon) sd[1] = 0.0f;
+	if(fabsf(sd[2]) < epsilon) sd[2] = 0.0f;
+
+	if(sd[0] > 0) {
+		if(sd[1] > 0) {
+			if(sd[2] > 0) {
+				// +++
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// ++-
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				vec_add_dir(q2, v1, v2, (sd[1]/(sd[1]-sd[2])));
+				vec_add_dir(q3, v0, v2, (sd[0]/(sd[0]-sd[2])));
+			}
+			else {
+				// ++0
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+		}
+		else if(sd[1] < 0) {
+			if(sd[2] > 0) {
+				// +-+
+				copy_v3_v3(q0, v0);
+				vec_add_dir(q1, v0, v1, (sd[0]/(sd[0]-sd[1])));
+				vec_add_dir(q2, v1, v2, (sd[1]/(sd[1]-sd[2])));
+				copy_v3_v3(q3, v2);
+			}
+			else if(sd[2] < 0) {
+				// +--
+				copy_v3_v3(q0, v0);
+				vec_add_dir(q1, v0, v1, (sd[0]/(sd[0]-sd[1])));
+				vec_add_dir(q2, v0, v2, (sd[0]/(sd[0]-sd[2])));
+				copy_v3_v3(q3, q2);
+			}
+			else {
+				// +-0
+				copy_v3_v3(q0, v0);
+				vec_add_dir(q1, v0, v1, (sd[0]/(sd[0]-sd[1])));
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+		}
+		else {
+			if(sd[2] > 0) {
+				// +0+
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// +0-
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				vec_add_dir(q2, v0, v2, (sd[0]/(sd[0]-sd[2])));
+				copy_v3_v3(q3, q2);
+			}
+			else {
+				// +00
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+		}
+	}
+	else if(sd[0] < 0) {
+		if(sd[1] > 0) {
+			if(sd[2] > 0) {
+				// -++
+				vec_add_dir(q0, v0, v1, (sd[0]/(sd[0]-sd[1])));
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				vec_add_dir(q3, v0, v2, (sd[0]/(sd[0]-sd[2])));
+			}
+			else if(sd[2] < 0) {
+				// -+-
+				vec_add_dir(q0, v0, v1, (sd[0]/(sd[0]-sd[1])));
+				copy_v3_v3(q1, v1);
+				vec_add_dir(q2, v1, v2, (sd[1]/(sd[1]-sd[2])));
+				copy_v3_v3(q3, q2);
+			}
+			else {
+				// -+0
+				vec_add_dir(q0, v0, v1, (sd[0]/(sd[0]-sd[1])));
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+		}
+		else if(sd[1] < 0) {
+			if(sd[2] > 0) {
+				// --+
+				vec_add_dir(q0, v0, v2, (sd[0]/(sd[0]-sd[2])));
+				vec_add_dir(q1, v1, v2, (sd[1]/(sd[1]-sd[2])));
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// ---
+				return 0;
+			}
+			else {
+				// --0
+				return 0;
+			}
+		}
+		else {
+			if(sd[2] > 0) {
+				// -0+
+				vec_add_dir(q0, v0, v2, (sd[0]/(sd[0]-sd[2])));
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// -0-
+				return 0;
+			}
+			else {
+				// -00
+				return 0;
+			}
+		}
+	}
+	else {
+		if(sd[1] > 0) {
+			if(sd[2] > 0) {
+				// 0++
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// 0+-
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				vec_add_dir(q2, v1, v2, (sd[1]/(sd[1]-sd[2])));
+				copy_v3_v3(q3, q2);
+			}
+			else {
+				// 0+0
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+		}
+		else if(sd[1] < 0) {
+			if(sd[2] > 0) {
+				// 0-+
+				copy_v3_v3(q0, v0);
+				vec_add_dir(q1, v1, v2, (sd[1]/(sd[1]-sd[2])));
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// 0--
+				return 0;
+			}
+			else {
+				// 0-0
+				return 0;
+			}
+		}
+		else {
+			if(sd[2] > 0) {
+				// 00+
+				copy_v3_v3(q0, v0);
+				copy_v3_v3(q1, v1);
+				copy_v3_v3(q2, v2);
+				copy_v3_v3(q3, q2);
+			}
+			else if(sd[2] < 0) {
+				// 00-
+				return 0;
+			}
+			else {
+				// 000
+				return 0;
+			}
+		}
+	}
+
+	return 1;
+}
+
+/* altivec optimization, this works, but is unused */
+
+#if 0
+#include <Accelerate/Accelerate.h>
+
+typedef union {
+	vFloat v;
+	float f[4];
+} vFloatResult;
+
+static vFloat vec_splat_float(float val)
+{
+	return (vFloat){val, val, val, val};
+}
+
+static float ff_quad_form_factor(float *p, float *n, float *q0, float *q1, float *q2, float *q3)
+{
+	vFloat vcos, rlen, vrx, vry, vrz, vsrx, vsry, vsrz, gx, gy, gz, vangle;
+	vUInt8 rotate = (vUInt8){4,5,6,7,8,9,10,11,12,13,14,15,0,1,2,3};
+	vFloatResult vresult;
+	float result;
+
+	/* compute r* */
+	vrx = (vFloat){q0[0], q1[0], q2[0], q3[0]} - vec_splat_float(p[0]);
+	vry = (vFloat){q0[1], q1[1], q2[1], q3[1]} - vec_splat_float(p[1]);
+	vrz = (vFloat){q0[2], q1[2], q2[2], q3[2]} - vec_splat_float(p[2]);
+
+	/* normalize r* */
+	rlen = vec_rsqrte(vrx*vrx + vry*vry + vrz*vrz + vec_splat_float(1e-16f));
+	vrx = vrx*rlen;
+	vry = vry*rlen;
+	vrz = vrz*rlen;
+
+	/* rotate r* for cross and dot */
+	vsrx= vec_perm(vrx, vrx, rotate);
+	vsry= vec_perm(vry, vry, rotate);
+	vsrz= vec_perm(vrz, vrz, rotate);
+
+	/* cross product */
+	gx = vsry*vrz - vsrz*vry;
+	gy = vsrz*vrx - vsrx*vrz;
+	gz = vsrx*vry - vsry*vrx;
+
+	/* normalize */
+	rlen = vec_rsqrte(gx*gx + gy*gy + gz*gz + vec_splat_float(1e-16f));
+	gx = gx*rlen;
+	gy = gy*rlen;
+	gz = gz*rlen;
+
+	/* angle */
+	vcos = vrx*vsrx + vry*vsry + vrz*vsrz;
+	vcos= vec_max(vec_min(vcos, vec_splat_float(1.0f)), vec_splat_float(-1.0f));
+	vangle= vacosf(vcos);
+
+	/* dot */
+	vresult.v = (vec_splat_float(n[0])*gx +
+	             vec_splat_float(n[1])*gy +
+	             vec_splat_float(n[2])*gz)*vangle;
+
+	result= (vresult.f[0] + vresult.f[1] + vresult.f[2] + vresult.f[3])*(0.5f/(float)M_PI);
+	result= MAX2(result, 0.0f);
+
+	return result;
+}
+
+#endif
+
+/* SSE optimization, acos code doesn't work */
+
+#if 0
+
+#include <xmmintrin.h>
+
+static __m128 sse_approx_acos(__m128 x)
+{
+	/* needs a better approximation than taylor expansion of acos, since that
+	 * gives big erros for near 1.0 values, sqrt(2*x)*acos(1-x) should work
+	 * better, see http://www.tom.womack.net/projects/sse-fast-arctrig.html */
+
+	return _mm_set_ps1(1.0f);
+}
+
+static float ff_quad_form_factor(float *p, float *n, float *q0, float *q1, float *q2, float *q3)
+{
+	float r0[3], r1[3], r2[3], r3[3], g0[3], g1[3], g2[3], g3[3];
+	float a1, a2, a3, a4, dot1, dot2, dot3, dot4, result;
+	float fresult[4] __attribute__((aligned(16)));
+	__m128 qx, qy, qz, rx, ry, rz, rlen, srx, sry, srz, gx, gy, gz, glen, rcos, angle, aresult;
+
+	/* compute r */
+	qx = _mm_set_ps(q3[0], q2[0], q1[0], q0[0]);
+	qy = _mm_set_ps(q3[1], q2[1], q1[1], q0[1]);
+	qz = _mm_set_ps(q3[2], q2[2], q1[2], q0[2]);
+
+	rx = qx - _mm_set_ps1(p[0]);
+	ry = qy - _mm_set_ps1(p[1]);
+	rz = qz - _mm_set_ps1(p[2]);
+
+	/* normalize r */
+	rlen = _mm_rsqrt_ps(rx*rx + ry*ry + rz*rz + _mm_set_ps1(1e-16f));
+	rx = rx*rlen;
+	ry = ry*rlen;
+	rz = rz*rlen;
+
+	/* cross product */
+	srx = _mm_shuffle_ps(rx, rx, _MM_SHUFFLE(0,3,2,1));
+	sry = _mm_shuffle_ps(ry, ry, _MM_SHUFFLE(0,3,2,1));
+	srz = _mm_shuffle_ps(rz, rz, _MM_SHUFFLE(0,3,2,1));
+
+	gx = sry*rz - srz*ry;
+	gy = srz*rx - srx*rz;
+	gz = srx*ry - sry*rx;
+
+	/* normalize g */
+	glen = _mm_rsqrt_ps(gx*gx + gy*gy + gz*gz + _mm_set_ps1(1e-16f));
+	gx = gx*glen;
+	gy = gy*glen;
+	gz = gz*glen;
+
+	/* compute angle */
+	rcos = rx*srx + ry*sry + rz*srz;
+	rcos= _mm_max_ps(_mm_min_ps(rcos, _mm_set_ps1(1.0f)), _mm_set_ps1(-1.0f));
+
+	angle = sse_approx_cos(rcos);
+	aresult = (_mm_set_ps1(n[0])*gx + _mm_set_ps1(n[1])*gy + _mm_set_ps1(n[2])*gz)*angle;
+
+	/* sum together */
+    result= (fresult[0] + fresult[1] + fresult[2] + fresult[3])*(0.5f/(float)M_PI);
+	result= MAX2(result, 0.0f);
+
+	return result;
+}
+
+#endif
+
+static void ff_normalize(float n[3])
+{
+	float d;
+	
+	d= dot_v3v3(n, n);
+
+	if(d > 1.0e-35F) {
+		d= 1.0f/sqrtf(d);
+
+		n[0] *= d; 
+		n[1] *= d; 
+		n[2] *= d;
+	} 
+}
+
+static float ff_quad_form_factor(float *p, float *n, float *q0, float *q1, float *q2, float *q3)
+{
+	float r0[3], r1[3], r2[3], r3[3], g0[3], g1[3], g2[3], g3[3];
+	float a1, a2, a3, a4, dot1, dot2, dot3, dot4, result;
+
+	sub_v3_v3v3(r0, q0, p);
+	sub_v3_v3v3(r1, q1, p);
+	sub_v3_v3v3(r2, q2, p);
+	sub_v3_v3v3(r3, q3, p);
+
+	ff_normalize(r0);
+	ff_normalize(r1);
+	ff_normalize(r2);
+	ff_normalize(r3);
+
+	cross_v3_v3v3(g0, r1, r0); ff_normalize(g0);
+	cross_v3_v3v3(g1, r2, r1); ff_normalize(g1);
+	cross_v3_v3v3(g2, r3, r2); ff_normalize(g2);
+	cross_v3_v3v3(g3, r0, r3); ff_normalize(g3);
+
+	a1= saacosf(dot_v3v3(r0, r1));
+	a2= saacosf(dot_v3v3(r1, r2));
+	a3= saacosf(dot_v3v3(r2, r3));
+	a4= saacosf(dot_v3v3(r3, r0));
+
+	dot1= dot_v3v3(n, g0);
+	dot2= dot_v3v3(n, g1);
+	dot3= dot_v3v3(n, g2);
+	dot4= dot_v3v3(n, g3);
+
+	result= (a1*dot1 + a2*dot2 + a3*dot3 + a4*dot4)*0.5f/(float)M_PI;
+	result= MAX2(result, 0.0f);
+
+	return result;
+}
+
+float form_factor_hemi_poly(float p[3], float n[3], float v1[3], float v2[3], float v3[3], float v4[3])
+{

@@ Diff output truncated at 10240 characters. @@




More information about the Bf-blender-cvs mailing list