1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/src/dbgray.cc	Fri Sep 10 16:47:00 2010 +0100
     1.3 @@ -0,0 +1,342 @@
     1.4 +#include <string.h>
     1.5 +#include <assert.h>
     1.6 +#include "rt.h"
     1.7 +#include "ogl.h"
     1.8 +#include "vector.h"
     1.9 +#include "timer.h"
    1.10 +
    1.11 +struct SurfPoint {
    1.12 +	float t;
    1.13 +	Vector3 pos, norm;
    1.14 +	const Face *face;
    1.15 +};
    1.16 +
    1.17 +static void trace_ray(float *pixel, const Ray &ray, int iter, float energy = 1.0f);
    1.18 +static void shade(float *pixel, const Ray &ray, const SurfPoint &sp, int iter, float energy = 1.0f);
    1.19 +static bool find_intersection(const Ray &ray, const Scene *scn, const KDNode *kd, SurfPoint *spret);
    1.20 +static bool ray_aabb_test(const Ray &ray, const AABBox &aabb);
    1.21 +static bool ray_triangle_test(const Ray &ray, const Face *face, SurfPoint *sp);
    1.22 +static Vector3 calc_bary(const Vector3 &pt, const Face *face, const Vector3 &norm);
    1.23 +static void transform(float *res, const float *v, const float *xform);
    1.24 +static void transform_ray(Ray *ray, const float *xform, const float *invtrans_xform);
    1.25 +
    1.26 +static int xsz, ysz;
    1.27 +static float *fb;
    1.28 +static unsigned int tex;
    1.29 +static Scene *scn;
    1.30 +static const Ray *prim_rays;
    1.31 +
    1.32 +
    1.33 +bool init_dbg_renderer(int width, int height, Scene *scene, unsigned int texid)
    1.34 +{
    1.35 +	try {
    1.36 +		fb = new float[3 * width * height];
    1.37 +	}
    1.38 +	catch(...) {
    1.39 +		return false;
    1.40 +	}
    1.41 +
    1.42 +	xsz = width;
    1.43 +	ysz = height;
    1.44 +	tex = texid;
    1.45 +	scn = scene;
    1.46 +	return true;
    1.47 +}
    1.48 +
    1.49 +void destroy_dbg_renderer()
    1.50 +{
    1.51 +	delete [] fb;
    1.52 +	delete [] prim_rays;
    1.53 +}
    1.54 +
    1.55 +void dbg_set_primary_rays(const Ray *rays)
    1.56 +{
    1.57 +	prim_rays = rays;
    1.58 +}
    1.59 +
    1.60 +void dbg_render(const float *xform, const float *invtrans_xform, int num_threads)
    1.61 +{
    1.62 +	unsigned long t0 = get_msec();
    1.63 +
    1.64 +	int iter = get_render_option_int(ROPT_ITER);
    1.65 +
    1.66 +	int offs = 0;
    1.67 +	for(int i=0; i<ysz; i++) {
    1.68 +		for(int j=0; j<xsz; j++) {
    1.69 +			Ray ray = prim_rays[offs];
    1.70 +			transform_ray(&ray, xform, invtrans_xform);
    1.71 +
    1.72 +			trace_ray(fb + offs * 3, ray, iter, 1.0);
    1.73 +			offs++;
    1.74 +		}
    1.75 +	}
    1.76 +
    1.77 +	glPushAttrib(GL_TEXTURE_BIT);
    1.78 +	glBindTexture(GL_TEXTURE_2D, tex);
    1.79 +	glTexSubImage2D(GL_TEXTURE_2D, 0, 0, 0, xsz, ysz, GL_RGB, GL_FLOAT, fb);
    1.80 +	glPopAttrib();
    1.81 +
    1.82 +	printf("rendered in %lu msec\n", get_msec() - t0);
    1.83 +}
    1.84 +
    1.85 +static void trace_ray(float *pixel, const Ray &ray, int iter, float energy)
    1.86 +{
    1.87 +	SurfPoint sp;
    1.88 +
    1.89 +	if(find_intersection(ray, scn, scn->kdtree, &sp)) {
    1.90 +		shade(pixel, ray, sp, iter, energy);
    1.91 +	} else {
    1.92 +		pixel[0] = pixel[1] = pixel[2] = 0.05f;
    1.93 +	}
    1.94 +}
    1.95 +
    1.96 +#define MAX(a, b)	((a) > (b) ? (a) : (b))
    1.97 +
    1.98 +static void shade(float *pixel, const Ray &ray, const SurfPoint &sp, int iter, float energy)
    1.99 +{
   1.100 +	const Material *mat = scn->get_materials() + sp.face->matid;
   1.101 +
   1.102 +	bool cast_shadows = get_render_option_bool(ROPT_SHAD);
   1.103 +	Vector3 raydir(ray.dir);
   1.104 +	Vector3 norm = sp.norm;
   1.105 +
   1.106 +	if(dot(raydir, norm) >= 0.0) {
   1.107 +		norm = -norm;
   1.108 +	}
   1.109 +
   1.110 +	float dcol[3] = {0, 0, 0};
   1.111 +	float scol[3] = {0, 0, 0};
   1.112 +
   1.113 +	for(int i=0; i<scn->get_num_lights(); i++) {
   1.114 +		Vector3 lpos(scn->lights[i].pos);
   1.115 +		Vector3 ldir = lpos - sp.pos;
   1.116 +
   1.117 +		Ray shadowray;
   1.118 +		shadowray.origin[0] = sp.pos.x;
   1.119 +		shadowray.origin[1] = sp.pos.y;
   1.120 +		shadowray.origin[2] = sp.pos.z;
   1.121 +		shadowray.dir[0] = ldir.x;
   1.122 +		shadowray.dir[1] = ldir.y;
   1.123 +		shadowray.dir[2] = ldir.z;
   1.124 +
   1.125 +		if(!cast_shadows || !find_intersection(shadowray, scn, scn->kdtree, 0)) {
   1.126 +
   1.127 +			ldir.normalize();
   1.128 +
   1.129 +			Vector3 vdir = -raydir / RAY_MAG;
   1.130 +			Vector3 vref = reflect(vdir, norm);
   1.131 +
   1.132 +			float ndotl = dot(ldir, norm);
   1.133 +			float diff = MAX(ndotl, 0.0f);
   1.134 +
   1.135 +			dcol[0] += mat->kd[0] * diff;
   1.136 +			dcol[1] += mat->kd[1] * diff;
   1.137 +			dcol[2] += mat->kd[2] * diff;
   1.138 +
   1.139 +			float ldotvr = dot(ldir, vref);
   1.140 +			float spec = pow(MAX(ldotvr, 0.0f), mat->spow);
   1.141 +
   1.142 +			scol[0] += mat->ks[0] * spec;
   1.143 +			scol[1] += mat->ks[1] * spec;
   1.144 +			scol[2] += mat->ks[2] * spec;
   1.145 +		}
   1.146 +	}
   1.147 +
   1.148 +	float refl_color[3];
   1.149 +	refl_color[0] = mat->ks[0] * mat->kr;
   1.150 +	refl_color[1] = mat->ks[1] * mat->kr;
   1.151 +	refl_color[2] = mat->ks[2] * mat->kr;
   1.152 +
   1.153 +	energy *= (refl_color[0] + refl_color[1] + refl_color[2]) / 3.0;
   1.154 +	if(iter >= 0 && energy > MIN_ENERGY) {
   1.155 +		Vector3 rdir = reflect(-raydir, norm);
   1.156 +
   1.157 +		Ray refl;
   1.158 +		refl.origin[0] = sp.pos.x;
   1.159 +		refl.origin[1] = sp.pos.y;
   1.160 +		refl.origin[2] = sp.pos.z;
   1.161 +		refl.dir[0] = rdir.x;
   1.162 +		refl.dir[1] = rdir.y;
   1.163 +		refl.dir[2] = rdir.z;
   1.164 +
   1.165 +		float rcol[3];
   1.166 +		trace_ray(rcol, refl, iter - 1, energy);
   1.167 +		scol[0] += rcol[0] * mat->ks[0] * mat->kr;
   1.168 +		scol[1] += rcol[1] * mat->ks[1] * mat->kr;
   1.169 +		scol[2] += rcol[2] * mat->ks[2] * mat->kr;
   1.170 +	}
   1.171 +
   1.172 +	pixel[0] = dcol[0] + scol[0];
   1.173 +	pixel[1] = dcol[1] + scol[1];
   1.174 +	pixel[2] = dcol[2] + scol[2];
   1.175 +}
   1.176 +
   1.177 +static bool find_intersection(const Ray &ray, const Scene *scn, const KDNode *kd, SurfPoint *spret)
   1.178 +{
   1.179 +	if(!ray_aabb_test(ray, kd->aabb)) {
   1.180 +		return false;
   1.181 +	}
   1.182 +
   1.183 +	SurfPoint sp, sptmp;
   1.184 +	if(!spret) {
   1.185 +		spret = &sptmp;
   1.186 +	}
   1.187 +
   1.188 +	spret->t = RAY_MAG;
   1.189 +	spret->face = 0;
   1.190 +
   1.191 +	if(kd->left) {
   1.192 +		assert(kd->right);
   1.193 +
   1.194 +		bool found = find_intersection(ray, scn, kd->left, spret);
   1.195 +		if(find_intersection(ray, scn, kd->right, &sp)) {
   1.196 +			if(!found || sp.t < spret->t) {
   1.197 +				*spret = sp;
   1.198 +			}
   1.199 +			found = true;
   1.200 +		}
   1.201 +		return found;
   1.202 +	}
   1.203 +
   1.204 +	const Face *faces = scn->get_face_buffer();
   1.205 +
   1.206 +	for(size_t i=0; i<kd->face_idx.size(); i++) {
   1.207 +		if(ray_triangle_test(ray, faces + kd->face_idx[i], &sp) && sp.t < spret->t) {
   1.208 +			*spret = sp;
   1.209 +		}
   1.210 +	}
   1.211 +	return spret->face != 0;
   1.212 +}
   1.213 +
   1.214 +static bool ray_aabb_test(const Ray &ray, const AABBox &aabb)
   1.215 +{
   1.216 +	if(ray.origin[0] >= aabb.min[0] && ray.origin[1] >= aabb.min[1] && ray.origin[2] >= aabb.min[2] &&
   1.217 +			ray.origin[0] < aabb.max[0] && ray.origin[1] < aabb.max[1] && ray.origin[2] < aabb.max[2]) {
   1.218 +		return true;
   1.219 +	}
   1.220 +
   1.221 +	float bbox[][3] = {
   1.222 +		{aabb.min[0], aabb.min[1], aabb.min[2]},
   1.223 +		{aabb.max[0], aabb.max[1], aabb.max[2]}
   1.224 +	};
   1.225 +
   1.226 +	int xsign = (int)(ray.dir[0] < 0.0);
   1.227 +	float invdirx = 1.0 / ray.dir[0];
   1.228 +	float tmin = (bbox[xsign][0] - ray.origin[0]) * invdirx;
   1.229 +	float tmax = (bbox[1 - xsign][0] - ray.origin[0]) * invdirx;
   1.230 +
   1.231 +	int ysign = (int)(ray.dir[1] < 0.0);
   1.232 +	float invdiry = 1.0 / ray.dir[1];
   1.233 +	float tymin = (bbox[ysign][1] - ray.origin[1]) * invdiry;
   1.234 +	float tymax = (bbox[1 - ysign][1] - ray.origin[1]) * invdiry;
   1.235 +
   1.236 +	if(tmin > tymax || tymin > tmax) {
   1.237 +		return false;
   1.238 +	}
   1.239 +
   1.240 +	if(tymin > tmin) tmin = tymin;
   1.241 +	if(tymax < tmax) tmax = tymax;
   1.242 +
   1.243 +	int zsign = (int)(ray.dir[2] < 0.0);
   1.244 +	float invdirz = 1.0 / ray.dir[2];
   1.245 +	float tzmin = (bbox[zsign][2] - ray.origin[2]) * invdirz;
   1.246 +	float tzmax = (bbox[1 - zsign][2] - ray.origin[2]) * invdirz;
   1.247 +
   1.248 +	if(tmin > tzmax || tzmin > tmax) {
   1.249 +		return false;
   1.250 +	}
   1.251 +
   1.252 +	return tmin < 1.0 && tmax > 0.0;
   1.253 +
   1.254 +}
   1.255 +
   1.256 +static bool ray_triangle_test(const Ray &ray, const Face *face, SurfPoint *sp)
   1.257 +{
   1.258 +	Vector3 origin = ray.origin;
   1.259 +	Vector3 dir = ray.dir;
   1.260 +	Vector3 norm = face->normal;
   1.261 +
   1.262 +	float ndotdir = dot(dir, norm);
   1.263 +
   1.264 +	if(fabs(ndotdir) <= EPSILON) {
   1.265 +		return false;
   1.266 +	}
   1.267 +
   1.268 +	Vector3 pt = face->v[0].pos;
   1.269 +	Vector3 vec = pt - origin;
   1.270 +
   1.271 +	float ndotvec = dot(norm, vec);
   1.272 +	float t = ndotvec / ndotdir;
   1.273 +
   1.274 +	if(t < EPSILON || t > 1.0) {
   1.275 +		return false;
   1.276 +	}
   1.277 +	pt = origin + dir * t;
   1.278 +
   1.279 +
   1.280 +	Vector3 bc = calc_bary(pt, face, norm);
   1.281 +	float bc_sum = bc.x + bc.y + bc.z;
   1.282 +
   1.283 +	if(bc_sum < 1.0 - EPSILON || bc_sum > 1.0 + EPSILON) {
   1.284 +		return false;
   1.285 +	}
   1.286 +
   1.287 +	Vector3 n0(face->v[0].normal);
   1.288 +	Vector3 n1(face->v[1].normal);
   1.289 +	Vector3 n2(face->v[2].normal);
   1.290 +
   1.291 +	sp->t = t;
   1.292 +	sp->pos = pt;
   1.293 +	sp->norm = n0 * bc.x + n1 * bc.y + n2 * bc.z;
   1.294 +	sp->norm.normalize();
   1.295 +	sp->face = face;
   1.296 +	return true;
   1.297 +}
   1.298 +
   1.299 +static Vector3 calc_bary(const Vector3 &pt, const Face *face, const Vector3 &norm)
   1.300 +{
   1.301 +	Vector3 bc(0.0f, 0.0f, 0.0f);
   1.302 +
   1.303 +	Vector3 v1 = Vector3(face->v[1].pos) - Vector3(face->v[0].pos);
   1.304 +	Vector3 v2 = Vector3(face->v[2].pos) - Vector3(face->v[0].pos);
   1.305 +	Vector3 xv1v2 = cross(v1, v2);
   1.306 +
   1.307 +	float area = fabs(dot(xv1v2, norm)) * 0.5;
   1.308 +	if(area < EPSILON) {
   1.309 +		return bc;
   1.310 +	}
   1.311 +
   1.312 +	Vector3 pv0 = face->v[0].pos - pt;
   1.313 +	Vector3 pv1 = face->v[1].pos - pt;
   1.314 +	Vector3 pv2 = face->v[2].pos - pt;
   1.315 +
   1.316 +	// calculate the area of each sub-triangle
   1.317 +	Vector3 x12 = cross(pv1, pv2);
   1.318 +	Vector3 x20 = cross(pv2, pv0);
   1.319 +	Vector3 x01 = cross(pv0, pv1);
   1.320 +
   1.321 +	float a0 = fabs(dot(x12, norm)) * 0.5;
   1.322 +	float a1 = fabs(dot(x20, norm)) * 0.5;
   1.323 +	float a2 = fabs(dot(x01, norm)) * 0.5;
   1.324 +
   1.325 +	bc.x = a0 / area;
   1.326 +	bc.y = a1 / area;
   1.327 +	bc.z = a2 / area;
   1.328 +	return bc;
   1.329 +
   1.330 +}
   1.331 +
   1.332 +static void transform(float *res, const float *v, const float *xform)
   1.333 +{
   1.334 +	float tmp[3];
   1.335 +	tmp[0] = v[0] * xform[0] + v[1] * xform[4] + v[2] * xform[8] + xform[12];
   1.336 +	tmp[1] = v[0] * xform[1] + v[1] * xform[5] + v[2] * xform[9] + xform[13];
   1.337 +	tmp[2] = v[0] * xform[2] + v[1] * xform[6] + v[2] * xform[10] + xform[14];
   1.338 +	memcpy(res, tmp, sizeof tmp);
   1.339 +}
   1.340 +
   1.341 +static void transform_ray(Ray *ray, const float *xform, const float *invtrans_xform)
   1.342 +{
   1.343 +	transform(ray->origin, ray->origin, xform);
   1.344 +	transform(ray->dir, ray->dir, invtrans_xform);
   1.345 +}