workaround for divide by zero bug in octree pyrit
authorRadek Brich <radek.brich@devl.cz>
Fri, 28 Mar 2008 00:53:20 +0100 (2008-03-27)
branchpyrit
changeset 49 558fde7da82a
parent 48 a4913301c626
child 50 14a727b70d07
workaround for divide by zero bug in octree enable oversampling for DefaultSampler sampler and camera performance enhancements
TODO
config.mk
include/sampler.h
include/scene.h
src/octree.cc
src/raytracer.cc
src/sampler.cc
src/scene.cc
--- a/TODO	Wed Mar 26 17:03:38 2008 +0100
+++ b/TODO	Fri Mar 28 00:53:20 2008 +0100
@@ -4,13 +4,11 @@
    Inside-out transition is ignored as every second intersection is ignored to avoid duplicit
    intersect points. This solution should be replaced by moving ray origin a little forward, which
    should work as well and without side-effects.
- * bug in octree - black lines in center of image
 
 Future Plans
 ============
  * changing main ray tracing algoritm to more flexible architecture with a Sampler object:
-  - update all demos for new Sampler/Camera architecture
-  - enhance Sampler to support oversampling and subsampling as before
+  - enhance Sampler to support subsampling
   - rewrite pthreads
  * namespace
  * kd-tree:
@@ -25,6 +23,11 @@
  * implement efficient AABB-ray intersection using Plucker coordinates
  * generalization: Camera "shader" (ray generator), surface shader and maybe light & background shaders
 
+Performance Weak Points
+=======================
+ * Raytracer::render -- allocating and deallocating every sample
+
+
 New Classes?
 ============
 
--- a/config.mk	Wed Mar 26 17:03:38 2008 +0100
+++ b/config.mk	Fri Mar 28 00:53:20 2008 +0100
@@ -1,3 +1,5 @@
+CXX=ccache g++
+
 CCFLAGS=-g -Wall -Wno-write-strings -fno-strict-aliasing -I$(ROOT)/include
 LDFLAGS=
 
