Cycles render engine, initial commit. This is the engine itself, blender modifications and build instructions will follow later.

Cycles uses code from some great open source projects, many thanks them:

* BVH building and traversal code from NVidia's "Understanding the Efficiency of Ray Traversal on GPUs":
http://code.google.com/p/understanding-the-efficiency-of-ray-traversal-on-gpus/
* Open Shading Language for a large part of the shading system:
http://code.google.com/p/openshadinglanguage/
* Blender for procedural textures and a few other nodes.
* Approximate Catmull Clark subdivision from NVidia Mesh tools:
http://code.google.com/p/nvidia-mesh-tools/
* Sobol direction vectors from:
http://web.maths.unsw.edu.au/~fkuo/sobol/
* Film response functions from:
http://www.cs.columbia.edu/CAVE/software/softlib/dorf.php

10 files changed, 1908 insertions, 0 deletions

diff --git a/intern/cycles/bvh/CMakeLists.txt b/intern/cycles/bvh/CMakeLists.txt
new file mode 100644
index 00000000000..c934cded6da
--- /dev/null
+++ b/intern/cycles/bvh/CMakeLists.txt
@@ -0,0 +1,18 @@
+
+INCLUDE_DIRECTORIES(. ../kernel ../kernel/svm ../render ../util ../device)
+
+SET(sources
+	bvh.cpp
+	bvh_build.cpp
+	bvh_node.cpp
+	bvh_sort.cpp)
+
+SET(headers
+	bvh.h
+	bvh_build.h
+	bvh_node.h
+	bvh_params.h
+	bvh_sort.h)
+
+ADD_LIBRARY(bvh ${sources} ${headers})
+
diff --git a/intern/cycles/bvh/bvh.cpp b/intern/cycles/bvh/bvh.cpp
new file mode 100644
index 00000000000..ab230794774
--- /dev/null
+++ b/intern/cycles/bvh/bvh.cpp
@@ -0,0 +1,661 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "mesh.h"
+#include "object.h"
+#include "scene.h"
+
+#include "bvh.h"
+#include "bvh_build.h"
+#include "bvh_node.h"
+#include "bvh_params.h"
+
+#include "util_cache.h"
+#include "util_debug.h"
+#include "util_foreach.h"
+#include "util_progress.h"
+#include "util_types.h"
+
+CCL_NAMESPACE_BEGIN
+
+/* Pack Utility */
+
+struct BVHStackEntry
+{
+	const BVHNode *node;
+	int idx;
+
+	BVHStackEntry(const BVHNode* n = 0, int i = 0)
+	: node(n), idx(i)
+	{
+	}
+
+	int encodeIdx() const
+	{
+		return (node->is_leaf())? ~idx: idx;
+	}
+};
+
+/* BVH */
+
+BVH::BVH(const BVHParams& params_, const vector<Object*>& objects_)
+: params(params_), objects(objects_)
+{
+}
+
+BVH *BVH::create(const BVHParams& params, const vector<Object*>& objects)
+{
+	if(params.use_qbvh)
+		return new QBVH(params, objects);
+	else
+		return new RegularBVH(params, objects);
+}
+
+/* Cache */
+
+bool BVH::cache_read(CacheData& key)
+{
+	key.add(&params, sizeof(params));
+
+	foreach(Object *ob, objects) {
+		key.add(ob->mesh->verts);
+		key.add(ob->mesh->triangles);
+	}
+
+	CacheData value;
+
+	if(Cache::global.lookup(key, value)) {
+		value.read(pack.root_index);
+
+		value.read(pack.nodes);
+		value.read(pack.object_node);
+		value.read(pack.tri_woop);
+		value.read(pack.prim_index);
+		value.read(pack.prim_object);
+		value.read(pack.is_leaf);
+
+		return true;
+	}
+
+	return false;
+}
+
+void BVH::cache_write(CacheData& key)
+{
+	CacheData value;
+
+	value.add(pack.root_index);
+
+	value.add(pack.nodes);
+	value.add(pack.object_node);
+	value.add(pack.tri_woop);
+	value.add(pack.prim_index);
+	value.add(pack.prim_object);
+	value.add(pack.is_leaf);
+
+	Cache::global.insert(key, value);
+}
+
+/* Building */
+
+void BVH::build(Progress& progress)
+{
+	progress.set_substatus("Building BVH");
+
+	/* cache read */
+	CacheData key("bvh");
+
+	if(params.use_cache) {
+		progress.set_substatus("Looking in BVH cache");
+
+		if(cache_read(key))
+			return;
+	}
+
+	/* build nodes */
+	vector<int> prim_index;
+	vector<int> prim_object;
+
+	BVHBuild bvh_build(objects, prim_index, prim_object, params, progress);
+	BVHNode *root = bvh_build.run();
+
+	if(progress.get_cancel()) {
+		if(root) root->deleteSubtree();
+		return;
+	}
+
+	/* todo: get rid of this copy */
+	pack.prim_index = prim_index;
+	pack.prim_object = prim_object;
+
+	/* compute SAH */
+	if(!params.top_level)
+		pack.SAH = root->computeSubtreeSAHCost(params);
+
+	if(progress.get_cancel()) {
+		root->deleteSubtree();
+		return;
+	}
+
+	/* pack triangles */
+	progress.set_substatus("Packing BVH triangles");
+	pack_triangles();
+
+	if(progress.get_cancel()) {
+		root->deleteSubtree();
+		return;
+	}
+
+	/* pack nodes */
+	progress.set_substatus("Packing BVH nodes");
+	array<int> tmp_prim_object = pack.prim_object;
+	pack_nodes(tmp_prim_object, root);
+	
+	/* free build nodes */
+	root->deleteSubtree();
+
+	if(progress.get_cancel()) return;
+
+	/* cache write */
+	if(params.use_cache) {
+		progress.set_substatus("Writing BVH cache");
+		cache_write(key);
+	}
+}
+
+/* Refitting */
+
+void BVH::refit(Progress& progress)
+{
+	progress.set_substatus("Packing BVH triangles");
+	pack_triangles();
+
+	if(progress.get_cancel()) return;
+
+	progress.set_substatus("Refitting BVH nodes");
+	refit_nodes();
+}
+
+/* Triangles */
+
+void BVH::pack_triangle(int idx, float4 woop[3])
+{
+	/* create Woop triangle */
+	int tob = pack.prim_object[idx];
+	const Mesh *mesh = objects[tob]->mesh;
+	int tidx = pack.prim_index[idx];
+	const int *vidx = mesh->triangles[tidx].v;
+	const float3* vpos = &mesh->verts[0];
+	float3 v0 = vpos[vidx[0]];
+	float3 v1 = vpos[vidx[1]];
+	float3 v2 = vpos[vidx[2]];
+
+	float3 r0 = v0 - v2;
+	float3 r1 = v1 - v2;
+	float3 r2 = cross(r0, r1);
+
+	if(dot(r0, r0) == 0.0f || dot(r1, r1) == 0.0f || dot(r2, r2) == 0.0f) {
+		/* degenerate */
+		woop[0] = make_float4(0.0f, 0.0f, 0.0f, 0.0f);
+		woop[1] = make_float4(0.0f, 0.0f, 0.0f, 0.0f);
+		woop[2] = make_float4(0.0f, 0.0f, 0.0f, 0.0f);
+	}
+	else {
+		Transform t = make_transform(
+			r0.x, r1.x, r2.x, v2.x,
+			r0.y, r1.y, r2.y, v2.y,
+			r0.z, r1.z, r2.z, v2.z,
+			0.0f, 0.0f, 0.0f, 1.0f);
+
+		t = transform_inverse(t);
+
+		woop[0] = make_float4(t.z.x, t.z.y, t.z.z, -t.z.w);
+		woop[1] = make_float4(t.x.x, t.x.y, t.x.z, t.x.w);
+		woop[2] = make_float4(t.y.x, t.y.y, t.y.z, t.y.w);
+	}
+}
+
+void BVH::pack_triangles()
+{
+	int nsize = TRI_NODE_SIZE;
+	size_t tidx_size = pack.prim_index.size();
+
+	pack.tri_woop.clear();
+	pack.tri_woop.resize(tidx_size * nsize);
+
+	for(unsigned int i = 0; i < tidx_size; i++) {
+		if(pack.prim_index[i] != -1) {
+			float4 woop[3];
+
+			pack_triangle(i, woop);
+			memcpy(&pack.tri_woop[i * nsize], woop, sizeof(float4)*3);
+		}
+	}
+}
+
+/* Pack Instances */
+
+void BVH::pack_instances(size_t nodes_size)
+{
+	/* The BVH's for instances are built separately, but for traversal all
+	   BVH's are stored in global arrays. This function merges them into the
+	   top level BVH, adjusting indexes and offsets where appropriate. */
+	bool use_qbvh = params.use_qbvh;
+	size_t nsize = (use_qbvh)? BVH_QNODE_SIZE: BVH_NODE_SIZE;
+
+	/* adjust primitive index to point to the triangle in the global array, for
+	   meshes with transform applied and already in the top level BVH */
+	for(size_t i = 0; i < pack.prim_index.size(); i++)
+		if(pack.prim_index[i] != -1)
+			pack.prim_index[i] += objects[pack.prim_object[i]]->mesh->tri_offset;
+
+	/* track offsets of instanced BVH data in global array */
+	size_t tri_offset = pack.prim_index.size();
+	size_t nodes_offset = nodes_size;
+
+	/* clear array that gives the node indexes for instanced objects */
+	pack.object_node.clear();
+
+	/* reserve */
+	size_t prim_index_size = pack.prim_index.size();
+	size_t tri_woop_size = pack.tri_woop.size();
+
+	size_t pack_prim_index_offset = prim_index_size;
+	size_t pack_tri_woop_offset = tri_woop_size;
+	size_t pack_nodes_offset = nodes_size;
+	size_t object_offset = 0;
+
+	foreach(Object *ob, objects) {
+		Mesh *mesh = ob->mesh;
+		BVH *bvh = mesh->bvh;
+
+		if(!mesh->transform_applied) {
+			prim_index_size += bvh->pack.prim_index.size();
+			tri_woop_size += bvh->pack.tri_woop.size();
+			nodes_size += bvh->pack.nodes.size()*nsize;
+		}
+	}
+
+	pack.prim_index.resize(prim_index_size);
+	pack.prim_object.resize(prim_index_size);
+	pack.tri_woop.resize(tri_woop_size);
+	pack.nodes.resize(nodes_size);
+	pack.object_node.resize(objects.size());
+
+	int *pack_prim_index = &pack.prim_index[0];
+	int *pack_prim_object = &pack.prim_object[0];
+	float4 *pack_tri_woop = &pack.tri_woop[0];
+	int4 *pack_nodes = &pack.nodes[0];
+
+	/* merge */
+	foreach(Object *ob, objects) {
+		Mesh *mesh = ob->mesh;
+
+		/* if mesh transform is applied, that means it's already in the top
+		   level BVH, and we don't need to merge it in */
+		if(mesh->transform_applied) {
+			pack.object_node[object_offset++] = 0;
+			continue;
+		}
+
+		BVH *bvh = mesh->bvh;
+
+		int noffset = nodes_offset/nsize;
+		int mesh_tri_offset = mesh->tri_offset;
+
+		/* fill in node indexes for instances */
+		if(bvh->pack.is_leaf[0])
+			pack.object_node[object_offset++] = -noffset-1;
+		else
+			pack.object_node[object_offset++] = noffset;
+
+		/* merge primitive and object indexes */
+		{
+			size_t bvh_prim_index_size = bvh->pack.prim_index.size();
+			int *bvh_prim_index = &bvh->pack.prim_index[0];
+
+			for(size_t i = 0; i < bvh_prim_index_size; i++) {
+				pack_prim_index[pack_prim_index_offset] = bvh_prim_index[i] + mesh_tri_offset;
+				pack_prim_object[pack_prim_index_offset] = 0;  // unused for instances
+				pack_prim_index_offset++;
+			}
+		}
+
+		/* merge triangle intersection data */
+		{
+			memcpy(pack_tri_woop+pack_tri_woop_offset, &bvh->pack.tri_woop[0],
+				bvh->pack.tri_woop.size()*sizeof(float4));
+			pack_tri_woop_offset += bvh->pack.tri_woop.size();
+		}
+
+		/* merge nodes */
+		{
+			size_t nsize_bbox = (use_qbvh)? nsize-2: nsize-1;
+			int4 *bvh_nodes = &bvh->pack.nodes[0];
+			size_t bvh_nodes_size = bvh->pack.nodes.size(); 
+			int *bvh_is_leaf = &bvh->pack.is_leaf[0];
+
+			for(size_t i = 0, j = 0; i < bvh_nodes_size; i+=nsize, j++) {
+				memcpy(pack_nodes + pack_nodes_offset, bvh_nodes + i, nsize_bbox*sizeof(int4));
+
+				/* modify offsets into arrays */
+				int4 data = bvh_nodes[i + nsize_bbox];
+
+				if(bvh_is_leaf[j]) {
+					data.x += tri_offset;
+					data.y += tri_offset;
+				}
+				else {
+					data.x += (data.x < 0)? -noffset: noffset;
+					data.y += (data.y < 0)? -noffset: noffset;
+
+					if(use_qbvh) {
+						data.z += (data.z < 0)? -noffset: noffset;
+						data.w += (data.w < 0)? -noffset: noffset;
+					}
+				}
+
+				pack_nodes[pack_nodes_offset + nsize_bbox] = data;
+
+				if(use_qbvh)
+					pack_nodes[pack_nodes_offset + nsize_bbox+1] = bvh_nodes[i + nsize_bbox+1];
+
+				pack_nodes_offset += nsize;
+			}
+		}
+
+		nodes_offset += bvh->pack.nodes.size();
+		tri_offset += bvh->pack.prim_index.size();
+	}
+}
+
+/* Regular BVH */
+
+RegularBVH::RegularBVH(const BVHParams& params_, const vector<Object*>& objects_)
+: BVH(params_, objects_)
+{
+}
+
+void RegularBVH::pack_leaf(const BVHStackEntry& e, const LeafNode *leaf)
+{
+	if(leaf->num_triangles() == 1 && pack.prim_index[leaf->m_lo] == -1)
+		/* object */
+		pack_node(e.idx, leaf->m_bounds, leaf->m_bounds, ~(leaf->m_lo), 0);
+	else
+		/* triangle */
+		pack_node(e.idx, leaf->m_bounds, leaf->m_bounds, leaf->m_lo, leaf->m_hi);
+}
+
+void RegularBVH::pack_inner(const BVHStackEntry& e, const BVHStackEntry& e0, const BVHStackEntry& e1)
+{
+	pack_node(e.idx, e0.node->m_bounds, e1.node->m_bounds, e0.encodeIdx(), e1.encodeIdx());
+}
+
+void RegularBVH::pack_node(int idx, const BoundBox& b0, const BoundBox& b1, int c0, int c1)
+{
+	int4 data[BVH_NODE_SIZE] =
+	{
+		make_int4(__float_as_int(b0.min.x), __float_as_int(b0.max.x), __float_as_int(b0.min.y), __float_as_int(b0.max.y)),
+		make_int4(__float_as_int(b1.min.x), __float_as_int(b1.max.x), __float_as_int(b1.min.y), __float_as_int(b1.max.y)),
+		make_int4(__float_as_int(b0.min.z), __float_as_int(b0.max.z), __float_as_int(b1.min.z), __float_as_int(b1.max.z)),
+		make_int4(c0, c1, 0, 0)
+	};
+
+	memcpy(&pack.nodes[idx * BVH_NODE_SIZE], data, sizeof(int4)*BVH_NODE_SIZE);
+}
+
+void RegularBVH::pack_nodes(const array<int>& prims, const BVHNode *root)
+{
+	size_t node_size = root->getSubtreeSize(BVH_STAT_NODE_COUNT);
+
+	/* resize arrays */
+	pack.nodes.clear();
+	pack.is_leaf.clear();
+	pack.is_leaf.resize(node_size);
+
+	/* for top level BVH, first merge existing BVH's so we know the offsets */
+	if(params.top_level)
+		pack_instances(node_size*BVH_NODE_SIZE);
+	else
+		pack.nodes.resize(node_size*BVH_NODE_SIZE);
+
+	int nextNodeIdx = 0;
+
+	vector<BVHStackEntry> stack;
+	stack.push_back(BVHStackEntry(root, nextNodeIdx++));
+
+	while(stack.size()) {
+		BVHStackEntry e = stack.back();
+		stack.pop_back();
+
+		pack.is_leaf[e.idx] = e.node->is_leaf();
+
+		if(e.node->is_leaf()) {
+			/* leaf node */
+			const LeafNode* leaf = reinterpret_cast<const LeafNode*>(e.node);
+			pack_leaf(e, leaf);
+		}
+		else {
+			/* innner node */
+			stack.push_back(BVHStackEntry(e.node->get_child(0), nextNodeIdx++));
+			stack.push_back(BVHStackEntry(e.node->get_child(1), nextNodeIdx++));
+
+			pack_inner(e, stack[stack.size()-2], stack[stack.size()-1]);
+		}
+	}
+
+	/* root index to start traversal at, to handle case of single leaf node */
+	pack.root_index = (pack.is_leaf[0])? -1: 0;
+}
+
+void RegularBVH::refit_nodes()
+{
+	assert(!params.top_level);
+
+	BoundBox bbox;
+	refit_node(0, pack.is_leaf[0], bbox);
+}
+
+void RegularBVH::refit_node(int idx, bool leaf, BoundBox& bbox)
+{
+	int4 *data = &pack.nodes[idx*4];
+
+	int c0 = data[3].x;
+	int c1 = data[3].y;
+
+	if(leaf) {
+		/* refit leaf node */
+		for(int tri = c0; tri < c1; tri++) {
+			int tidx = pack.prim_index[tri];
+			int tob = pack.prim_object[tri];
+			Object *ob = objects[tob];
+
+			if(tidx == -1) {
+				/* object instance */
+				bbox.grow(ob->bounds);
+			}
+			else {
+				/* triangles */
+				const Mesh *mesh = ob->mesh;
+				int tri_offset = (params.top_level)? mesh->tri_offset: 0;
+				const int *vidx = mesh->triangles[tidx - tri_offset].v;
+				const float3 *vpos = &mesh->verts[0];
+
+				bbox.grow(vpos[vidx[0]]);
+				bbox.grow(vpos[vidx[1]]);
+				bbox.grow(vpos[vidx[2]]);
+			}
+		}
+
+		pack_node(idx, bbox, bbox, c0, c1);
+	}
+	else {
+		/* refit inner node, set bbox from children */
+		BoundBox bbox0, bbox1;
+
+		refit_node((c0 < 0)? -c0-1: c0, (c0 < 0), bbox0);
+		refit_node((c1 < 0)? -c1-1: c1, (c1 < 0), bbox1);
+
+		bbox.grow(bbox0);
+		bbox.grow(bbox1);
+
+		pack_node(idx, bbox0, bbox1, c0, c1);
+	}
+}
+
+/* QBVH */
+
+QBVH::QBVH(const BVHParams& params_, const vector<Object*>& objects_)
+: BVH(params_, objects_)
+{
+	params.use_qbvh = true;
+}
+
+void QBVH::pack_leaf(const BVHStackEntry& e, const LeafNode *leaf)
+{
+	float4 data[BVH_QNODE_SIZE];
+
+	memset(data, 0, sizeof(data));
+
+	if(leaf->num_triangles() == 1 && pack.prim_index[leaf->m_lo] == -1) {
+		/* object */
+		data[6].x = __int_as_float(~(leaf->m_lo));
+		data[6].y = __int_as_float(0);
+	}
+	else {
+		/* triangle */
+		data[6].x = __int_as_float(leaf->m_lo);
+		data[6].y = __int_as_float(leaf->m_hi);
+	}
+
+	memcpy(&pack.nodes[e.idx * BVH_QNODE_SIZE], data, sizeof(float4)*BVH_QNODE_SIZE);
+}
+
+void QBVH::pack_inner(const BVHStackEntry& e, const BVHStackEntry *en, int num)
+{
+	float4 data[BVH_QNODE_SIZE];
+
+	for(int i = 0; i < num; i++) {
+		float3 bb_min = en[i].node->m_bounds.min;
+		float3 bb_max = en[i].node->m_bounds.max;
+
+		data[0][i] = bb_min.x;
+		data[1][i] = bb_max.x;
+		data[2][i] = bb_min.y;
+		data[3][i] = bb_max.y;
+		data[4][i] = bb_min.z;
+		data[5][i] = bb_max.z;
+
+		data[6][i] = __int_as_float(en[i].encodeIdx());
+		data[7][i] = 0.0f;
+	}
+
+	for(int i = num; i < 4; i++) {
+		data[0][i] = 0.0f;
+		data[1][i] = 0.0f;
+		data[2][i] = 0.0f;
+
+		data[3][i] = 0.0f;
+		data[4][i] = 0.0f;
+		data[5][i] = 0.0f;
+
+		data[6][i] = __int_as_float(0);
+		data[7][i] = 0.0f;
+	}
+
+	memcpy(&pack.nodes[e.idx * BVH_QNODE_SIZE], data, sizeof(float4)*BVH_QNODE_SIZE);
+}
+
+/* Quad SIMD Nodes */
+
+void QBVH::pack_nodes(const array<int>& prims, const BVHNode *root)
+{
+	size_t node_size = root->getSubtreeSize(BVH_STAT_NODE_COUNT);
+
+	/* resize arrays */
+	pack.nodes.clear();
+	pack.is_leaf.clear();
+	pack.is_leaf.resize(node_size);
+
+	/* for top level BVH, first merge existing BVH's so we know the offsets */
+	if(params.top_level)
+		pack_instances(node_size*BVH_QNODE_SIZE);
+	else
+		pack.nodes.resize(node_size*BVH_QNODE_SIZE);
+
+	int nextNodeIdx = 0;
+
+	vector<BVHStackEntry> stack;
+	stack.push_back(BVHStackEntry(root, nextNodeIdx++));
+
+	while(stack.size()) {
+		BVHStackEntry e = stack.back();
+		stack.pop_back();
+
+		pack.is_leaf[e.idx] = e.node->is_leaf();
+
+		if(e.node->is_leaf()) {
+			/* leaf node */
+			const LeafNode* leaf = reinterpret_cast<const LeafNode*>(e.node);
+			pack_leaf(e, leaf);
+		}
+		else {
+			/* inner node */
+			const BVHNode *node = e.node;
+			const BVHNode *node0 = node->get_child(0);
+			const BVHNode *node1 = node->get_child(1);
+
+			/* collect nodes */
+			const BVHNode *nodes[4];
+			int numnodes = 0;
+
+			if(node0->is_leaf()) {
+				nodes[numnodes++] = node0;
+			}
+			else {
+				nodes[numnodes++] = node0->get_child(0);
+				nodes[numnodes++] = node0->get_child(1);
+			}
+
+			if(node1->is_leaf()) {
+				nodes[numnodes++] = node1;
+			}
+			else {
+				nodes[numnodes++] = node1->get_child(0);
+				nodes[numnodes++] = node1->get_child(1);
+			}
+
+			/* push entries on the stack */
+			for(int i = 0; i < numnodes; i++)
+				stack.push_back(BVHStackEntry(nodes[i], nextNodeIdx++));
+
+			/* set node */
+			pack_inner(e, &stack[stack.size()-numnodes], numnodes);
+		}
+	}
+
+	/* root index to start traversal at, to handle case of single leaf node */
+	pack.root_index = (pack.is_leaf[0])? -1: 0;
+}
+
+void QBVH::refit_nodes()
+{
+	assert(0); /* todo */
+}
+
+CCL_NAMESPACE_END
+
diff --git a/intern/cycles/bvh/bvh.h b/intern/cycles/bvh/bvh.h
new file mode 100644
index 00000000000..acc25291da3
--- /dev/null
+++ b/intern/cycles/bvh/bvh.h
@@ -0,0 +1,152 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef __BVH_H__
+#define __BVH_H__
+
+#include "bvh_params.h"
+
+#include "util_types.h"
+#include "util_vector.h"
+
+CCL_NAMESPACE_BEGIN
+
+class BVHNode;
+class BVHStackEntry;
+class BVHParams;
+class BoundBox;
+class CacheData;
+class LeafNode;
+class Object;
+class Progress;
+
+#define BVH_NODE_SIZE	4
+#define BVH_QNODE_SIZE	8
+#define BVH_ALIGN		4096
+#define TRI_NODE_SIZE	3
+
+/* Packed BVH
+ *
+ * BVH stored as it will be used for traversal on the rendering device. */
+
+struct PackedBVH {
+	/* BVH nodes storage, one node is 4x int4, and contains two bounding boxes,
+	   and child, triangle or object indexes dependening on the node type */
+	array<int4> nodes; 
+	/* object index to BVH node index mapping for instances */
+	array<int> object_node; 
+	/* precomputed triangle intersection data, one triangle is 4x float4 */
+	array<float4> tri_woop; 
+	/* mapping from BVH primitive index to true primitive index, as primitives
+	   may be duplicated due to spatial splits. -1 for instances. */
+	array<int> prim_index;
+	/* mapping from BVH primitive index, to the object id of that primitive. */
+	array<int> prim_object;
+	/* quick array to lookup if a node is a leaf, not used for traversal, only
+	   for instance BVH merging  */
+	array<int> is_leaf;
+
+	/* index of the root node. */
+	int root_index;
+
+	/* surface area heuristic, for building top level BVH */
+	float SAH;
+
+	PackedBVH()
+	{
+		root_index = 0;
+		SAH = 0.0f;
+	}
+};
+
+/* BVH */
+
+class BVH
+{
+public:
+	PackedBVH pack;
+	BVHParams params;
+	vector<Object*> objects;
+
+	static BVH *create(const BVHParams& params, const vector<Object*>& objects);
+
+	void build(Progress& progress);
+	void refit(Progress& progress);
+
+protected:
+	BVH(const BVHParams& params, const vector<Object*>& objects);
+
+	/* cache */
+	bool cache_read(CacheData& key);
+	void cache_write(CacheData& key);
+
+	/* triangles */
+	void pack_triangles();
+	void pack_triangle(int idx, float4 woop[3]);
+
+	/* merge instance BVH's */
+	void pack_instances(size_t nodes_size);
+
+	/* for subclasses to implement */
+	virtual void pack_nodes(const array<int>& prims, const BVHNode *root) = 0;
+	virtual void refit_nodes() = 0;
+};
+
+/* Regular BVH
+ *
+ * Typical BVH with each node having two children. */
+
+class RegularBVH : public BVH {
+protected:
+	/* constructor */
+	friend class BVH;
+	RegularBVH(const BVHParams& params, const vector<Object*>& objects);
+
+	/* pack */
+	void pack_nodes(const array<int>& prims, const BVHNode *root);
+	void pack_leaf(const BVHStackEntry& e, const LeafNode *leaf);
+	void pack_inner(const BVHStackEntry& e, const BVHStackEntry& e0, const BVHStackEntry& e1);
+	void pack_node(int idx, const BoundBox& b0, const BoundBox& b1, int c0, int c1);
+
+	/* refit */
+	void refit_nodes();
+	void refit_node(int idx, bool leaf, BoundBox& bbox);
+};
+
+/* QBVH
+ *
+ * Quad BVH, with each node having four children, to use with SIMD instructions. */
+
+class QBVH : public BVH {
+protected:
+	/* constructor */
+	friend class BVH;
+	QBVH(const BVHParams& params, const vector<Object*>& objects);
+
+	/* pack */
+	void pack_nodes(const array<int>& prims, const BVHNode *root);
+	void pack_leaf(const BVHStackEntry& e, const LeafNode *leaf);
+	void pack_inner(const BVHStackEntry& e, const BVHStackEntry *en, int num);
+
+	/* refit */
+	void refit_nodes();
+};
+
+CCL_NAMESPACE_END
+
+#endif /* __BVH_H__ */
+
diff --git a/intern/cycles/bvh/bvh_build.cpp b/intern/cycles/bvh/bvh_build.cpp
new file mode 100644
index 00000000000..6a9cc915f01
--- /dev/null
+++ b/intern/cycles/bvh/bvh_build.cpp
@@ -0,0 +1,545 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "bvh_build.h"
+#include "bvh_node.h"
+#include "bvh_params.h"
+#include "bvh_sort.h"
+
+#include "mesh.h"
+#include "object.h"
+#include "scene.h"
+
+#include "util_algorithm.h"
+#include "util_foreach.h"
+#include "util_progress.h"
+#include "util_time.h"
+
+CCL_NAMESPACE_BEGIN
+
+/* Constructor / Destructor */
+
+BVHBuild::BVHBuild(const vector<Object*>& objects_,
+	vector<int>& prim_index_, vector<int>& prim_object_,
+	const BVHParams& params_, Progress& progress_)
+: objects(objects_),
+  prim_index(prim_index_),
+  prim_object(prim_object_),
+  params(params_),
+  progress(progress_),
+  progress_start_time(0.0)
+{
+	spatial_min_overlap = 0.0f;
+	progress_num_duplicates = 0;
+}
+
+BVHBuild::~BVHBuild()
+{
+}
+
+/* Adding References */
+
+void BVHBuild::add_reference_mesh(NodeSpec& root, Mesh *mesh, int i)
+{
+	for(uint j = 0; j < mesh->triangles.size(); j++) {
+		Mesh::Triangle t = mesh->triangles[j];
+		Reference ref;
+
+		ref.prim_index = j;
+		ref.prim_object = i;
+
+		for(int k = 0; k < 3; k++) {
+			float3 pt = mesh->verts[t.v[k]];
+			ref.bounds.grow(pt);
+		}
+
+		references.push_back(ref);
+		root.bounds.grow(ref.bounds);
+	}
+}
+
+void BVHBuild::add_reference_object(NodeSpec& root, Object *ob, int i)
+{
+	Reference ref;
+
+	ref.prim_index = -1;
+	ref.prim_object = i;
+	ref.bounds = ob->bounds;
+
+	references.push_back(ref);
+	root.bounds.grow(ref.bounds);
+}
+
+void BVHBuild::add_references(NodeSpec& root)
+{
+	/* init root spec */
+	root.num = 0;
+	root.bounds = BoundBox();
+
+	/* add objects */
+	int i = 0;
+
+	foreach(Object *ob, objects) {
+		if(params.top_level) {
+			if(ob->mesh->transform_applied)
+				add_reference_mesh(root, ob->mesh, i);
+			else
+				add_reference_object(root, ob, i);
+		}
+		else
+			add_reference_mesh(root, ob->mesh, i);
+
+		i++;
+
+		if(progress.get_cancel()) return;
+	}
+
+	/* happens mostly on empty meshes */
+	if(!root.bounds.valid())
+		root.bounds.grow(make_float3(0.0f, 0.0f, 0.0f));
+
+	root.num = references.size();
+}
+
+/* Build */
+
+BVHNode* BVHBuild::run()
+{
+	NodeSpec root;
+
+	/* add references */
+	add_references(root);
+
+	if(progress.get_cancel()) return NULL;
+
+	/* init spatial splits */
+	if(params.top_level) /* todo: get rid of this */
+		params.use_spatial_split = false;
+
+	spatial_min_overlap = root.bounds.area() * params.spatial_split_alpha;
+	spatial_right_bounds.clear();
+	spatial_right_bounds.resize(max(root.num, (int)BVHParams::NUM_SPATIAL_BINS) - 1);
+
+	/* init progress updates */
+	progress_num_duplicates = 0;
+	progress_start_time = time_dt();
+
+	/* build recursively */
+	return build_node(root, 0, 0.0f, 1.0f);
+}
+
+void BVHBuild::progress_update(float progress_start, float progress_end)
+{
+	if(time_dt() - progress_start_time < 0.25f)
+		return;
+
+	float duplicates = (float)progress_num_duplicates/(float)references.size();
+	string msg = string_printf("Building BVH %.0f%%, duplicates %.0f%%",
+		progress_start*100.0f, duplicates*100.0f);
+
+	progress.set_substatus(msg);
+	progress_start_time = time_dt();
+}
+
+BVHNode* BVHBuild::build_node(const NodeSpec& spec, int level, float progress_start, float progress_end)
+{
+	/* progress update */
+	progress_update(progress_start, progress_end);
+	if(progress.get_cancel()) return NULL;
+
+	/* small enough or too deep => create leaf. */
+	if(spec.num <= params.min_leaf_size || level >= BVHParams::MAX_DEPTH)
+		return create_leaf_node(spec);
+
+	/* find split candidates. */
+	float area = spec.bounds.area();
+	float leafSAH = area * params.triangle_cost(spec.num);
+	float nodeSAH = area * params.node_cost(2);
+	ObjectSplit object = find_object_split(spec, nodeSAH);
+	SpatialSplit spatial;
+
+	if(params.use_spatial_split && level < BVHParams::MAX_SPATIAL_DEPTH) {
+		BoundBox overlap = object.left_bounds;
+		overlap.intersect(object.right_bounds);
+
+		if(overlap.area() >= spatial_min_overlap)
+			spatial = find_spatial_split(spec, nodeSAH);
+	}
+
+	/* leaf SAH is the lowest => create leaf. */
+	float minSAH = min(min(leafSAH, object.sah), spatial.sah);
+
+	if(minSAH == leafSAH && spec.num <= params.max_leaf_size)
+		return create_leaf_node(spec);
+
+	/* perform split. */
+	NodeSpec left, right;
+
+	if(params.use_spatial_split && minSAH == spatial.sah)
+		do_spatial_split(left, right, spec, spatial);
+	if(!left.num || !right.num)
+		do_object_split(left, right, spec, object);
+
+	/* create inner node. */
+	progress_num_duplicates += left.num + right.num - spec.num;
+
+	float progress_mid = lerp(progress_start, progress_end, (float)right.num / (float)(left.num + right.num));
+
+	BVHNode* rightNode = build_node(right, level + 1, progress_start, progress_mid);
+	if(progress.get_cancel()) {
+		if(rightNode) rightNode->deleteSubtree();
+		return NULL;
+	}
+
+	BVHNode* leftNode = build_node(left, level + 1, progress_mid, progress_end);
+	if(progress.get_cancel()) {
+		if(leftNode) leftNode->deleteSubtree();
+		return NULL;
+	}
+
+	return new InnerNode(spec.bounds, leftNode, rightNode);
+}
+
+BVHNode *BVHBuild::create_object_leaf_nodes(const Reference *ref, int num)
+{
+	if(num == 0) {
+		BoundBox bounds;
+		return new LeafNode(bounds, 0, 0);
+	}
+	else if(num == 1) {
+		prim_index.push_back(ref[0].prim_index);
+		prim_object.push_back(ref[0].prim_object);
+		return new LeafNode(ref[0].bounds, prim_index.size()-1, prim_index.size());
+	}
+	else {
+		int mid = num/2;
+		BVHNode *leaf0 = create_object_leaf_nodes(ref, mid); 
+		BVHNode *leaf1 = create_object_leaf_nodes(ref+mid, num-mid); 
+
+		BoundBox bounds;
+		bounds.grow(leaf0->m_bounds);
+		bounds.grow(leaf1->m_bounds);
+
+		return new InnerNode(bounds, leaf0, leaf1);
+	}
+}
+
+BVHNode* BVHBuild::create_leaf_node(const NodeSpec& spec)
+{
+	vector<int>& p_index = prim_index;
+	vector<int>& p_object = prim_object;
+	BoundBox bounds;
+	int num = 0;
+
+	for(int i = 0; i < spec.num; i++) {
+		if(references.back().prim_index != -1) {
+			p_index.push_back(references.back().prim_index);
+			p_object.push_back(references.back().prim_object);
+			bounds.grow(references.back().bounds);
+			references.pop_back();
+			num++;
+		}
+	}
+
+	BVHNode *leaf = NULL;
+	
+	if(num > 0) {
+		leaf = new LeafNode(bounds, p_index.size() - num, p_index.size());
+
+		if(num == spec.num)
+			return leaf;
+	}
+
+	/* while there may be multiple triangles in a leaf, for object primitives
+	 * we want them to be the only one, so we  */
+	int ob_num = spec.num - num;
+	BVHNode *oleaf = create_object_leaf_nodes(&references.back() - (ob_num - 1), ob_num);
+	for(int i = 0; i < ob_num; i++)
+		references.pop_back();
+	
+	if(leaf)
+		return new InnerNode(spec.bounds, leaf, oleaf);
+	else
+		return oleaf;
+}
+
+/* Object Split */
+
+BVHBuild::ObjectSplit BVHBuild::find_object_split(const NodeSpec& spec, float nodeSAH)
+{
+	ObjectSplit split;
+	const Reference *ref_ptr = &references[references.size() - spec.num];
+
+	for(int dim = 0; dim < 3; dim++) {
+		/* sort references */
+		bvh_reference_sort(references.size() - spec.num, references.size(), &references[0], dim);
+
+		/* sweep right to left and determine bounds. */
+		BoundBox right_bounds;
+
+		for(int i = spec.num - 1; i > 0; i--) {
+			right_bounds.grow(ref_ptr[i].bounds);
+			spatial_right_bounds[i - 1] = right_bounds;
+		}
+
+		/* sweep left to right and select lowest SAH. */
+		BoundBox left_bounds;
+
+		for(int i = 1; i < spec.num; i++) {
+			left_bounds.grow(ref_ptr[i - 1].bounds);
+			right_bounds = spatial_right_bounds[i - 1];
+
+			float sah = nodeSAH +
+				left_bounds.area() * params.triangle_cost(i) +
+				right_bounds.area() * params.triangle_cost(spec.num - i);
+
+			if(sah < split.sah) {
+				split.sah = sah;
+				split.dim = dim;
+				split.num_left = i;
+				split.left_bounds = left_bounds;
+				split.right_bounds = right_bounds;
+			}
+		}
+	}
+
+	return split;
+}
+
+void BVHBuild::do_object_split(NodeSpec& left, NodeSpec& right, const NodeSpec& spec, const ObjectSplit& split)
+{
+	/* sort references according to split */
+	int start = references.size() - spec.num;
+	int end = references.size(); /* todo: is this right? */
+
+	bvh_reference_sort(start, end, &references[0], split.dim);
+
+	/* split node specs */
+	left.num = split.num_left;
+	left.bounds = split.left_bounds;
+	right.num = spec.num - split.num_left;
+	right.bounds = split.right_bounds;
+}
+
+/* Spatial Split */
+
+BVHBuild::SpatialSplit BVHBuild::find_spatial_split(const NodeSpec& spec, float nodeSAH)
+{
+	/* initialize bins. */
+	float3 origin = spec.bounds.min;
+	float3 binSize = (spec.bounds.max - origin) * (1.0f / (float)BVHParams::NUM_SPATIAL_BINS);
+	float3 invBinSize = 1.0f / binSize;
+
+	for(int dim = 0; dim < 3; dim++) {
+		for(int i = 0; i < BVHParams::NUM_SPATIAL_BINS; i++) {
+			SpatialBin& bin = spatial_bins[dim][i];
+
+			bin.bounds = BoundBox();
+			bin.enter = 0;
+			bin.exit = 0;
+		}
+	}
+
+	/* chop references into bins. */
+	for(unsigned int refIdx = references.size() - spec.num; refIdx < references.size(); refIdx++) {
+		const Reference& ref = references[refIdx];
+		float3 firstBinf = (ref.bounds.min - origin) * invBinSize;
+		float3 lastBinf = (ref.bounds.max - origin) * invBinSize;
+		int3 firstBin = make_int3(firstBinf.x, firstBinf.y, firstBinf.z);
+		int3 lastBin = make_int3(lastBinf.x, lastBinf.y, lastBinf.z);
+
+		firstBin = clamp(firstBin, 0, BVHParams::NUM_SPATIAL_BINS - 1);
+		lastBin = clamp(lastBin, firstBin, BVHParams::NUM_SPATIAL_BINS - 1);
+
+		for(int dim = 0; dim < 3; dim++) {
+			Reference currRef = ref;
+
+			for(int i = firstBin[dim]; i < lastBin[dim]; i++) {
+				Reference leftRef, rightRef;
+
+				split_reference(leftRef, rightRef, currRef, dim, origin[dim] + binSize[dim] * (float)(i + 1));
+				spatial_bins[dim][i].bounds.grow(leftRef.bounds);
+				currRef = rightRef;
+			}
+
+			spatial_bins[dim][lastBin[dim]].bounds.grow(currRef.bounds);
+			spatial_bins[dim][firstBin[dim]].enter++;
+			spatial_bins[dim][lastBin[dim]].exit++;
+		}
+	}
+
+	/* select best split plane. */
+	SpatialSplit split;
+
+	for(int dim = 0; dim < 3; dim++) {
+		/* sweep right to left and determine bounds. */
+		BoundBox right_bounds;
+
+		for(int i = BVHParams::NUM_SPATIAL_BINS - 1; i > 0; i--) {
+			right_bounds.grow(spatial_bins[dim][i].bounds);
+			spatial_right_bounds[i - 1] = right_bounds;
+		}
+
+		/* sweep left to right and select lowest SAH. */
+		BoundBox left_bounds;
+		int leftNum = 0;
+		int rightNum = spec.num;
+
+		for(int i = 1; i < BVHParams::NUM_SPATIAL_BINS; i++) {
+			left_bounds.grow(spatial_bins[dim][i - 1].bounds);
+			leftNum += spatial_bins[dim][i - 1].enter;
+			rightNum -= spatial_bins[dim][i - 1].exit;
+
+			float sah = nodeSAH +
+				left_bounds.area() * params.triangle_cost(leftNum) +
+				spatial_right_bounds[i - 1].area() * params.triangle_cost(rightNum);
+
+			if(sah < split.sah) {
+				split.sah = sah;
+				split.dim = dim;
+				split.pos = origin[dim] + binSize[dim] * (float)i;
+			}
+		}
+	}
+
+	return split;
+}
+
+void BVHBuild::do_spatial_split(NodeSpec& left, NodeSpec& right, const NodeSpec& spec, const SpatialSplit& split)
+{
+	/* Categorize references and compute bounds.
+	 *
+	 * Left-hand side:			[left_start, left_end[
+	 * Uncategorized/split:		[left_end, right_start[
+	 * Right-hand side:			[right_start, refs.size()[ */
+
+	vector<Reference>& refs = references;
+	int left_start = refs.size() - spec.num;
+	int left_end = left_start;
+	int right_start = refs.size();
+
+	left.bounds = right.bounds = BoundBox();
+
+	for(int i = left_end; i < right_start; i++) {
+		if(refs[i].bounds.max[split.dim] <= split.pos) {
+			/* entirely on the left-hand side */
+			left.bounds.grow(refs[i].bounds);
+			swap(refs[i], refs[left_end++]);
+		}
+		else if(refs[i].bounds.min[split.dim] >= split.pos) {
+			/* entirely on the right-hand side */
+			right.bounds.grow(refs[i].bounds);
+			swap(refs[i--], refs[--right_start]);
+		}
+	}
+
+	/* duplicate or unsplit references intersecting both sides. */
+	while(left_end < right_start) {
+		/* split reference. */
+		Reference lref, rref;
+
+		split_reference(lref, rref, refs[left_end], split.dim, split.pos);
+
+		/* compute SAH for duplicate/unsplit candidates. */
+		BoundBox lub = left.bounds;		// Unsplit to left:		new left-hand bounds.
+		BoundBox rub = right.bounds;	// Unsplit to right:	new right-hand bounds.
+		BoundBox ldb = left.bounds;		// Duplicate:			new left-hand bounds.
+		BoundBox rdb = right.bounds;	// Duplicate:			new right-hand bounds.
+
+		lub.grow(refs[left_end].bounds);
+		rub.grow(refs[left_end].bounds);
+		ldb.grow(lref.bounds);
+		rdb.grow(rref.bounds);
+
+		float lac = params.triangle_cost(left_end - left_start);
+		float rac = params.triangle_cost(refs.size() - right_start);
+		float lbc = params.triangle_cost(left_end - left_start + 1);
+		float rbc = params.triangle_cost(refs.size() - right_start + 1);
+
+		float unsplitLeftSAH = lub.area() * lbc + right.bounds.area() * rac;
+		float unsplitRightSAH = left.bounds.area() * lac + rub.area() * rbc;
+		float duplicateSAH = ldb.area() * lbc + rdb.area() * rbc;
+		float minSAH = min(min(unsplitLeftSAH, unsplitRightSAH), duplicateSAH);
+
+		if(minSAH == unsplitLeftSAH) {
+			/* unsplit to left */
+			left.bounds = lub;
+			left_end++;
+		}
+		else if(minSAH == unsplitRightSAH) {
+			/* unsplit to right */
+			right.bounds = rub;
+			swap(refs[left_end], refs[--right_start]);
+		}
+		else {
+			/* duplicate */
+			left.bounds = ldb;
+			right.bounds = rdb;
+			refs[left_end++] = lref;
+			refs.push_back(rref);
+		}
+	}
+
+	left.num = left_end - left_start;
+	right.num = refs.size() - right_start;
+}
+
+void BVHBuild::split_reference(Reference& left, Reference& right, const Reference& ref, int dim, float pos)
+{
+	/* initialize references. */
+	left.prim_index = right.prim_index = ref.prim_index;
+	left.prim_object = right.prim_object = ref.prim_object;
+	left.bounds = right.bounds = BoundBox();
+
+	/* loop over vertices/edges. */
+	Object *ob = objects[ref.prim_object];
+	const Mesh *mesh = ob->mesh;
+	const int *inds = mesh->triangles[ref.prim_index].v;
+	const float3 *verts = &mesh->verts[0];
+	const float3* v1 = &verts[inds[2]];
+
+	for(int i = 0; i < 3; i++) {
+		const float3* v0 = v1;
+		int vindex = inds[i];
+		v1 = &verts[vindex];
+		float v0p = (*v0)[dim];
+		float v1p = (*v1)[dim];
+
+		/* insert vertex to the boxes it belongs to. */
+		if(v0p <= pos)
+			left.bounds.grow(*v0);
+
+		if(v0p >= pos)
+			right.bounds.grow(*v0);
+
+		/* edge intersects the plane => insert intersection to both boxes. */
+		if((v0p < pos && v1p > pos) || (v0p > pos && v1p < pos)) {
+			float3 t = lerp(*v0, *v1, clamp((pos - v0p) / (v1p - v0p), 0.0f, 1.0f));
+			left.bounds.grow(t);
+			right.bounds.grow(t);
+		}
+	}
+
+	/* intersect with original bounds. */
+	left.bounds.max[dim] = pos;
+	right.bounds.min[dim] = pos;
+	left.bounds.intersect(ref.bounds);
+	right.bounds.intersect(ref.bounds);
+}
+
+CCL_NAMESPACE_END
+
diff --git a/intern/cycles/bvh/bvh_build.h b/intern/cycles/bvh/bvh_build.h
new file mode 100644
index 00000000000..1fa1951d7f2
--- /dev/null
+++ b/intern/cycles/bvh/bvh_build.h
@@ -0,0 +1,152 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef __BVH_BUILD_H__
+#define __BVH_BUILD_H__
+
+#include <float.h>
+
+#include "bvh.h"
+
+#include "util_boundbox.h"
+#include "util_vector.h"
+
+CCL_NAMESPACE_BEGIN
+
+class BVHParams;
+class Mesh;
+class Object;
+class Progress;
+
+/* BVH Builder */
+
+class BVHBuild
+{
+public:
+	struct Reference
+	{
+		int prim_index;
+		int prim_object;
+		BoundBox bounds;
+
+		Reference()
+		{
+		}
+	};
+
+	struct NodeSpec
+	{
+		int num;
+		BoundBox bounds;
+
+		NodeSpec()
+		{
+			num = 0;
+		}
+	};
+
+	BVHBuild(
+		const vector<Object*>& objects,
+		vector<int>& prim_index,
+		vector<int>& prim_object,
+		const BVHParams& params,
+		Progress& progress);
+	~BVHBuild();
+
+	BVHNode *run();
+
+protected:
+	/* adding references */
+	void add_reference_mesh(NodeSpec& root, Mesh *mesh, int i);
+	void add_reference_object(NodeSpec& root, Object *ob, int i);
+	void add_references(NodeSpec& root);
+
+	/* building */
+	BVHNode *build_node(const NodeSpec& spec, int level, float progress_start, float progress_end);
+	BVHNode *create_leaf_node(const NodeSpec& spec);
+	BVHNode *create_object_leaf_nodes(const Reference *ref, int num);
+
+	void progress_update(float progress_start, float progress_end);
+
+	/* object splits */
+	struct ObjectSplit
+	{
+		float sah;
+		int dim;
+		int num_left;
+		BoundBox left_bounds;
+		BoundBox right_bounds;
+
+		ObjectSplit()
+		: sah(FLT_MAX), dim(0), num_left(0)
+		{
+		}
+	};
+
+	ObjectSplit find_object_split(const NodeSpec& spec, float nodeSAH);
+	void do_object_split(NodeSpec& left, NodeSpec& right, const NodeSpec& spec, const ObjectSplit& split);
+
+	/* spatial splits */
+	struct SpatialSplit
+	{
+		float sah;
+		int dim;
+		float pos;
+
+		SpatialSplit()
+		: sah(FLT_MAX), dim(0), pos(0.0f)
+		{
+		}
+	};
+
+	struct SpatialBin
+	{
+		BoundBox bounds;
+		int enter;
+		int exit;
+	};
+
+	SpatialSplit find_spatial_split(const NodeSpec& spec, float nodeSAH);
+	void do_spatial_split(NodeSpec& left, NodeSpec& right, const NodeSpec& spec, const SpatialSplit& split);
+	void split_reference(Reference& left, Reference& right, const Reference& ref, int dim, float pos);
+
+	/* objects and primitive references */
+	vector<Object*> objects;
+	vector<Reference> references;
+
+	/* output primitive indexes and objects */
+	vector<int>& prim_index;
+	vector<int>& prim_object;
+
+	/* build parameters */
+	BVHParams params;
+
+	/* progress reporting */
+	Progress& progress;
+	double progress_start_time;
+	int progress_num_duplicates;
+
+	/* spatial splitting */
+	float spatial_min_overlap;
+	vector<BoundBox> spatial_right_bounds;
+	SpatialBin spatial_bins[3][BVHParams::NUM_SPATIAL_BINS];
+};
+
+CCL_NAMESPACE_END
+
+#endif /* __BVH_BUILD_H__ */
+
diff --git a/intern/cycles/bvh/bvh_node.cpp b/intern/cycles/bvh/bvh_node.cpp
new file mode 100644
index 00000000000..63683bae4a3
--- /dev/null
+++ b/intern/cycles/bvh/bvh_node.cpp
@@ -0,0 +1,101 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "bvh.h"
+#include "bvh_build.h"
+#include "bvh_node.h"
+
+#include "util_debug.h"
+#include "util_vector.h"
+
+CCL_NAMESPACE_BEGIN
+
+int BVHNode::getSubtreeSize(BVH_STAT stat) const
+{
+	int cnt = 0;
+
+	switch(stat)
+	{
+		case BVH_STAT_NODE_COUNT:
+			cnt = 1;
+			break;
+		case BVH_STAT_LEAF_COUNT:
+			cnt = is_leaf() ? 1 : 0;
+			break;
+		case BVH_STAT_INNER_COUNT:
+			cnt = is_leaf() ? 0 : 1;
+			break;
+		case BVH_STAT_TRIANGLE_COUNT:
+			cnt = is_leaf() ? reinterpret_cast<const LeafNode*>(this)->num_triangles() : 0;
+			break;
+		case BVH_STAT_CHILDNODE_COUNT:
+			cnt = num_children();
+			break;
+		default:
+			assert(0); /* unknown mode */
+	}
+
+	if(!is_leaf())
+		for(int i=0;i<num_children();i++)
+			cnt += get_child(i)->getSubtreeSize(stat);
+
+	return cnt;
+}
+
+void BVHNode::deleteSubtree()
+{
+	for(int i=0;i<num_children();i++)
+		get_child(i)->deleteSubtree();
+
+	delete this;
+}
+
+float BVHNode::computeSubtreeSAHCost(const BVHParams& p, float probability) const
+{
+	float SAH = probability * p.cost(num_children(), num_triangles());
+
+	for(int i=0;i<num_children();i++) {
+		BVHNode *child = get_child(i);
+		SAH += child->computeSubtreeSAHCost(p, probability * child->m_bounds.area()/m_bounds.area());
+	}
+
+	return SAH;
+}
+
+void InnerNode::print(int depth) const
+{
+	for(int i = 0; i < depth; i++)
+		printf("  ");
+	
+	printf("inner node %p\n", (void*)this);
+
+	if(children[0])
+		children[0]->print(depth+1);
+	if(children[1])
+		children[1]->print(depth+1);
+}
+
+void LeafNode::print(int depth) const
+{
+	for(int i = 0; i < depth; i++)
+		printf("  ");
+	
+	printf("leaf node %d to %d\n", m_lo, m_hi);
+}
+
+CCL_NAMESPACE_END
+
diff --git a/intern/cycles/bvh/bvh_node.h b/intern/cycles/bvh/bvh_node.h
new file mode 100644
index 00000000000..d83c006b93d
--- /dev/null
+++ b/intern/cycles/bvh/bvh_node.h
@@ -0,0 +1,108 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef __BVH_NODE_H__
+#define __BVH_NODE_H__
+
+#include "util_boundbox.h"
+#include "util_debug.h"
+#include "util_types.h"
+
+CCL_NAMESPACE_BEGIN
+
+enum BVH_STAT
+{
+	BVH_STAT_NODE_COUNT,
+	BVH_STAT_INNER_COUNT,
+	BVH_STAT_LEAF_COUNT,
+	BVH_STAT_TRIANGLE_COUNT,
+	BVH_STAT_CHILDNODE_COUNT
+};
+
+class BVHParams;
+
+class BVHNode
+{
+public:
+	BVHNode()
+	{
+	}
+
+	virtual bool is_leaf() const = 0;
+	virtual int num_children() const = 0;
+	virtual BVHNode *get_child(int i) const = 0;
+	virtual int num_triangles() const { return 0; }
+	virtual void print(int depth = 0) const = 0;
+
+	float getArea() const { return m_bounds.area(); }
+
+	BoundBox m_bounds;
+
+	// Subtree functions
+	int getSubtreeSize(BVH_STAT stat=BVH_STAT_NODE_COUNT) const;
+	float computeSubtreeSAHCost(const BVHParams& p, float probability = 1.0f) const;	
+	void deleteSubtree();
+};
+
+class InnerNode : public BVHNode
+{
+public:
+	InnerNode(const BoundBox& bounds, BVHNode* child0, BVHNode* child1)
+	{
+		m_bounds = bounds;
+		children[0] = child0;
+		children[1] = child1;
+	}
+
+	bool is_leaf() const { return false; }
+	int num_children() const { return 2; }
+	BVHNode *get_child(int i) const{ assert(i>=0 && i<2); return children[i]; }
+	void print(int depth) const;
+
+	BVHNode *children[2];
+};
+
+class LeafNode : public BVHNode
+{
+public:
+	LeafNode(const BoundBox& bounds, int lo, int hi) 
+	{
+		m_bounds = bounds;
+		m_lo = lo;
+		m_hi = hi;
+	}
+
+	LeafNode(const LeafNode& s)
+	: BVHNode()
+	{
+		*this = s;
+	}
+
+	bool is_leaf() const { return true; }
+	int num_children() const { return 0; }
+	BVHNode *get_child(int) const { return NULL; }
+	int num_triangles() const { return m_hi - m_lo; }
+	void print(int depth) const;
+
+	int m_lo;
+	int m_hi;
+};
+
+CCL_NAMESPACE_END
+
+#endif /* __BVH_NODE_H__ */
+
diff --git a/intern/cycles/bvh/bvh_params.h b/intern/cycles/bvh/bvh_params.h
new file mode 100644
index 00000000000..b38e40cfbda
--- /dev/null
+++ b/intern/cycles/bvh/bvh_params.h
@@ -0,0 +1,86 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef __BVH_PARAMS_H__
+#define __BVH_PARAMS_H__
+
+CCL_NAMESPACE_BEGIN
+
+/* BVH Parameters */
+
+class BVHParams
+{
+public:
+	/* spatial split area threshold */
+	bool use_spatial_split;
+	float spatial_split_alpha;
+
+	/* SAH costs */
+	float sah_node_cost;
+	float sah_triangle_cost;
+
+	/* number of triangles in leaf */
+	int min_leaf_size;
+	int max_leaf_size;
+
+	/* object or mesh level bvh */
+	bool top_level;
+
+	/* disk cache */
+	bool use_cache;
+
+	/* QBVH */
+	bool use_qbvh;
+
+	/* fixed parameters */
+	enum {
+		MAX_DEPTH = 64,
+		MAX_SPATIAL_DEPTH = 48,
+		NUM_SPATIAL_BINS = 32
+	};
+
+	BVHParams()
+	{
+		use_spatial_split = true;
+		spatial_split_alpha = 1e-5f;
+
+		sah_node_cost = 1.0f;
+		sah_triangle_cost = 1.0f;
+
+		min_leaf_size = 1;
+		max_leaf_size = 0x7FFFFFF;
+
+		top_level = false;
+		use_cache = false;
+		use_qbvh = false;
+	}
+
+	/* SAH costs */
+	float cost(int num_nodes, int num_tris) const
+	{ return node_cost(num_nodes) + triangle_cost(num_tris); }
+
+	float triangle_cost(int n) const
+	{ return n*sah_triangle_cost; }
+
+	float node_cost(int n) const
+	{ return n*sah_node_cost; }
+};
+
+CCL_NAMESPACE_END
+
+#endif /* __BVH_PARAMS_H__ */
+
diff --git a/intern/cycles/bvh/bvh_sort.cpp b/intern/cycles/bvh/bvh_sort.cpp
new file mode 100644
index 00000000000..ee4531a4843
--- /dev/null
+++ b/intern/cycles/bvh/bvh_sort.cpp
@@ -0,0 +1,57 @@
+/*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+ 
+#include "bvh_build.h"
+#include "bvh_sort.h"
+
+#include "util_algorithm.h"
+#include "util_debug.h"
+
+CCL_NAMESPACE_BEGIN
+
+struct BVHReferenceCompare {
+public:
+	int dim;
+
+	BVHReferenceCompare(int dim_)
+	{
+		dim = dim_;
+	}
+
+	bool operator()(const BVHBuild::Reference& ra, const BVHBuild::Reference& rb)
+	{
+		float ca = ra.bounds.min[dim] + ra.bounds.max[dim];
+		float cb = rb.bounds.min[dim] + rb.bounds.max[dim];
+
+		if(ca < cb) return true;
+		else if(ca > cb) return false;
+		else if(ra.prim_object < rb.prim_object) return true;
+		else if(ra.prim_object > rb.prim_object) return false;
+		else if(ra.prim_index < rb.prim_index) return true;
+		else if(ra.prim_index > rb.prim_index) return false;
+
+		return false;
+	}
+};
+
+void bvh_reference_sort(int start, int end, BVHBuild::Reference *data, int dim)
+{
+	sort(data+start, data+end, BVHReferenceCompare(dim));
+}
+
+CCL_NAMESPACE_END
+
diff --git a/intern/cycles/bvh/bvh_sort.h b/intern/cycles/bvh/bvh_sort.h
new file mode 100644
index 00000000000..f0676948146
--- /dev/null
+++ b/intern/cycles/bvh/bvh_sort.h
@@ -0,0 +1,28 @@
+ /*
+ * Adapted from code copyright 2009-2010 NVIDIA Corporation
+ * Modifications Copyright 2011, Blender Foundation.
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#ifndef __BVH_SORT_H__
+#define __BVH_SORT_H__
+
+CCL_NAMESPACE_BEGIN
+
+void bvh_reference_sort(int start, int end, BVHBuild::Reference *data, int dim);
+
+CCL_NAMESPACE_END
+
+#endif /* __BVH_SORT_H__ */
+