nuclear@54: #include <string.h>
nuclear@54: #include <assert.h>
nuclear@55: #include <limits.h>
nuclear@54: #include "rt.h"
nuclear@54: #include "ogl.h"
nuclear@54: #include "vector.h"
nuclear@54: #include "timer.h"
nuclear@54: 
nuclear@54: struct SurfPoint {
nuclear@54: 	float t;
nuclear@54: 	Vector3 pos, norm;
nuclear@54: 	const Face *face;
nuclear@54: };
nuclear@54: 
nuclear@54: static void trace_ray(float *pixel, const Ray &ray, int iter, float energy = 1.0f);
nuclear@54: static void shade(float *pixel, const Ray &ray, const SurfPoint &sp, int iter, float energy = 1.0f);
nuclear@54: static bool find_intersection(const Ray &ray, const Scene *scn, const KDNode *kd, SurfPoint *spret);
nuclear@54: static bool ray_aabb_test(const Ray &ray, const AABBox &aabb);
nuclear@54: static bool ray_triangle_test(const Ray &ray, const Face *face, SurfPoint *sp);
nuclear@54: static Vector3 calc_bary(const Vector3 &pt, const Face *face, const Vector3 &norm);
nuclear@54: static void transform(float *res, const float *v, const float *xform);
nuclear@54: static void transform_ray(Ray *ray, const float *xform, const float *invtrans_xform);
nuclear@54: 
nuclear@54: static int xsz, ysz;
nuclear@54: static float *fb;
nuclear@54: static unsigned int tex;
nuclear@54: static Scene *scn;
nuclear@54: static const Ray *prim_rays;
nuclear@55: static int max_iter;
nuclear@54: 
nuclear@55: static RenderStats *rstat;
nuclear@55: static int cur_ray_aabb_tests, cur_ray_triangle_tests;
nuclear@54: 
nuclear@54: bool init_dbg_renderer(int width, int height, Scene *scene, unsigned int texid)
nuclear@54: {
nuclear@54: 	try {
nuclear@54: 		fb = new float[3 * width * height];
nuclear@54: 	}
nuclear@54: 	catch(...) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	xsz = width;
nuclear@54: 	ysz = height;
nuclear@54: 	tex = texid;
nuclear@54: 	scn = scene;
nuclear@55: 
nuclear@55: 	rstat = (RenderStats*)get_render_stats();
nuclear@55: 
nuclear@54: 	return true;
nuclear@54: }
nuclear@54: 
nuclear@54: void destroy_dbg_renderer()
nuclear@54: {
nuclear@54: 	delete [] fb;
nuclear@54: 	delete [] prim_rays;
nuclear@54: }
nuclear@54: 
nuclear@54: void dbg_set_primary_rays(const Ray *rays)
nuclear@54: {
nuclear@54: 	prim_rays = rays;
nuclear@54: }
nuclear@54: 
nuclear@54: void dbg_render(const float *xform, const float *invtrans_xform, int num_threads)
nuclear@54: {
nuclear@54: 	unsigned long t0 = get_msec();
nuclear@54: 
nuclear@55: 	max_iter = get_render_option_int(ROPT_ITER);
nuclear@55: 
nuclear@55: 	// initialize render-stats
nuclear@55: 	memset(rstat, 0, sizeof *rstat);
nuclear@55: 	rstat->min_aabb_tests = rstat->min_triangle_tests = INT_MAX;
nuclear@55: 	rstat->max_aabb_tests = rstat->max_triangle_tests = 0;
nuclear@54: 
nuclear@54: 	int offs = 0;
nuclear@54: 	for(int i=0; i<ysz; i++) {
nuclear@54: 		for(int j=0; j<xsz; j++) {
nuclear@54: 			Ray ray = prim_rays[offs];
nuclear@54: 			transform_ray(&ray, xform, invtrans_xform);
nuclear@54: 
nuclear@55: 			cur_ray_aabb_tests = cur_ray_triangle_tests = 0;
nuclear@55: 
nuclear@55: 			trace_ray(fb + offs * 3, ray, max_iter, 1.0);
nuclear@54: 			offs++;
nuclear@55: 
nuclear@55: 			// update stats as needed
nuclear@55: 			if(cur_ray_aabb_tests < rstat->min_aabb_tests) {
nuclear@55: 				rstat->min_aabb_tests = cur_ray_aabb_tests;
nuclear@55: 			}
nuclear@55: 			if(cur_ray_aabb_tests > rstat->max_aabb_tests) {
nuclear@55: 				rstat->max_aabb_tests = cur_ray_aabb_tests;
nuclear@55: 			}
nuclear@55: 			if(cur_ray_triangle_tests < rstat->min_triangle_tests) {
nuclear@55: 				rstat->min_triangle_tests = cur_ray_triangle_tests;
nuclear@55: 			}
nuclear@55: 			if(cur_ray_triangle_tests > rstat->max_triangle_tests) {
nuclear@55: 				rstat->max_triangle_tests = cur_ray_triangle_tests;
nuclear@55: 			}
nuclear@55: 			rstat->prim_rays++;
nuclear@55: 			rstat->aabb_tests += cur_ray_aabb_tests;
nuclear@55: 			rstat->triangle_tests += cur_ray_triangle_tests;
nuclear@54: 		}
nuclear@54: 	}
nuclear@54: 
nuclear@55: 	unsigned long t1 = get_msec();
nuclear@55: 
nuclear@54: 	glPushAttrib(GL_TEXTURE_BIT);
nuclear@54: 	glBindTexture(GL_TEXTURE_2D, tex);
nuclear@54: 	glTexSubImage2D(GL_TEXTURE_2D, 0, 0, 0, xsz, ysz, GL_RGB, GL_FLOAT, fb);
nuclear@54: 	glPopAttrib();
nuclear@55: 	glFinish();
nuclear@54: 
nuclear@55: 	rstat->render_time = t1 - t0;
nuclear@55: 	rstat->tex_update_time = get_msec() - t1;
nuclear@55: 
nuclear@55: 	rstat->rays_cast = rstat->prim_rays + rstat->refl_rays + rstat->shadow_rays;
nuclear@55: 	rstat->rays_per_sec = 1000 * rstat->rays_cast / rstat->render_time;
nuclear@55: 	rstat->avg_aabb_tests = (float)rstat->aabb_tests / (float)rstat->rays_cast;
nuclear@55: 	rstat->avg_triangle_tests = (float)rstat->triangle_tests / (float)rstat->rays_cast;
nuclear@54: }
nuclear@54: 
nuclear@54: static void trace_ray(float *pixel, const Ray &ray, int iter, float energy)
nuclear@54: {
nuclear@54: 	SurfPoint sp;
nuclear@54: 
nuclear@54: 	if(find_intersection(ray, scn, scn->kdtree, &sp)) {
nuclear@54: 		shade(pixel, ray, sp, iter, energy);
nuclear@54: 	} else {
nuclear@54: 		pixel[0] = pixel[1] = pixel[2] = 0.05f;
nuclear@54: 	}
nuclear@54: }
nuclear@54: 
nuclear@54: #define MAX(a, b)	((a) > (b) ? (a) : (b))
nuclear@54: 
nuclear@54: static void shade(float *pixel, const Ray &ray, const SurfPoint &sp, int iter, float energy)
nuclear@54: {
nuclear@54: 	const Material *mat = scn->get_materials() + sp.face->matid;
nuclear@54: 
nuclear@54: 	bool cast_shadows = get_render_option_bool(ROPT_SHAD);
nuclear@54: 	Vector3 raydir(ray.dir);
nuclear@54: 	Vector3 norm = sp.norm;
nuclear@54: 
nuclear@54: 	if(dot(raydir, norm) >= 0.0) {
nuclear@54: 		norm = -norm;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	float dcol[3] = {0, 0, 0};
nuclear@54: 	float scol[3] = {0, 0, 0};
nuclear@54: 
nuclear@54: 	for(int i=0; i<scn->get_num_lights(); i++) {
nuclear@54: 		Vector3 lpos(scn->lights[i].pos);
nuclear@54: 		Vector3 ldir = lpos - sp.pos;
nuclear@54: 
nuclear@54: 		Ray shadowray;
nuclear@54: 		shadowray.origin[0] = sp.pos.x;
nuclear@54: 		shadowray.origin[1] = sp.pos.y;
nuclear@54: 		shadowray.origin[2] = sp.pos.z;
nuclear@54: 		shadowray.dir[0] = ldir.x;
nuclear@54: 		shadowray.dir[1] = ldir.y;
nuclear@54: 		shadowray.dir[2] = ldir.z;
nuclear@54: 
nuclear@54: 		if(!cast_shadows || !find_intersection(shadowray, scn, scn->kdtree, 0)) {
nuclear@55: 			rstat->brdf_evals++;
nuclear@54: 
nuclear@54: 			ldir.normalize();
nuclear@54: 
nuclear@54: 			Vector3 vdir = -raydir / RAY_MAG;
nuclear@54: 			Vector3 vref = reflect(vdir, norm);
nuclear@54: 
nuclear@54: 			float ndotl = dot(ldir, norm);
nuclear@54: 			float diff = MAX(ndotl, 0.0f);
nuclear@54: 
nuclear@54: 			dcol[0] += mat->kd[0] * diff;
nuclear@54: 			dcol[1] += mat->kd[1] * diff;
nuclear@54: 			dcol[2] += mat->kd[2] * diff;
nuclear@54: 
nuclear@54: 			float ldotvr = dot(ldir, vref);
nuclear@54: 			float spec = pow(MAX(ldotvr, 0.0f), mat->spow);
nuclear@54: 
nuclear@54: 			scol[0] += mat->ks[0] * spec;
nuclear@54: 			scol[1] += mat->ks[1] * spec;
nuclear@54: 			scol[2] += mat->ks[2] * spec;
nuclear@54: 		}
nuclear@55: 
nuclear@55: 		if(cast_shadows) {
nuclear@55: 			rstat->shadow_rays++;
nuclear@55: 		}
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	float refl_color[3];
nuclear@54: 	refl_color[0] = mat->ks[0] * mat->kr;
nuclear@54: 	refl_color[1] = mat->ks[1] * mat->kr;
nuclear@54: 	refl_color[2] = mat->ks[2] * mat->kr;
nuclear@54: 
nuclear@54: 	energy *= (refl_color[0] + refl_color[1] + refl_color[2]) / 3.0;
nuclear@54: 	if(iter >= 0 && energy > MIN_ENERGY) {
nuclear@54: 		Vector3 rdir = reflect(-raydir, norm);
nuclear@54: 
nuclear@54: 		Ray refl;
nuclear@54: 		refl.origin[0] = sp.pos.x;
nuclear@54: 		refl.origin[1] = sp.pos.y;
nuclear@54: 		refl.origin[2] = sp.pos.z;
nuclear@54: 		refl.dir[0] = rdir.x;
nuclear@54: 		refl.dir[1] = rdir.y;
nuclear@54: 		refl.dir[2] = rdir.z;
nuclear@54: 
nuclear@54: 		float rcol[3];
nuclear@54: 		trace_ray(rcol, refl, iter - 1, energy);
nuclear@54: 		scol[0] += rcol[0] * mat->ks[0] * mat->kr;
nuclear@54: 		scol[1] += rcol[1] * mat->ks[1] * mat->kr;
nuclear@54: 		scol[2] += rcol[2] * mat->ks[2] * mat->kr;
nuclear@55: 
nuclear@55: 		rstat->refl_rays++;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	pixel[0] = dcol[0] + scol[0];
nuclear@54: 	pixel[1] = dcol[1] + scol[1];
nuclear@54: 	pixel[2] = dcol[2] + scol[2];
nuclear@54: }
nuclear@54: 
nuclear@54: static bool find_intersection(const Ray &ray, const Scene *scn, const KDNode *kd, SurfPoint *spret)
nuclear@54: {
nuclear@54: 	if(!ray_aabb_test(ray, kd->aabb)) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	SurfPoint sp, sptmp;
nuclear@54: 	if(!spret) {
nuclear@54: 		spret = &sptmp;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	spret->t = RAY_MAG;
nuclear@54: 	spret->face = 0;
nuclear@54: 
nuclear@54: 	if(kd->left) {
nuclear@54: 		assert(kd->right);
nuclear@54: 
nuclear@54: 		bool found = find_intersection(ray, scn, kd->left, spret);
nuclear@54: 		if(find_intersection(ray, scn, kd->right, &sp)) {
nuclear@54: 			if(!found || sp.t < spret->t) {
nuclear@54: 				*spret = sp;
nuclear@54: 			}
nuclear@54: 			found = true;
nuclear@54: 		}
nuclear@54: 		return found;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	const Face *faces = scn->get_face_buffer();
nuclear@54: 
nuclear@54: 	for(size_t i=0; i<kd->face_idx.size(); i++) {
nuclear@54: 		if(ray_triangle_test(ray, faces + kd->face_idx[i], &sp) && sp.t < spret->t) {
nuclear@54: 			*spret = sp;
nuclear@54: 		}
nuclear@54: 	}
nuclear@54: 	return spret->face != 0;
nuclear@54: }
nuclear@54: 
nuclear@54: static bool ray_aabb_test(const Ray &ray, const AABBox &aabb)
nuclear@54: {
nuclear@55: 	cur_ray_aabb_tests++;
nuclear@55: 
nuclear@54: 	if(ray.origin[0] >= aabb.min[0] && ray.origin[1] >= aabb.min[1] && ray.origin[2] >= aabb.min[2] &&
nuclear@54: 			ray.origin[0] < aabb.max[0] && ray.origin[1] < aabb.max[1] && ray.origin[2] < aabb.max[2]) {
nuclear@54: 		return true;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	float bbox[][3] = {
nuclear@54: 		{aabb.min[0], aabb.min[1], aabb.min[2]},
nuclear@54: 		{aabb.max[0], aabb.max[1], aabb.max[2]}
nuclear@54: 	};
nuclear@54: 
nuclear@54: 	int xsign = (int)(ray.dir[0] < 0.0);
nuclear@54: 	float invdirx = 1.0 / ray.dir[0];
nuclear@54: 	float tmin = (bbox[xsign][0] - ray.origin[0]) * invdirx;
nuclear@54: 	float tmax = (bbox[1 - xsign][0] - ray.origin[0]) * invdirx;
nuclear@54: 
nuclear@54: 	int ysign = (int)(ray.dir[1] < 0.0);
nuclear@54: 	float invdiry = 1.0 / ray.dir[1];
nuclear@54: 	float tymin = (bbox[ysign][1] - ray.origin[1]) * invdiry;
nuclear@54: 	float tymax = (bbox[1 - ysign][1] - ray.origin[1]) * invdiry;
nuclear@54: 
nuclear@54: 	if(tmin > tymax || tymin > tmax) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	if(tymin > tmin) tmin = tymin;
nuclear@54: 	if(tymax < tmax) tmax = tymax;
nuclear@54: 
nuclear@54: 	int zsign = (int)(ray.dir[2] < 0.0);
nuclear@54: 	float invdirz = 1.0 / ray.dir[2];
nuclear@54: 	float tzmin = (bbox[zsign][2] - ray.origin[2]) * invdirz;
nuclear@54: 	float tzmax = (bbox[1 - zsign][2] - ray.origin[2]) * invdirz;
nuclear@54: 
nuclear@54: 	if(tmin > tzmax || tzmin > tmax) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	return tmin < 1.0 && tmax > 0.0;
nuclear@54: 
nuclear@54: }
nuclear@54: 
nuclear@54: static bool ray_triangle_test(const Ray &ray, const Face *face, SurfPoint *sp)
nuclear@54: {
nuclear@55: 	cur_ray_triangle_tests++;
nuclear@55: 
nuclear@54: 	Vector3 origin = ray.origin;
nuclear@54: 	Vector3 dir = ray.dir;
nuclear@54: 	Vector3 norm = face->normal;
nuclear@54: 
nuclear@54: 	float ndotdir = dot(dir, norm);
nuclear@54: 
nuclear@54: 	if(fabs(ndotdir) <= EPSILON) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	Vector3 pt = face->v[0].pos;
nuclear@54: 	Vector3 vec = pt - origin;
nuclear@54: 
nuclear@54: 	float ndotvec = dot(norm, vec);
nuclear@54: 	float t = ndotvec / ndotdir;
nuclear@54: 
nuclear@54: 	if(t < EPSILON || t > 1.0) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 	pt = origin + dir * t;
nuclear@54: 
nuclear@54: 
nuclear@54: 	Vector3 bc = calc_bary(pt, face, norm);
nuclear@54: 	float bc_sum = bc.x + bc.y + bc.z;
nuclear@54: 
nuclear@54: 	if(bc_sum < 1.0 - EPSILON || bc_sum > 1.0 + EPSILON) {
nuclear@54: 		return false;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	Vector3 n0(face->v[0].normal);
nuclear@54: 	Vector3 n1(face->v[1].normal);
nuclear@54: 	Vector3 n2(face->v[2].normal);
nuclear@54: 
nuclear@54: 	sp->t = t;
nuclear@54: 	sp->pos = pt;
nuclear@54: 	sp->norm = n0 * bc.x + n1 * bc.y + n2 * bc.z;
nuclear@54: 	sp->norm.normalize();
nuclear@54: 	sp->face = face;
nuclear@54: 	return true;
nuclear@54: }
nuclear@54: 
nuclear@54: static Vector3 calc_bary(const Vector3 &pt, const Face *face, const Vector3 &norm)
nuclear@54: {
nuclear@54: 	Vector3 bc(0.0f, 0.0f, 0.0f);
nuclear@54: 
nuclear@54: 	Vector3 v1 = Vector3(face->v[1].pos) - Vector3(face->v[0].pos);
nuclear@54: 	Vector3 v2 = Vector3(face->v[2].pos) - Vector3(face->v[0].pos);
nuclear@54: 	Vector3 xv1v2 = cross(v1, v2);
nuclear@54: 
nuclear@54: 	float area = fabs(dot(xv1v2, norm)) * 0.5;
nuclear@54: 	if(area < EPSILON) {
nuclear@54: 		return bc;
nuclear@54: 	}
nuclear@54: 
nuclear@54: 	Vector3 pv0 = face->v[0].pos - pt;
nuclear@54: 	Vector3 pv1 = face->v[1].pos - pt;
nuclear@54: 	Vector3 pv2 = face->v[2].pos - pt;
nuclear@54: 
nuclear@54: 	// calculate the area of each sub-triangle
nuclear@54: 	Vector3 x12 = cross(pv1, pv2);
nuclear@54: 	Vector3 x20 = cross(pv2, pv0);
nuclear@54: 	Vector3 x01 = cross(pv0, pv1);
nuclear@54: 
nuclear@54: 	float a0 = fabs(dot(x12, norm)) * 0.5;
nuclear@54: 	float a1 = fabs(dot(x20, norm)) * 0.5;
nuclear@54: 	float a2 = fabs(dot(x01, norm)) * 0.5;
nuclear@54: 
nuclear@54: 	bc.x = a0 / area;
nuclear@54: 	bc.y = a1 / area;
nuclear@54: 	bc.z = a2 / area;
nuclear@54: 	return bc;
nuclear@54: 
nuclear@54: }
nuclear@54: 
nuclear@54: static void transform(float *res, const float *v, const float *xform)
nuclear@54: {
nuclear@54: 	float tmp[3];
nuclear@54: 	tmp[0] = v[0] * xform[0] + v[1] * xform[4] + v[2] * xform[8] + xform[12];
nuclear@54: 	tmp[1] = v[0] * xform[1] + v[1] * xform[5] + v[2] * xform[9] + xform[13];
nuclear@54: 	tmp[2] = v[0] * xform[2] + v[1] * xform[6] + v[2] * xform[10] + xform[14];
nuclear@54: 	memcpy(res, tmp, sizeof tmp);
nuclear@54: }
nuclear@54: 
nuclear@54: static void transform_ray(Ray *ray, const float *xform, const float *invtrans_xform)
nuclear@54: {
nuclear@54: 	transform(ray->origin, ray->origin, xform);
nuclear@54: 	transform(ray->dir, ray->dir, invtrans_xform);
nuclear@54: }