--- a/include/sampler.h	Wed Mar 26 17:03:38 2008 +0100
+++ b/include/sampler.h	Fri Mar 28 00:53:20 2008 +0100
@@ -72,7 +72,7 @@
 class DefaultSample: public Sample
 {
 	friend class DefaultSampler;
-	int sx,sy;
+	int sx,sy,osa_samp;
 };
 
 /**
--- a/include/scene.h	Wed Mar 26 17:03:38 2008 +0100
+++ b/include/scene.h	Fri Mar 28 00:53:20 2008 +0100
@@ -63,13 +63,13 @@
 {
 public:
 	Vector3 eye, p, u, v;
-	Float f;
+	Float f,F;
 
-	Camera(): eye(0,0,10), p(0,0,-1), u(-1,0,0), v(0,1,0), f(3.14/4.0) {};
+	Camera(): eye(0,0,10), p(0,0,-1), u(-1,0,0), v(0,1,0), f(3.14/4.0), F(0.5/f) {};
 	Camera(const Vector3 &C, const Vector3 &ap, const Vector3 &au, const Vector3 &av):
-		eye(C), p(ap), u(au), v(av), f(3.14/4.0) {};
+		eye(C), p(ap), u(au), v(av), f(3.14/4.0), F(0.5/f) {};
 	void setEye(const Vector3 &aeye) { eye = aeye; };
-	void setFocalLength(const Float af) { f = af; };
+	void setFocalLength(const Float af) { f = af; F = 0.5/f; };
 	void rotate(const Quaternion &q);
 	void move(const Float fw, const Float left, const Float up);
 
--- a/src/octree.cc	Wed Mar 26 17:03:38 2008 +0100
+++ b/src/octree.cc	Fri Mar 28 00:53:20 2008 +0100
@@ -171,7 +171,7 @@
 
 	int a = 0;
 	Vector3 ro(ray.o);
-	Vector3 rdir(1.0/ray.dir.x, 1.0/ray.dir.y, 1.0/ray.dir.z);
+	Vector3 rdir(ray.dir);
 
 	if (rdir.x < 0.0)
 	{
@@ -192,6 +192,13 @@
 		a |= 1;
 	}
 
+	if (rdir.x == 0.0) rdir.x = Eps;
+	if (rdir.y == 0.0) rdir.y = Eps;
+	if (rdir.z == 0.0) rdir.z = Eps;
+	rdir.x = 1.0/rdir.x;
+	rdir.y = 1.0/rdir.y;
+	rdir.z = 1.0/rdir.z;
+
 	tx0 = (bbox.L.x - ro.x) * rdir.x;
 	tx1 = (bbox.H.x - ro.x) * rdir.x;
 	ty0 = (bbox.L.y - ro.y) * rdir.y;
@@ -199,7 +206,7 @@
 	tz0 = (bbox.L.z - ro.z) * rdir.z;
 	tz1 = (bbox.H.z - ro.z) * rdir.z;
 
-	if (max3(tx0,ty0,tz0) > min3(tx1,ty1,tz1))
+	if (max3(tx0,ty0,tz0) >= min3(tx1,ty1,tz1))
 		return NULL;
 
 	node = root;
--- a/src/raytracer.cc	Wed Mar 26 17:03:38 2008 +0100
+++ b/src/raytracer.cc	Fri Mar 28 00:53:20 2008 +0100
@@ -238,54 +238,6 @@
 	}
 }
 
-static inline void samplepixel(Colour &c, Vector3 &dir, RenderrowData* d, const int &oversample)
-{
-	if (oversample <= 0)
-	{
-		// no oversampling
-		dir.normalize();
-		Ray ray(d->eye, dir);
-		c = d->rt->raytrace(ray, 0, NULL);
-	}
-	else
-	if (oversample <= 3)
-	{
-		// grid oversampling
-		c = Colour(0,0,0);
-		static const int gridsamples[] = {5,9,16};
-		static const Float osa5x[] = {0.0, -0.4, +0.4, +0.4, -0.4};
-		static const Float osa5y[] = {0.0, -0.4, -0.4, +0.4, +0.4};
-		static const Float osa9x[] = {-0.34,  0.00, +0.34,
-			-0.34,  0.00, +0.34, -0.34,  0.00, +0.34};
-		static const Float osa9y[] = {-0.34, -0.34, -0.34,
-				0.00,  0.00,  0.00, +0.34, +0.34, +0.34};
-		static const Float osa16x[] = {-0.375, -0.125, +0.125, +0.375,
-			-0.375, -0.125, +0.125, +0.375, -0.375, -0.125, +0.125, +0.375,
-			-0.375, -0.125, +0.125, +0.375};
-		static const Float osa16y[] = {-0.375, -0.375, -0.375, -0.375,
-			-0.125, -0.125, -0.125, -0.125, +0.125, +0.125, +0.125, +0.125,
-			+0.375, +0.375, +0.375, +0.375};
-		static const Float *osaSx[] = {osa5x, osa9x, osa16x};
-		static const Float *osaSy[] = {osa5y, osa9y, osa16y};
-		const int samples = gridsamples[oversample-1];
-		const Float *osax = osaSx[oversample-1];
-		const Float *osay = osaSy[oversample-1];
-		for (int i = 0; i < samples; i++)
-		{
-			Vector3 tmpdir = dir + osax[i]*d->dx + osay[i]*d->dy;
-			tmpdir.normalize();
-			Ray ray(d->eye, tmpdir);
-			c += d->rt->raytrace(ray, 0, NULL);
-		}
-		c = c * (1.0/samples);
-	}
-	else
-	{
-		// stochastic oversampling
-		// ...todo
-	}
-}
-
 #if 0
 static void *renderrow(void *data)
 {
@@ -468,7 +420,11 @@
 
 			delete prev;
 			prev = sample;
+			sampnum--;
+			if ((sampnum % 1000) == 0)
+				dbgmsg(2, "\b\b\b\b\b\b\b\b%8d", sampnum);
 		}
+		dbgmsg(2, "\n");
 	}
 
 	// wait for workers
--- a/src/sampler.cc	Wed Mar 26 17:03:38 2008 +0100
+++ b/src/sampler.cc	Fri Mar 28 00:53:20 2008 +0100
@@ -36,10 +36,22 @@
 
 int DefaultSampler::initSampleSet()
 {
+	static const int gridsamples[] = {1,5,9,16};
+	const int samples = gridsamples[oversample];
 	if ( phase == 0 )
 	{
+		cout << "phase 1" << endl;
 		phase++;
-		return w*h;
+		return w*h*samples;
+	}
+	else if ( phase == 1 && oversample )
+	{
+		cout << "phase 2" << endl;
+		phase = -1;
+		Float *buf;
+		for (buf = buffer; buf != buffer + w*h*3; buf++)
+			*buf = *buf * (1.0/samples);
+		return 0;
 	}
 	else
 	{
@@ -51,31 +63,79 @@
 Sample* DefaultSampler::nextSample(Sample *prev)
 {
 	DefaultSample *s = new DefaultSample;
-	if (prev)
+
+	/* grid oversampling */
+	static const int gridsamples[] = {1,5,9,16};
+	static const Float osa5x[] = {0.0, -0.4, +0.4, +0.4, -0.4};
+	static const Float osa5y[] = {0.0, -0.4, -0.4, +0.4, +0.4};
+	static const Float osa9x[] = {-0.34,  0.00, +0.34,
+		-0.34,  0.00, +0.34, -0.34,  0.00, +0.34};
+	static const Float osa9y[] = {-0.34, -0.34, -0.34,
+			0.00,  0.00,  0.00, +0.34, +0.34, +0.34};
+	static const Float osa16x[] = {-0.375, -0.125, +0.125, +0.375,
+		-0.375, -0.125, +0.125, +0.375, -0.375, -0.125, +0.125, +0.375,
+		-0.375, -0.125, +0.125, +0.375};
+	static const Float osa16y[] = {-0.375, -0.375, -0.375, -0.375,
+		-0.125, -0.125, -0.125, -0.125, +0.125, +0.125, +0.125, +0.125,
+		+0.375, +0.375, +0.375, +0.375};
+	static const Float *osaSx[] = {NULL, osa5x, osa9x, osa16x};
+	static const Float *osaSy[] = {NULL, osa5y, osa9y, osa16y};
+	const int samples = gridsamples[oversample];
+	const Float *osax = osaSx[oversample];
+	const Float *osay = osaSy[oversample];
+
+	if (!prev)
 	{
-		DefaultSample *sp = static_cast<DefaultSample*>(prev);
-		s->sx = sp->sx + 1;
-		s->sy = sp->sy;
-		if (s->sx >= w)
-		{
-			s->sx = 0;
-			s->sy++;
-		}
-		if (s->sy >= h)
-		{
-			delete s;
-			return NULL;
-		}
-		s->x = (Float)s->sx/h - (Float)w/h/2.0;
-		s->y = (Float)s->sy/h - 0.5;
-	}
-	else
-	{
+		// first sample
 		s->x = -(Float)w/h/2.0;
 		s->y = -0.5;
 		s->sx = 0;
 		s->sy = 0;
+		s->osa_samp = 0;
 	}
+	else
+	{
+		DefaultSample *sp = static_cast<DefaultSample*>(prev);
+
+		s->osa_samp = sp->osa_samp + 1;
+
+		if (oversample && oversample <= 3 && s->osa_samp < samples)
+		{
+			s->sx = sp->sx;
+			s->sy = sp->sy;
+			s->x = osax[s->osa_samp]/h + (Float)s->sx/h - (Float)w/h/2.0;
+			s->y = osay[s->osa_samp]/h + (Float)s->sy/h - 0.5;
+		}
+		else
+		{
+			s->sx = sp->sx + 1;
+			s->sy = sp->sy;
+			if (s->sx >= w)
+			{
+				s->sx = 0;
+				s->sy++;
+			}
+			if (s->sy >= h)
+			{
+				delete s;
+				return NULL;
+			}
+			s->x = (Float)s->sx/h - (Float)w/h/2.0;
+			s->y = (Float)s->sy/h - 0.5;
+			s->osa_samp = 0;
+		}
+	}
+
+	if (s->osa_samp == 0 && oversample && oversample <= 3)
+	{
+		s->x += osax[0]/h;
+		s->y += osay[0]/h;
+		Float *buf = buffer + 3*(s->sy*w + s->sx);
+		*(buf++) = 0;
+		*(buf++) = 0;
+		*(buf++) = 0;
+	}
+
 	return s;
 }
 
@@ -83,8 +143,17 @@
 {
 	DefaultSample *sp = static_cast<DefaultSample*>(samp);
 	Float *buf = buffer + 3*(sp->sy*w + sp->sx);
-	*(buf++) = col.r;
-	*(buf++) = col.g;
-	*(buf++) = col.b;
+	if (oversample)
+	{
+		*(buf+0) += col.r;
+		*(buf+1) += col.g;
+		*(buf+2) += col.b;
+	}
+	else
+	{
+		*(buf++) = col.r;
+		*(buf++) = col.g;
+		*(buf++) = col.b;
+	}
 }
 
--- a/src/scene.cc	Wed Mar 26 17:03:38 2008 +0100
+++ b/src/scene.cc	Fri Mar 28 00:53:20 2008 +0100
@@ -77,7 +77,7 @@
 
 Ray Camera::makeRay(Sample *samp)
 {
-	Vector3 dir = p + u*(0.5/f)*samp->x + v*(-0.5/f)*samp->y;
+	Vector3 dir = p + (u*samp->x - v*samp->y)*F;
 	dir.normalize();
 	return Ray(eye, dir);
 }