Merge remote-tracking branch 'origin/master' into openvdb

This node was already disabled. All vdb data should eventually be loaded
through the Cycles sync procedure as object data.

Specialized nodes for accessing voxel data may eventually be added again
for convenience, but these would just be advanced attribute nodes.

CMakeLists.txt

@@ -249,8 +249,8 @@ option(WITH_COMPOSITOR         "Enable the tile based nodal compositor" ON)
 
 option(WITH_OPENSUBDIV    "Enable OpenSubdiv for surface subdivision" _init_OPENSUBDIV)
 
-option(WITH_OPENVDB       "Enable features relying on OpenVDB" OFF)
-option(WITH_OPENVDB_BLOSC "Enable blosc compression for OpenVDB, only enable if OpenVDB was built with blosc support" OFF)
+option(WITH_OPENVDB       "Enable features relying on OpenVDB" ON)
+option(WITH_OPENVDB_BLOSC "Enable blosc compression for OpenVDB, only enable if OpenVDB was built with blosc support" ON)
 
 # GHOST Windowing Library Options
 option(WITH_GHOST_DEBUG   "Enable debugging output for the GHOST library" OFF)

intern/cycles/CMakeLists.txt

@@ -163,10 +163,20 @@ if(WITH_CYCLES_OPENSUBDIV)
 	)
 endif()
 
+if(WITH_OPENVDB)
+	add_definitions(-DWITH_OPENVDB)
+	add_definitions(-DDWREAL_IS_DOUBLE=0)
+	add_definitions(-DTBB_USE_EXCEPTIONS=0)
+	include_directories(
+		${OPENVDB_INCLUDE_DIRS}
+	)
+endif()
+
 set(WITH_CYCLES_DEVICE_OPENCL TRUE)
 set(WITH_CYCLES_DEVICE_CUDA TRUE)
 set(WITH_CYCLES_DEVICE_MULTI TRUE)
 
+
 if(CYCLES_STANDALONE_REPOSITORY)
 	TEST_UNORDERED_MAP_SUPPORT()
 endif()

intern/cycles/app/cycles_xml.cpp

@@ -223,7 +223,6 @@ static void xml_read_shader_graph(XMLReadState& state, Shader *shader, pugi::xml
 	for(pugi::xml_node node = graph_node.first_child(); node; node = node.next_sibling()) {
 		ustring node_name(node.name());
 
-		if(node_name == "connect") {
 			/* connect nodes */
 			vector<string> from_tokens, to_tokens;
 

intern/cycles/blender/CMakeLists.txt

@@ -28,6 +28,7 @@ set(SRC
 	blender_shader.cpp
 	blender_sync.cpp
 	blender_texture.cpp
+	blender_volume.cpp
 
 	CCL_api.h
 	blender_object_cull.h

intern/cycles/blender/blender_mesh.cpp

@@ -1088,7 +1088,7 @@ Mesh *BlenderSync::sync_mesh(BL::Object& b_ob,
 				else
 					create_mesh(scene, mesh, b_mesh, used_shaders, false);
 
-				create_mesh_volume_attributes(scene, b_ob, mesh, b_scene.frame_current());
+				//create_mesh_volume_attributes(scene, b_ob, mesh, b_scene.frame_current());
 			}
 
 			if(render_layer.use_hair && mesh->subdivision_type == Mesh::SUBDIVISION_NONE)

intern/cycles/blender/blender_object.cpp

@@ -24,6 +24,7 @@
 #include "render/nodes.h"
 #include "render/particles.h"
 #include "render/shader.h"
+#include "volume.h"
 
 #include "blender/blender_object_cull.h"
 #include "blender/blender_sync.h"
@@ -92,6 +93,21 @@ bool BlenderSync::object_is_light(BL::Object& b_ob)
 	return (b_ob_data && b_ob_data.is_a(&RNA_Lamp));
 }
 
+bool BlenderSync::object_has_sparse_volume(BL::Object& b_ob)
+{
+	BL::SmokeDomainSettings b_domain = object_smoke_domain_find(b_ob);
+	if (b_domain) {
+		BL::PointCache b_ptcache = b_domain.point_cache();
+		if (b_ptcache.is_baked() && b_domain.cache_file_format() == BL::SmokeDomainSettings::cache_file_format_OPENVDB) {
+			char filename[1024];
+			SmokeDomainSettings_cache_filename_get(&b_domain.ptr, filename);
+			return strcmp(filename, "");
+		}
+	}
+	
+	return false;
+}
+
 static uint object_ray_visibility(BL::Object& b_ob)
 {
 	PointerRNA cvisibility = RNA_pointer_get(&b_ob.ptr, "cycles_visibility");
@@ -349,8 +365,13 @@ Object *BlenderSync::sync_object(BL::Object& b_parent,
 	if(object_map.sync(&object, b_ob, b_parent, key))
 		object_updated = true;
 	
-	/* mesh sync */
-	object->mesh = sync_mesh(b_ob, object_updated, hide_tris);
+	if(object_has_sparse_volume(b_ob)) {
+		sync_volume(b_ob);
+	}
+	/*else*/ {
+		/* mesh sync */
+		object->mesh = sync_mesh(b_ob, object_updated, hide_tris);
+	}
 
 	/* special case not tracked by object update flags */
 
@@ -513,6 +534,7 @@ void BlenderSync::sync_objects(float motion_time)
 		mesh_map.pre_sync();
 		object_map.pre_sync();
 		particle_system_map.pre_sync();
+		volume_map.pre_sync();
 		motion_times.clear();
 	}
 	else {
@@ -633,6 +655,8 @@ void BlenderSync::sync_objects(float motion_time)
 			scene->object_manager->tag_update(scene);
 		if(particle_system_map.post_sync())
 			scene->particle_system_manager->tag_update(scene);
+		if(volume_map.post_sync())
+			scene->volume_manager->tag_update(scene);
 	}
 
 	if(motion)

intern/cycles/blender/blender_sync.cpp

@@ -56,6 +56,7 @@ BlenderSync::BlenderSync(BL::RenderEngine& b_engine,
   mesh_map(&scene->meshes),
   light_map(&scene->lights),
   particle_system_map(&scene->particle_systems),
+  volume_map(&scene->volumes),
   world_map(NULL),
   world_recalc(false),
   scene(scene),
@@ -148,6 +149,10 @@ bool BlenderSync::sync_recalc()
 			for(b_ob->particle_systems.begin(b_psys); b_psys != b_ob->particle_systems.end(); ++b_psys)
 				particle_system_map.set_recalc(*b_ob);
 		}
+		
+		if(b_ob->is_updated()) {
+			volume_map.set_recalc(*b_ob);
+		}
 	}
 
 	BL::BlendData::meshes_iterator b_mesh;
@@ -182,6 +187,7 @@ bool BlenderSync::sync_recalc()
 		light_map.has_recalc() ||
 		mesh_map.has_recalc() ||
 		particle_system_map.has_recalc() ||
+		volume_map.has_recalc() ||
 		BlendDataObjects_is_updated_get(&b_data.ptr) ||
 		world_recalc;
 

intern/cycles/blender/blender_sync.h

@@ -46,6 +46,7 @@ class Scene;
 class Shader;
 class ShaderGraph;
 class ShaderNode;
+class Volume;
 
 class BlenderSync {
 public:
@@ -119,6 +120,7 @@ private:
 	                 BL::Object& b_ob,
 	                 bool motion,
 	                 int time_index = 0);
+	Volume *sync_volume(BL::Object& b_ob);
 	Object *sync_object(BL::Object& b_parent,
 	                    int persistent_id[OBJECT_PERSISTENT_ID_SIZE],
 	                    BL::DupliObject& b_dupli_ob,
@@ -155,6 +157,7 @@ private:
 	bool BKE_object_is_modified(BL::Object& b_ob);
 	bool object_is_mesh(BL::Object& b_ob);
 	bool object_is_light(BL::Object& b_ob);
+	bool object_has_sparse_volume(BL::Object& b_ob);
 
 	/* variables */
 	BL::RenderEngine b_engine;
@@ -166,8 +169,10 @@ private:
 	id_map<void*, Mesh> mesh_map;
 	id_map<ObjectKey, Light> light_map;
 	id_map<ParticleSystemKey, ParticleSystem> particle_system_map;
+	id_map<VolumeKey, Volume> volume_map;
 	set<Mesh*> mesh_synced;
 	set<Mesh*> mesh_motion_synced;
+	set<Volume*> volume_synced;
 	set<float> motion_times;
 	void *world_map;
 	bool world_recalc;

intern/cycles/blender/blender_util.h

@@ -813,6 +813,26 @@ protected:
 	set< std::pair<int, int> > edges_;
 };
 
+/* Volume Key */
+
+/* XXX For now we just use Object ID as a volume key;
+ * Volumes may become a true ID block in Blender later,
+ * or the key can be augmented to distinguish multiple volumes inside the same object.
+ */
+struct VolumeKey {
+	void *ob;
+
+	VolumeKey(void *ob_)
+	: ob(ob_)
+	{
+	}
+
+	bool operator<(const VolumeKey& k) const
+	{
+		return ob < k.ob;
+	}
+};
+
 CCL_NAMESPACE_END
 
 #endif /* __BLENDER_UTIL_H__ */

intern/cycles/blender/blender_volume.cpp

@@ -0,0 +1,185 @@
+/*
+ * Copyright 2011-2016 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 "blender_sync.h"
+
+#include "attribute.h"
+#include "../render/volume.h"
+
+#include "util_foreach.h"
+
+#include <openvdb/openvdb.h>
+
+CCL_NAMESPACE_BEGIN
+
+static Attribute *create_openvdb_attribute(Volume *volume,
+                                           const string& filename,
+                                           const ustring& name)
+{
+	Attribute *attr = NULL;
+
+	openvdb::initialize();
+
+	openvdb::io::File file(filename);
+	file.open();
+
+	openvdb::GridBase::ConstPtr grid = file.readGrid(name.string());
+
+	openvdb::Name value_type = grid->valueType();
+
+	if(value_type == "float") {
+		attr = volume->attributes.add(name, TypeDesc::TypeFloat, ATTR_ELEMENT_VOXEL);
+	}
+	else if(value_type == "vec3s") {
+		if (grid->getMetadata< openvdb::TypedMetadata<bool> >("is_color")) {
+			attr = volume->attributes.add(name, TypeDesc::TypeColor, ATTR_ELEMENT_VOXEL);
+		}
+		else {
+			attr = volume->attributes.add(name, TypeDesc::TypeVector, ATTR_ELEMENT_VOXEL);
+		}
+	}
+	else {
+		fprintf(stderr, "Skipping volume attribute: %s\n", name.string().c_str());
+	}
+
+	return attr;
+}
+
+static Attribute *create_smoke_attribute(BL::Object& b_ob,
+                                         Volume *volume,
+                                         const ustring& name,
+                                         float /*frame*/,
+                                         string *filename)
+{
+	BL::SmokeDomainSettings b_domain = object_smoke_domain_find(b_ob);
+	if(b_domain) {
+		char filename_buf[1024];
+		SmokeDomainSettings_cache_filename_get(&b_domain.ptr, filename_buf);
+		*filename = string(filename_buf);
+		
+		return create_openvdb_attribute(volume, *filename, name);
+	}
+	
+	return NULL;
+}
+
+static bool is_volume_attribute(AttributeStandard std) {
+	return std == ATTR_STD_VOLUME_DENSITY
+	        || std == ATTR_STD_VOLUME_COLOR
+	        || std == ATTR_STD_VOLUME_FLAME
+	        || std == ATTR_STD_VOLUME_HEAT
+	        || std == ATTR_STD_VOLUME_VELOCITY;
+}
+
+static void create_volume_attributes(Scene *scene,
+                                     BL::Object& b_ob,
+                                     Volume *volume,
+                                     float frame)
+{
+	foreach(Shader *shader, volume->used_shaders) {
+		foreach(AttributeRequest req, shader->attributes.requests) {
+			ustring name;
+			if (is_volume_attribute(req.std)) {
+				name = ustring(Attribute::standard_name(req.std));
+			}
+			else {
+				name = req.name;
+			}
+			
+			if (!name.empty()) {
+				string filename;
+				Attribute *attr = create_smoke_attribute(b_ob, volume, name, frame, &filename);
+				
+				if (attr) {
+					VoxelAttribute *volume_data = attr->data_voxel();
+					assert(volume_data && "Failed to create volume data!\n");
+					
+					// TODO(kevin): add volume fields to the Volume*
+					//volume_data->manager = volume_manager;
+					volume_data->slot = scene->volume_manager->add_volume(volume, filename, name.string());
+				}
+			}
+		}
+	}
+}
+
+Volume *BlenderSync::sync_volume(BL::Object &b_ob)
+{
+	BL::ID key = b_ob;
+	BL::Material material_override = render_layer.material_override;
+
+	/* find shader indices */
+	vector<Shader*> used_shaders;
+	BL::ID b_ob_data = b_ob.data();
+
+	BL::Object::material_slots_iterator slot;
+	for(b_ob.material_slots.begin(slot); slot != b_ob.material_slots.end(); ++slot) {
+		if(material_override) {
+			find_shader(material_override, used_shaders, scene->default_volume);
+		}
+		else {
+			BL::ID b_material(slot->material());
+			find_shader(b_material, used_shaders, scene->default_volume);
+		}
+	}
+
+	if(used_shaders.size() == 0) {
+		if(material_override)
+			find_shader(material_override, used_shaders, scene->default_volume);
+		else
+			used_shaders.push_back(scene->default_volume);
+	}
+
+	Volume *volume;
+
+	if(!volume_map.sync(&volume, key)) {
+		/* test if shaders changed, these can be object level so mesh
+		 * does not get tagged for recalc */
+		if(volume->used_shaders != used_shaders);
+		else {
+			/* even if not tagged for recalc, we may need to sync anyway
+			 * because the shader needs different volume attributes */
+			bool attribute_recalc = false;
+
+			foreach(Shader *shader, volume->used_shaders)
+				if(shader->need_update_attributes)
+					attribute_recalc = true;
+
+			if(!attribute_recalc)
+				return volume;
+		}
+	}
+
+	/* ensure we only sync instanced meshes once */
+	if(volume_synced.find(volume) != volume_synced.end())
+		return volume;
+	
+	volume_synced.insert(volume);
+	
+	volume->used_shaders = used_shaders;
+	volume->name = ustring(b_ob_data.name().c_str());
+
+	create_volume_attributes(scene, b_ob, volume, b_scene.frame_current());
+
+	/* tag update */
+	bool rebuild = false;
+
+	volume->tag_update(scene, rebuild);
+
+	return volume;
+}
+
+CCL_NAMESPACE_END

intern/cycles/device/CMakeLists.txt

@@ -1,6 +1,7 @@
 
 set(INC
 	..
+	../kernel/openvdb
 	../../glew-mx
 )
 

intern/cycles/device/device.h

@@ -35,6 +35,8 @@ CCL_NAMESPACE_BEGIN
 class Progress;
 class RenderTile;
 
+struct OpenVDBGlobals;
+
 /* Device Types */
 
 enum DeviceType {
@@ -310,6 +312,9 @@ public:
 	/* open shading language, only for CPU device */
 	virtual void *osl_memory() { return NULL; }
 
+	/* OpenVDB data */
+	virtual OpenVDBGlobals *vdb_memory() { return NULL; }
+
 	/* load/compile kernels, must be called before adding tasks */ 
 	virtual bool load_kernels(
 	        const DeviceRequestedFeatures& /*requested_features*/)

intern/cycles/device/device_cpu.cpp

@@ -40,6 +40,12 @@
 #include "kernel/osl/osl_shader.h"
 #include "kernel/osl/osl_globals.h"
 
+#ifdef WITH_OPENVDB
+#include "vdb_globals.h"
+#include "vdb_thread.h"
+#endif
+
+
 #include "render/buffers.h"
 
 #include "util/util_debug.h"
@@ -170,6 +176,10 @@ public:
 	OSLGlobals osl_globals;
 #endif
 
+#ifdef WITH_OPENVDB
+	OpenVDBGlobals vdb_globals;
+#endif
+
 	bool use_split_kernel;
 
 	DeviceRequestedFeatures requested_features;
@@ -413,6 +423,15 @@ public:
 #endif
 	}
 
+	OpenVDBGlobals *vdb_memory()
+	{
+#ifdef WITH_OPENVDB
+		return &vdb_globals;
+#else
+		return NULL;
+#endif
+	}
+
 	void thread_run(DeviceTask *task)
 	{
 		if(task->type == DeviceTask::RENDER) {
@@ -800,6 +819,10 @@ public:
 #ifdef WITH_OSL
 		OSLShader::thread_init(&kg, &kernel_globals, &osl_globals);
 #endif
+#ifdef WITH_OPENVDB
+		vdb_thread_init(&kg, &kernel_globals, &vdb_globals);
+#endif
+
 		for(int sample = 0; sample < task.num_samples; sample++) {
 			for(int x = task.shader_x; x < task.shader_x + task.shader_w; x++)
 				shader_kernel()(&kg,
@@ -820,6 +843,9 @@ public:
 
 #ifdef WITH_OSL
 		OSLShader::thread_free(&kg);
+#endif
+#ifdef WITH_OPENVDB
+		vdb_thread_free(&kg);
 #endif
 	}
 
@@ -871,6 +897,9 @@ protected:
 		kg.decoupled_volume_steps_index = 0;
 #ifdef WITH_OSL
 		OSLShader::thread_init(&kg, &kernel_globals, &osl_globals);
+#endif
+#ifdef WITH_OPENVDB
+		vdb_thread_init(&kg, &kernel_globals, &vdb_globals);
 #endif
 		return kg;
 	}
@@ -893,6 +922,9 @@ protected:
 		}
 #ifdef WITH_OSL
 		OSLShader::thread_free(kg);
+#endif
+#ifdef WITH_OPENVDB
+		vdb_thread_free(kg);
 #endif
 	}
 

intern/cycles/kernel/CMakeLists.txt

@@ -430,6 +430,18 @@ if(WITH_CYCLES_OSL)
 	add_subdirectory(shaders)
 endif()
 
+# OpenVDB module
+
+list(APPEND SRC
+	openvdb/vdb_thread.cpp
+	)
+
+list(APPEND SRC_HEADERS
+	openvdb/vdb_globals.h
+	openvdb/vdb_thread.h
+	openvdb/vdb_intern.h
+	)
+
 # CPU module
 
 include_directories(${INC})

intern/cycles/kernel/bvh/bvh_volume.h

@@ -97,6 +97,24 @@ bool BVH_FUNCTION_FULL_NAME(BVH)(KernelGlobals *kg,
 	gen_idirsplat_swap(pn, shuf_identity, shuf_swap, idir, idirsplat, shufflexyz);
 #endif
 
+#if 1
+	/* try to intersect with VDB volumes */
+	int num_volumes = kernel_data.tables.num_volumes;
+
+	for(int i = 0; i < num_volumes; i++) {
+		float t;
+
+		if(vdb_volume_intersect(kg->vdb_tdata, i, ray, &t)) {
+			isect->type = PRIMITIVE_VOLUME;
+			isect->prim = i;
+			isect->t = t;
+			isect->u = 1.0f;
+			isect->v = 1.0f;
+			return true;
+		}
+	}
+#endif
+
 	/* traversal loop */
 	do {
 		do {

intern/cycles/kernel/bvh/bvh_volume_all.h

@@ -101,6 +101,33 @@ uint BVH_FUNCTION_FULL_NAME(BVH)(KernelGlobals *kg,
 	gen_idirsplat_swap(pn, shuf_identity, shuf_swap, idir, idirsplat, shufflexyz);
 #endif  /* __KERNEL_SSE2__ */
 
+#if 1
+	/* try to intersect with VDB volumes */
+	int num_volumes = kernel_data.tables.num_volumes;
+
+	for(int i = 0; i < num_volumes; i++) {
+		float t;
+
+		if(vdb_volume_intersect(kg->vdb_tdata, i, ray, &t)) {
+			isect_array->type = PRIMITIVE_VOLUME;
+			isect_array->prim = i;
+			isect_array->t = t;
+			isect_array->u = 1.0f;
+			isect_array->v = 1.0f;
+			isect_array++;
+			num_hits++;
+
+			if(num_hits == max_hits) {
+				return num_hits;
+			}
+		}
+	}
+
+	if(num_hits > 0) {
+		return num_hits;
+	}
+#endif
+
 	/* traversal loop */
 	do {
 		do {

intern/cycles/kernel/bvh/qbvh_volume.h

@@ -91,6 +91,24 @@ ccl_device bool BVH_FUNCTION_FULL_NAME(QBVH)(KernelGlobals *kg,
 	                       &near_x, &near_y, &near_z,
 	                       &far_x, &far_y, &far_z);
 
+#if 1
+	/* try to intersect with VDB volumes */
+	int num_volumes = kernel_data.tables.num_volumes;
+
+	for(int i = 0; i < num_volumes; i++) {
+		float t;
+
+		if(vdb_volume_intersect(kg->vdb_tdata, i, ray, &t)) {
+			isect->type = PRIMITIVE_VOLUME;
+			isect->prim = i;
+			isect->t = t;
+			isect->u = 1.0f;
+			isect->v = 1.0f;
+			return true;
+		}
+	}
+#endif
+
 	/* Traversal loop. */
 	do {
 		do {

intern/cycles/kernel/bvh/qbvh_volume_all.h

@@ -95,6 +95,33 @@ ccl_device uint BVH_FUNCTION_FULL_NAME(QBVH)(KernelGlobals *kg,
 	                       &near_x, &near_y, &near_z,
 	                       &far_x, &far_y, &far_z);
 
+#if 1
+	/* try to intersect with VDB volumes */
+	int num_volumes = kernel_data.tables.num_volumes;
+
+	for(int i = 0; i < num_volumes; i++) {
+		float t;
+
+		if(vdb_volume_intersect(kg->vdb_tdata, i, ray, &t)) {
+			isect_array->type = PRIMITIVE_VOLUME;
+			isect_array->prim = i;
+			isect_array->t = t;
+			isect_array->u = 1.0f;
+			isect_array->v = 1.0f;
+			isect_array++;
+			num_hits++;
+
+			if(num_hits == max_hits) {
+				return num_hits;
+			}
+		}
+	}
+
+	if(num_hits > 0) {
+		return num_hits;
+	}
+#endif
+
 	/* Traversal loop. */
 	do {
 		do {

intern/cycles/kernel/geom/geom_primitive.h

@@ -40,7 +40,7 @@ ccl_device_inline float primitive_attribute_float(KernelGlobals *kg,
 	}
 #endif
 #ifdef __VOLUME__
-	else if(sd->object != OBJECT_NONE && desc.element == ATTR_ELEMENT_VOXEL) {
+	else if(sd->type & PRIMITIVE_ALL_VOLUME) {
 		return volume_attribute_float(kg, sd, desc, dx, dy);
 	}
 #endif
@@ -68,7 +68,7 @@ ccl_device_inline float3 primitive_attribute_float3(KernelGlobals *kg,
 	}
 #endif
 #ifdef __VOLUME__
-	else if(sd->object != OBJECT_NONE && desc.element == ATTR_ELEMENT_VOXEL) {
+	else if(sd->type & PRIMITIVE_ALL_VOLUME) {
 		return volume_attribute_float3(kg, sd, desc, dx, dy);
 	}
 #endif

intern/cycles/kernel/geom/geom_volume.h

@@ -14,6 +14,10 @@
  * limitations under the License.
  */
 
+#ifdef WITH_OPENVDB
+#include "../openvdb/vdb_thread.h"
+#endif
+
 /* Volume Primitive
  *
  * Volumes are just regions inside meshes with the mesh surface as boundaries.
@@ -49,26 +53,38 @@ ccl_device_inline float3 volume_normalized_position(KernelGlobals *kg,
 
 ccl_device float volume_attribute_float(KernelGlobals *kg, const ShaderData *sd, const AttributeDescriptor desc, float *dx, float *dy)
 {
+#ifdef __OPENVDB__
+	float3 P = sd->P;
+	/* XXX OpenVDB does not support cubic interpolation - lukas_t */
+	float r = kernel_tex_voxel_float(desc.offset, P.x, P.y, P.z, OPENVDB_SAMPLE_BOX);
+#else
 	float3 P = volume_normalized_position(kg, sd, sd->P);
 	InterpolationType interp = (sd->flag & SD_VOLUME_CUBIC)? INTERPOLATION_CUBIC: INTERPOLATION_NONE;
-	float4 r = kernel_tex_image_interp_3d(kg, desc.offset, P.x, P.y, P.z, interp);
+	float r = average(float4_to_float3(kernel_tex_image_interp_3d_float(kg, desc.offset, P.x, P.y, P.z, interp)));
+#endif
 
 	if(dx) *dx = 0.0f;
 	if(dy) *dy = 0.0f;
 
-	return average(float4_to_float3(r));
+	return r;
 }
 
 ccl_device float3 volume_attribute_float3(KernelGlobals *kg, const ShaderData *sd, const AttributeDescriptor desc, float3 *dx, float3 *dy)
 {
+#ifdef __OPENVDB__
+	float3 P = sd->P;
+	/* XXX OpenVDB does not support cubic interpolation - lukas_t */
+	float3 r = kernel_tex_voxel_float3(desc.offset, P.x, P.y, P.z, OPENVDB_SAMPLE_BOX);
+#else
 	float3 P = volume_normalized_position(kg, sd, sd->P);
 	InterpolationType interp = (sd->flag & SD_VOLUME_CUBIC)? INTERPOLATION_CUBIC: INTERPOLATION_NONE;
-	float4 r = kernel_tex_image_interp_3d(kg, desc.offset, P.x, P.y, P.z, interp);
+	float3 r = float4_to_float3(kernel_tex_image_interp_3d(kg, desc.offset, P.x, P.y, P.z, interp));
+#endif
 
 	if(dx) *dx = make_float3(0.0f, 0.0f, 0.0f);
 	if(dy) *dy = make_float3(0.0f, 0.0f, 0.0f);
 
-	return float4_to_float3(r);
+	return r;
 }
 
 #endif

intern/cycles/kernel/kernel_compat_cpu.h

@@ -119,6 +119,8 @@ template<typename T> struct texture  {
 #define kernel_tex_fetch_ssef(tex, index) (kg->tex.fetch_ssef(index))
 #define kernel_tex_fetch_ssei(tex, index) (kg->tex.fetch_ssei(index))
 #define kernel_tex_lookup(tex, t, offset, size) (kg->tex.lookup(t, offset, size))
+#define kernel_tex_voxel_float(tex, x, y, z, sampling) (vdb_volume_sample_scalar(kg->vdb, kg->vdb_tdata, tex, x, y, z, sampling))
+#define kernel_tex_voxel_float3(tex, x, y, z, sampling) (vdb_volume_sample_vector(kg->vdb, kg->vdb_tdata, tex, x, y, z, sampling))
 
 #define kernel_data (kg->__data)
 

intern/cycles/kernel/kernel_globals.h

@@ -42,9 +42,16 @@ struct OSLThreadData;
 struct OSLShadingSystem;
 #  endif
 
+#  ifdef WITH_OPENVDB
+struct OpenVDBGlobals;
+struct OpenVDBThreadData;
+#  endif
+
 struct Intersection;
 struct VolumeStep;
 
+#  define MAX_VOLUME        1024
+
 typedef struct KernelGlobals {
 #  define KERNEL_TEX(type, name) texture<type> name;
 #  define KERNEL_IMAGE_TEX(type, ttype, name)
@@ -75,6 +82,12 @@ typedef struct KernelGlobals {
 
 	int2 global_size;
 	int2 global_id;
+
+#  ifdef WITH_OPENVDB
+	/* OpenVDB */
+	OpenVDBGlobals *vdb;
+	OpenVDBThreadData *vdb_tdata;
+#  endif
 } KernelGlobals;
 
 #endif  /* __KERNEL_CPU__ */

intern/cycles/kernel/kernel_shader.h

@@ -69,7 +69,12 @@ ccl_device_noinline void shader_setup_from_ray(KernelGlobals *kg,
 #endif
 	sd->time = ray->time;
 
-	sd->prim = kernel_tex_fetch(__prim_index, isect->prim);
+	if (sd->type & PRIMITIVE_VOLUME) {
+		sd->prim = isect->prim;
+	}
+	else {
+		sd->prim = kernel_tex_fetch(__prim_index, isect->prim);
+	}
 	sd->ray_length = isect->t;
 
 #ifdef __UV__
@@ -106,6 +111,9 @@ ccl_device_noinline void shader_setup_from_ray(KernelGlobals *kg,
 		triangle_dPdudv(kg, sd->prim, &sd->dPdu, &sd->dPdv);
 #endif
 	}
+	else if(sd->type & PRIMITIVE_VOLUME) {
+		sd->shader = kernel_tex_fetch(__vol_shader, sd->prim);
+	}
 	else {
 		/* motion triangle */
 		motion_triangle_shader_setup(kg, sd, isect, ray, false);
@@ -432,7 +440,7 @@ ccl_device_inline void shader_setup_from_volume(KernelGlobals *kg, ShaderData *s
 	sd->object = PRIM_NONE; /* todo: fill this for texture coordinates */
 #endif
 	sd->prim = PRIM_NONE;
-	sd->type = PRIMITIVE_NONE;
+	sd->type = PRIMITIVE_VOLUME;
 
 #ifdef __UV__
 	sd->u = 0.0f;

intern/cycles/kernel/kernel_textures.h

@@ -78,6 +78,9 @@ KERNEL_TEX(float, __lookup_table)
 /* sobol */
 KERNEL_TEX(uint, __sobol_directions)
 
+/* volume */
+KERNEL_TEX(uint, __vol_shader)
+
 #if !defined(__KERNEL_CUDA__) || __CUDA_ARCH__ >= 300
 /* image textures */
 KERNEL_TEX(TextureInfo, __texture_info)

intern/cycles/kernel/kernel_types.h

@@ -87,6 +87,9 @@ CCL_NAMESPACE_BEGIN
 #  ifdef WITH_OSL
 #    define __OSL__
 #  endif
+#  ifdef WITH_OPENVDB
+#    define __OPENVDB__
+#  endif
 #  define __PRINCIPLED__
 #  define __SUBSURFACE__
 #  define __CMJ__
@@ -686,11 +689,13 @@ typedef enum PrimitiveType {
 	 * since it is no real traceable primitive.
 	 */
 	PRIMITIVE_LAMP            = (1 << 4),
+	PRIMITIVE_VOLUME          = (1 << 5),
 
 	PRIMITIVE_ALL_TRIANGLE = (PRIMITIVE_TRIANGLE|PRIMITIVE_MOTION_TRIANGLE),
 	PRIMITIVE_ALL_CURVE = (PRIMITIVE_CURVE|PRIMITIVE_MOTION_CURVE),
 	PRIMITIVE_ALL_MOTION = (PRIMITIVE_MOTION_TRIANGLE|PRIMITIVE_MOTION_CURVE),
-	PRIMITIVE_ALL = (PRIMITIVE_ALL_TRIANGLE|PRIMITIVE_ALL_CURVE),
+	PRIMITIVE_ALL_VOLUME = (PRIMITIVE_VOLUME),
+	PRIMITIVE_ALL = (PRIMITIVE_ALL_TRIANGLE|PRIMITIVE_ALL_CURVE|PRIMITIVE_ALL_VOLUME),
 
 	/* Total number of different traceable primitives.
 	 * NOTE: This is an actual value, not a bitflag.
@@ -996,6 +1001,7 @@ typedef ccl_addr_space struct ShaderData {
 typedef struct VolumeStack {
 	int object;
 	int shader;
+	int volume;
 } VolumeStack;
 #endif
 
@@ -1326,7 +1332,8 @@ static_assert_align(KernelCurves, 16);
 
 typedef struct KernelTables {
 	int beckmann_offset;
-	int pad1, pad2, pad3;
+	int num_volumes;
+	int density_index, pad2;
 } KernelTables;
 static_assert_align(KernelTables, 16);
 

intern/cycles/kernel/kernel_volume.h

@@ -14,6 +14,8 @@
  * limitations under the License.
  */
 
+#include "openvdb/vdb_thread.h"
+
 CCL_NAMESPACE_BEGIN
 
 /* Events for probalistic scattering */
@@ -21,7 +23,8 @@ CCL_NAMESPACE_BEGIN
 typedef enum VolumeIntegrateResult {
 	VOLUME_PATH_SCATTERED = 0,
 	VOLUME_PATH_ATTENUATED = 1,
-	VOLUME_PATH_MISSED = 2
+	VOLUME_PATH_MISSED = 2,
+	VOLUME_PATH_CONTINUE = 3,
 } VolumeIntegrateResult;
 
 /* Volume shader properties
@@ -173,6 +176,38 @@ ccl_device void kernel_volume_shadow_homogeneous(KernelGlobals *kg,
 		*throughput *= volume_color_transmittance(sigma_t, ray->t);
 }
 
+ccl_device_inline bool kernel_volume_integrate_shadow_ray(
+        KernelGlobals *kg, PathState *state, Ray *ray, ShaderData *sd,
+        float3 *tp, float t, float new_t, float random_jitter_offset,
+        float3 *sum, float tp_eps, int i)
+{
+	float dt = new_t - t;
+
+	/* use random position inside this segment to sample shader */
+	if(new_t == ray->t)
+		random_jitter_offset = lcg_step_float(&state->rng_congruential) * dt;
+
+	float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
+	float3 sigma_t;
+
+	/* compute attenuation over segment */
+	if(volume_shader_extinction_sample(kg, sd, state, new_P, &sigma_t)) {
+		/* Compute expf() only for every Nth step, to save some calculations
+		 * because exp(a)*exp(b) = exp(a+b), also do a quick tp_eps check then. */
+
+		*sum += (-sigma_t * (new_t - t));
+		if((i & 0x07) == 0) { /* ToDo: Other interval? */
+			*tp = *tp * make_float3(expf(sum->x), expf(sum->y), expf(sum->z));
+
+			/* stop if nearly all light is blocked */
+			if(tp->x < tp_eps && tp->y < tp_eps && tp->z < tp_eps)
+				return true;
+		}
+	}
+
+	return false;
+}
+
 /* heterogeneous volume: integrate stepping through the volume until we
  * reach the end, get absorbed entirely, or run out of iterations */
 ccl_device void kernel_volume_shadow_heterogeneous(KernelGlobals *kg,
@@ -194,42 +229,75 @@ ccl_device void kernel_volume_shadow_heterogeneous(KernelGlobals *kg,
 
 	float3 sum = make_float3(0.0f, 0.0f, 0.0f);
 
-	for(int i = 0; i < max_steps; i++) {
-		/* advance to new position */
-		float new_t = min(ray->t, (i+1) * step);
-		float dt = new_t - t;
+#ifdef __OPENVDB__
+//	int density_index = kernel_data.tables.density_index;
+	int num_volumes = kernel_data.tables.num_volumes;
+	bool has_vdb_volume = num_volumes > 0;
+	float t1 = 0.0f;
+	int v = 0;
 
-		/* use random position inside this segment to sample shader */
-		if(new_t == ray->t)
-			random_jitter_offset = lcg_step_float_addrspace(&state->rng_congruential) * dt;
-
-		float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
-		float3 sigma_t;
-
-		/* compute attenuation over segment */
-		if(volume_shader_extinction_sample(kg, sd, state, new_P, &sigma_t)) {
-			/* Compute expf() only for every Nth step, to save some calculations
-			 * because exp(a)*exp(b) = exp(a+b), also do a quick tp_eps check then. */
-
-			sum += (-sigma_t * (new_t - t));
-			if((i & 0x07) == 0) { /* ToDo: Other interval? */
-				tp = *throughput * make_float3(expf(sum.x), expf(sum.y), expf(sum.z));
-
-				/* stop if nearly all light is blocked */
-				if(tp.x < tp_eps && tp.y < tp_eps && tp.z < tp_eps)
-					break;
-			}
-		}
-
-		/* stop if at the end of the volume */
-		t = new_t;
-		if(t == ray->t) {
-			/* Update throughput in case we haven't done it above */
-			tp = *throughput * make_float3(expf(sum.x), expf(sum.y), expf(sum.z));
+	float isec_t = 0.0f;
+	for(; v < num_volumes; v++) {
+		if(vdb_volume_intersect(kg->vdb_tdata, v, ray, &isec_t)) {
 			break;
 		}
 	}
 
+	if(has_vdb_volume && v < num_volumes && vdb_volume_scalar_has_uniform_voxels(kg->vdb, v)) {
+		/* TODO(kevin): this call should be moved out of here, all it does is
+		 * checking if we have an intersection with the boundbox of the volumue
+		 * which in most cases corresponds to the boundbox of the object that has
+		 * this volume. Also it initializes the rays for the ray marching. */
+		//if(!vdb_volume_intersect(kg->vdb_tdata, density_index, ray, NULL)) {
+		//	return;
+		//}
+
+		/* t and t1 represent the entry and exit points for each leaf node or tile
+		 * containing active voxels. If we don't have any active node in the current
+		 * ray path (i.e. empty space) the ray march loop is not executed,
+		 * otherwise we loop through all leaves until the end of the volume. */
+		while(vdb_volume_march(kg->vdb_tdata, v, &t, &t1)) {
+			int i = 0;
+
+			/* Perform small steps through the current leaf or tile. */
+			for(float new_t = step * ceilf(t / step); new_t <= t1; new_t += step) {
+				bool ok = kernel_volume_integrate_shadow_ray(
+				              kg, state, ray, sd, &tp, t, new_t, random_jitter_offset,
+				              &sum, tp_eps, i);
+
+				if (ok) {
+					*throughput = tp;
+					return;
+				}
+
+				/* stop if at the end of the volume */
+				t = new_t;
+				i++;
+			}
+		}
+	}
+	else
+#endif
+	{
+		for(int i = 0; i < max_steps; i++) {
+			/* advance to new position */
+			float new_t = min(ray->t, (i+1) * step);
+
+			bool ok = kernel_volume_integrate_shadow_ray(
+			              kg, state, ray, sd, &tp, t, new_t, random_jitter_offset,
+			              &sum, tp_eps, i);
+
+			/* stop if at the end of the volume */
+			t = new_t;
+			if(ok || t == ray->t) {
+				break;
+			}
+		}
+	}
+
+	/* Update throughput in case we haven't done it above */
+	tp = *throughput * make_float3(expf(sum.x), expf(sum.y), expf(sum.z));
+
 	*throughput = tp;
 }
 
@@ -456,19 +524,130 @@ ccl_device VolumeIntegrateResult kernel_volume_integrate_homogeneous(
 	return VOLUME_PATH_ATTENUATED;
 }
 
+ccl_device_inline VolumeIntegrateResult kernel_volume_integrate_ray(
+        KernelGlobals *kg,
+        PathState *state,
+        Ray *ray,
+        ShaderData *sd,
+        PathRadiance *L,
+        float3 *throughput,
+        float t,
+        float new_t,
+        float random_jitter_offset,
+        bool has_scatter,
+        float3 *accum_transmittance,
+        int channel,
+        const float tp_eps,
+        float *xi)
+{
+	float dt = new_t - t;
+	float3 tp = *throughput;
+
+	/* use random position inside this segment to sample shader */
+	if(new_t == ray->t)
+		random_jitter_offset = lcg_step_float_addrspace(&state->rng_congruential) * dt;
+
+	float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
+	VolumeShaderCoefficients coeff;
+
+	/* compute segment */
+	if(volume_shader_sample(kg, sd, state, new_P, &coeff)) {
+		int closure_flag = sd->flag;
+		float3 new_tp;
+		float3 transmittance = make_float3(1.0f, 1.0f, 1.0f);
+		bool scatter = false;
+
+		/* distance sampling */
+#ifdef __VOLUME_SCATTER__
+		if((closure_flag & SD_SCATTER) || (has_scatter && (closure_flag & SD_ABSORPTION))) {
+			has_scatter = true;
+
+			float3 sigma_t = coeff.sigma_a + coeff.sigma_s;
+			float3 sigma_s = coeff.sigma_s;
+
+			/* compute transmittance over full step */
+			transmittance = volume_color_transmittance(sigma_t, dt);
+
+			/* decide if we will scatter or continue */
+			float sample_transmittance = kernel_volume_channel_get(transmittance, channel);
+
+			if(1.0f - *xi >= sample_transmittance) {
+				/* compute sampling distance */
+				float sample_sigma_t = kernel_volume_channel_get(sigma_t, channel);
+				float new_dt = -logf(1.0f - *xi)/sample_sigma_t;
+				new_t = t + new_dt;
+
+				/* transmittance and pdf */
+				float3 new_transmittance = volume_color_transmittance(sigma_t, new_dt);
+				float3 pdf = sigma_t * new_transmittance;
+
+				/* throughput */
+				new_tp = tp * sigma_s * new_transmittance / average(pdf);
+				scatter = true;
+			}
+			else {
+				/* throughput */
+				float pdf = average(transmittance);
+				new_tp = tp * transmittance / pdf;
+
+				/* remap xi so we can reuse it and keep thing stratified */
+				*xi = 1.0f - (1.0f - *xi)/sample_transmittance;
+			}
+		}
+		else
+#endif
+			if(closure_flag & SD_ABSORPTION) {
+				/* absorption only, no sampling needed */
+				float3 sigma_a = coeff.sigma_a;
+
+				transmittance = volume_color_transmittance(sigma_a, dt);
+				new_tp = tp * transmittance;
+			}
+
+		/* integrate emission attenuated by absorption */
+		if(L && (closure_flag & SD_EMISSION)) {
+			float3 emission = kernel_volume_emission_integrate(&coeff, closure_flag, transmittance, dt);
+			path_radiance_accum_emission(L, state, tp, emission);
+		}
+
+		/* modify throughput */
+		if(closure_flag & (SD_ABSORPTION|SD_SCATTER)) {
+			tp = new_tp;
+
+			/* stop if nearly all light blocked */
+			if(tp.x < tp_eps && tp.y < tp_eps && tp.z < tp_eps) {
+				tp = make_float3(0.0f, 0.0f, 0.0f);
+				*throughput = tp;
+				return VOLUME_PATH_ATTENUATED;
+			}
+		}
+
+		/* prepare to scatter to new direction */
+		if(scatter) {
+			/* adjust throughput and move to new location */
+			sd->P = ray->P + new_t*ray->D;
+			*throughput = tp;
+
+			return VOLUME_PATH_SCATTERED;
+		}
+		else {
+			/* accumulate transmittance */
+			*accum_transmittance *= transmittance;
+		}
+	}
+
+	*throughput = tp;
+	return VOLUME_PATH_CONTINUE;
+}
+
 /* heterogeneous volume distance sampling: integrate stepping through the
  * volume until we reach the end, get absorbed entirely, or run out of
  * iterations. this does probabilistically scatter or get transmitted through
  * for path tracing where we don't want to branch. */
-ccl_device VolumeIntegrateResult kernel_volume_integrate_heterogeneous_distance(
-    KernelGlobals *kg,
-    ccl_addr_space PathState *state,
-    Ray *ray,
-    ShaderData *sd,
-    PathRadiance *L,
-    ccl_addr_space float3 *throughput)
+ccl_device VolumeIntegrateResult kernel_volume_integrate_heterogeneous_distance(KernelGlobals *kg,
+	PathState *state, Ray *ray, ShaderData *sd, PathRadiance *L, float3 *throughput)
 {
-	float3 tp = *throughput;
+	VolumeIntegrateResult result = VOLUME_PATH_MISSED;
 	const float tp_eps = 1e-6f; /* todo: this is likely not the right value */
 
 	/* prepare for stepping */
@@ -485,114 +664,77 @@ ccl_device VolumeIntegrateResult kernel_volume_integrate_heterogeneous_distance(
 	float xi = path_state_rng_1D(kg, state, PRNG_SCATTER_DISTANCE);
 	float rphase = path_state_rng_1D(kg, state, PRNG_PHASE_CHANNEL);
 	int channel = (int)(rphase*3.0f);
+	sd->randb_closure = rphase*3.0f - channel;
 	bool has_scatter = false;
+	bool path_missed = true;
 
-	for(int i = 0; i < max_steps; i++) {
-		/* advance to new position */
-		float new_t = min(ray->t, (i+1) * step_size);
-		float dt = new_t - t;
+#ifdef __OPENVDB__
+//	int density_index = kernel_data.tables.density_index;
+	int num_volumes = kernel_data.tables.num_volumes;
+	bool has_vdb_volume = num_volumes > 0;
+	float t1 = 0.0f;
+	int i;
 
-		/* use random position inside this segment to sample shader */
-		if(new_t == ray->t)
-			random_jitter_offset = lcg_step_float_addrspace(&state->rng_congruential) * dt;
-
-		float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
-		VolumeShaderCoefficients coeff;
-
-		/* compute segment */
-		if(volume_shader_sample(kg, sd, state, new_P, &coeff)) {
-			int closure_flag = sd->flag;
-			float3 new_tp;
-			float3 transmittance;
-			bool scatter = false;
-
-			/* distance sampling */
-#ifdef __VOLUME_SCATTER__
-			if((closure_flag & SD_SCATTER) || (has_scatter && (closure_flag & SD_ABSORPTION))) {
-				has_scatter = true;
-
-				float3 sigma_t = coeff.sigma_a + coeff.sigma_s;
-				float3 sigma_s = coeff.sigma_s;
-
-				/* compute transmittance over full step */
-				transmittance = volume_color_transmittance(sigma_t, dt);
-
-				/* decide if we will scatter or continue */
-				float sample_transmittance = kernel_volume_channel_get(transmittance, channel);
-
-				if(1.0f - xi >= sample_transmittance) {
-					/* compute sampling distance */
-					float sample_sigma_t = kernel_volume_channel_get(sigma_t, channel);
-					float new_dt = -logf(1.0f - xi)/sample_sigma_t;
-					new_t = t + new_dt;
-
-					/* transmittance and pdf */
-					float3 new_transmittance = volume_color_transmittance(sigma_t, new_dt);
-					float3 pdf = sigma_t * new_transmittance;
-
-					/* throughput */
-					new_tp = tp * sigma_s * new_transmittance / average(pdf);
-					scatter = true;
-				}
-				else {
-					/* throughput */
-					float pdf = average(transmittance);
-					new_tp = tp * transmittance / pdf;
-
-					/* remap xi so we can reuse it and keep thing stratified */
-					xi = 1.0f - (1.0f - xi)/sample_transmittance;
-				}
-			}
-			else 
-#endif
-			if(closure_flag & SD_ABSORPTION) {
-				/* absorption only, no sampling needed */
-				float3 sigma_a = coeff.sigma_a;
-
-				transmittance = volume_color_transmittance(sigma_a, dt);
-				new_tp = tp * transmittance;
-			}
-
-			/* integrate emission attenuated by absorption */
-			if(L && (closure_flag & SD_EMISSION)) {
-				float3 emission = kernel_volume_emission_integrate(&coeff, closure_flag, transmittance, dt);
-				path_radiance_accum_emission(L, state, tp, emission);
-			}
-
-			/* modify throughput */
-			if(closure_flag & (SD_ABSORPTION|SD_SCATTER)) {
-				tp = new_tp;
-
-				/* stop if nearly all light blocked */
-				if(tp.x < tp_eps && tp.y < tp_eps && tp.z < tp_eps) {
-					tp = make_float3(0.0f, 0.0f, 0.0f);
-					break;
-				}
-			}
-
-			/* prepare to scatter to new direction */
-			if(scatter) {
-				/* adjust throughput and move to new location */
-				sd->P = ray->P + new_t*ray->D;
-				*throughput = tp;
-
-				return VOLUME_PATH_SCATTERED;
-			}
-			else {
-				/* accumulate transmittance */
-				accum_transmittance *= transmittance;
-			}
-		}
-
-		/* stop if at the end of the volume */
-		t = new_t;
-		if(t == ray->t)
+	for(i = 0; i < num_volumes; i++) {
+		float isec_t = 0.0f;
+		if(vdb_volume_intersect(kg->vdb_tdata, i, ray, &isec_t)) {
 			break;
+		}
 	}
 
-	*throughput = tp;
+	if(has_vdb_volume /*&& i >= 0*/ && vdb_volume_scalar_has_uniform_voxels(kg->vdb, i)) {
+		/* TODO(kevin): this call should be moved out of here, all it does is
+		 * checking if we have an intersection with the boundbox of the volumue
+		 * which in most cases corresponds to the boundbox of the object that has
+		 * this volume. Also it initializes the rays for the ray marching. */
+		//if(!vdb_volume_intersect(kg->vdb_tdata, density_index, ray, NULL)) {
+		//	return VOLUME_PATH_MISSED;
+		//}
 
-	return VOLUME_PATH_ATTENUATED;
+		/* t and t1 represent the entry and exit points for each leaf node or tile
+		 * containing active voxels. If we don't have any active node in the current
+		 * ray path (i.e. empty space) the ray march loop is not executed,
+		 * otherwise we loop through all leaves until the end of the volume. */
+		while(vdb_volume_march(kg->vdb_tdata, i, &t, &t1)) {
+			path_missed = false;
+
+			/* Perform small steps through the current leaf or tile. */
+			for(float new_t = step_size * ceilf(t / step_size); new_t <= t1; new_t += step_size) {
+				result = kernel_volume_integrate_ray(kg, state, ray, sd, L, throughput, t, new_t,
+				                                     random_jitter_offset, has_scatter,
+				                                     &accum_transmittance, channel, tp_eps, &xi);
+
+				if(result != VOLUME_PATH_CONTINUE)
+					return result;
+
+				t = new_t;
+			}
+		}
+	}
+	else
+#endif
+	{
+		path_missed = false;
+
+		for(int i = 0; i < max_steps; i++) {
+			/* advance to new position */
+			float new_t = min(ray->t, (i+1) * step_size);
+
+			result = kernel_volume_integrate_ray(kg, state, ray, sd, L, throughput, t, new_t,
+			                                     random_jitter_offset, has_scatter,
+			                                     &accum_transmittance, channel, tp_eps, &xi);
+
+			if(result != VOLUME_PATH_CONTINUE)
+				return result;
+
+			/* stop if at the end of the volume */
+			t = new_t;
+			if(t == ray->t)
+				break;
+		}
+	}
+
+	return (path_missed) ? VOLUME_PATH_MISSED : VOLUME_PATH_ATTENUATED;
 }
 
 /* get the volume attenuation and emission over line segment defined by
@@ -719,78 +861,175 @@ ccl_device void kernel_volume_decoupled_record(KernelGlobals *kg, PathState *sta
 	bool is_last_step_empty = false;
 
 	VolumeStep *step = segment->steps;
+#ifdef __OPENVDB__
+	int density_index = kernel_data.tables.density_index;
+	bool has_vdb_volume = kernel_data.tables.num_volumes > 0;
+	float t1 = 0.0f;
 
-	for(int i = 0; i < max_steps; i++, step++) {
-		/* advance to new position */
-		float new_t = min(ray->t, (i+1) * step_size);
-		float dt = new_t - t;
-
-		/* use random position inside this segment to sample shader */
-		if(heterogeneous && new_t == ray->t)
-			random_jitter_offset = lcg_step_float(&state->rng_congruential) * dt;
-
-		float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
-		VolumeShaderCoefficients coeff;
-
-		/* compute segment */
-		if(volume_shader_sample(kg, sd, state, new_P, &coeff)) {
-			int closure_flag = sd->flag;
-			float3 sigma_t = coeff.sigma_a + coeff.sigma_s;
-
-			/* compute accumulated transmittance */
-			float3 transmittance = volume_color_transmittance(sigma_t, dt);
-
-			/* compute emission attenuated by absorption */
-			if(closure_flag & SD_EMISSION) {
-				float3 emission = kernel_volume_emission_integrate(&coeff, closure_flag, transmittance, dt);
-				accum_emission += accum_transmittance * emission;
-			}
-
-			accum_transmittance *= transmittance;
-
-			/* compute pdf for distance sampling */
-			float3 pdf_distance = dt * accum_transmittance * coeff.sigma_s;
-			cdf_distance = cdf_distance + pdf_distance;
-
-			/* write step data */
-			step->sigma_t = sigma_t;
-			step->sigma_s = coeff.sigma_s;
-			step->closure_flag = closure_flag;
-
-			segment->closure_flag |= closure_flag;
-
-			is_last_step_empty = false;
-			segment->numsteps++;
+	if(has_vdb_volume && vdb_volume_scalar_has_uniform_voxels(kg->vdb, density_index)) {
+		/* TODO(kevin): this call should be moved out of here, all it does is
+		 * checking if we have an intersection with the boundbox of the volumue
+		 * which in most cases corresponds to the boundbox of the object that has
+		 * this volume. Also it initializes the rays for the ray marching. */
+		float isect_t = 0.0f;
+		if(!vdb_volume_intersect(kg->vdb_tdata, density_index, ray, &isect_t)) {
+			return;
 		}
-		else {
-			if(is_last_step_empty) {
-				/* consecutive empty step, merge */
-				step--;
+
+		/* t and t1 represent the entry and exit points for each leaf node or tile
+		 * containing active voxels. If we don't have any active node in the current
+		 * ray path (i.e. empty space) the ray march loop is not executed,
+		 * otherwise we loop through all leaves until the end of the volume. */
+		while(vdb_volume_march(kg->vdb_tdata, density_index, &t, &t1)) {
+
+			/* Perform small steps through the current leaf or tile. */
+			for(float new_t = step_size * ceilf(t / step_size); new_t <= t1; new_t += step_size, step++) {
+				float dt = new_t - t;
+
+				/* use random position inside this segment to sample shader */
+				if(heterogeneous && new_t == ray->t)
+					random_jitter_offset = lcg_step_float(&state->rng_congruential) * dt;
+
+				float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
+				VolumeShaderCoefficients coeff;
+
+				/* compute segment */
+				if(volume_shader_sample(kg, sd, state, new_P, &coeff)) {
+					int closure_flag = sd->flag;
+					float3 sigma_t = coeff.sigma_a + coeff.sigma_s;
+
+					/* compute accumulated transmittance */
+					float3 transmittance = volume_color_transmittance(sigma_t, dt);
+
+					/* compute emission attenuated by absorption */
+					if(closure_flag & SD_EMISSION) {
+						float3 emission = kernel_volume_emission_integrate(&coeff, closure_flag, transmittance, dt);
+						accum_emission += accum_transmittance * emission;
+					}
+
+					accum_transmittance *= transmittance;
+
+					/* compute pdf for distance sampling */
+					float3 pdf_distance = dt * accum_transmittance * coeff.sigma_s;
+					cdf_distance = cdf_distance + pdf_distance;
+
+					/* write step data */
+					step->sigma_t = sigma_t;
+					step->sigma_s = coeff.sigma_s;
+					step->closure_flag = closure_flag;
+
+					segment->closure_flag |= closure_flag;
+
+					is_last_step_empty = false;
+					segment->numsteps++;
+				}
+				else {
+					if(is_last_step_empty) {
+						/* consecutive empty step, merge */
+						step--;
+					}
+					else {
+						/* store empty step */
+						step->sigma_t = make_float3(0.0f, 0.0f, 0.0f);
+						step->sigma_s = make_float3(0.0f, 0.0f, 0.0f);
+						step->closure_flag = 0;
+
+						segment->numsteps++;
+						is_last_step_empty = true;
+					}
+				}
+
+				step->accum_transmittance = accum_transmittance;
+				step->cdf_distance = cdf_distance;
+				step->t = new_t;
+				step->shade_t = t + random_jitter_offset;
+
+				/* stop if at the end of the volume */
+				t = new_t;
+				if(t == ray->t)
+					break;
+
+				/* stop if nearly all light blocked */
+				if(accum_transmittance.x < tp_eps && accum_transmittance.y < tp_eps && accum_transmittance.z < tp_eps)
+					break;
+			}
+		}
+	}
+	else
+#endif
+	{
+		for(int i = 0; i < max_steps; i++, step++) {
+			/* advance to new position */
+			float new_t = min(ray->t, (i+1) * step_size);
+			float dt = new_t - t;
+
+			/* use random position inside this segment to sample shader */
+			if(heterogeneous && new_t == ray->t)
+				random_jitter_offset = lcg_step_float(&state->rng_congruential) * dt;
+
+			float3 new_P = ray->P + ray->D * (t + random_jitter_offset);
+			VolumeShaderCoefficients coeff;
+
+			/* compute segment */
+			if(volume_shader_sample(kg, sd, state, new_P, &coeff)) {
+				int closure_flag = sd->flag;
+				float3 sigma_t = coeff.sigma_a + coeff.sigma_s;
+
+				/* compute accumulated transmittance */
+				float3 transmittance = volume_color_transmittance(sigma_t, dt);
+
+				/* compute emission attenuated by absorption */
+				if(closure_flag & SD_EMISSION) {
+					float3 emission = kernel_volume_emission_integrate(&coeff, closure_flag, transmittance, dt);
+					accum_emission += accum_transmittance * emission;
+				}
+
+				accum_transmittance *= transmittance;
+
+				/* compute pdf for distance sampling */
+				float3 pdf_distance = dt * accum_transmittance * coeff.sigma_s;
+				cdf_distance = cdf_distance + pdf_distance;
+
+				/* write step data */
+				step->sigma_t = sigma_t;
+				step->sigma_s = coeff.sigma_s;
+				step->closure_flag = closure_flag;
+
+				segment->closure_flag |= closure_flag;
+
+				is_last_step_empty = false;
+				segment->numsteps++;
 			}
 			else {
-				/* store empty step */
-				step->sigma_t = make_float3(0.0f, 0.0f, 0.0f);
-				step->sigma_s = make_float3(0.0f, 0.0f, 0.0f);
-				step->closure_flag = 0;
+				if(is_last_step_empty) {
+					/* consecutive empty step, merge */
+					step--;
+				}
+				else {
+					/* store empty step */
+					step->sigma_t = make_float3(0.0f, 0.0f, 0.0f);
+					step->sigma_s = make_float3(0.0f, 0.0f, 0.0f);
+					step->closure_flag = 0;
 
-				segment->numsteps++;
-				is_last_step_empty = true;
+					segment->numsteps++;
+					is_last_step_empty = true;
+				}
 			}
+
+			step->accum_transmittance = accum_transmittance;
+			step->cdf_distance = cdf_distance;
+			step->t = new_t;
+			step->shade_t = t + random_jitter_offset;
+
+			/* stop if at the end of the volume */
+			t = new_t;
+			if(t == ray->t)
+				break;
+
+			/* stop if nearly all light blocked */
+			if(accum_transmittance.x < tp_eps && accum_transmittance.y < tp_eps && accum_transmittance.z < tp_eps)
+				break;
 		}
-
-		step->accum_transmittance = accum_transmittance;
-		step->cdf_distance = cdf_distance;
-		step->t = new_t;
-		step->shade_t = t + random_jitter_offset;
-
-		/* stop if at the end of the volume */
-		t = new_t;
-		if(t == ray->t)
-			break;
-
-		/* stop if nearly all light blocked */
-		if(accum_transmittance.x < tp_eps && accum_transmittance.y < tp_eps && accum_transmittance.z < tp_eps)
-			break;
 	}
 
 	/* store total emission and transmittance */
@@ -798,15 +1037,19 @@ ccl_device void kernel_volume_decoupled_record(KernelGlobals *kg, PathState *sta
 	segment->accum_transmittance = accum_transmittance;
 
 	/* normalize cumulative density function for distance sampling */
-	VolumeStep *last_step = segment->steps + segment->numsteps - 1;
+	int numsteps = segment->numsteps - 1;
 
-	if(!is_zero(last_step->cdf_distance)) {
-		VolumeStep *step = &segment->steps[0];
-		int numsteps = segment->numsteps;
-		float3 inv_cdf_distance_sum = safe_invert_color(last_step->cdf_distance);
+	if(numsteps >= 0) {
+		VolumeStep *last_step = segment->steps + segment->numsteps - 1;
 
-		for(int i = 0; i < numsteps; i++, step++)
-			step->cdf_distance *= inv_cdf_distance_sum;
+		if(!is_zero(last_step->cdf_distance)) {
+			VolumeStep *step = &segment->steps[0];
+			int numsteps = segment->numsteps;
+			float3 inv_cdf_distance_sum = safe_invert_color(last_step->cdf_distance);
+
+			for(int i = 0; i < numsteps; i++, step++)
+				step->cdf_distance *= inv_cdf_distance_sum;
+		}
 	}
 }
 
@@ -1166,6 +1409,7 @@ ccl_device void kernel_volume_stack_init(KernelGlobals *kg,
 			if(need_add) {
 				stack[stack_index].object = stack_sd->object;
 				stack[stack_index].shader = stack_sd->shader;
+				stack[stack_index].volume = stack_sd->prim;
 				++stack_index;
 			}
 		}

intern/cycles/kernel/openvdb/vdb_globals.h

@@ -0,0 +1,45 @@
+/*
+ * Copyright 2016 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 __VDB_GLOBALS_H__
+#define __VDB_GLOBALS_H__
+
+#include "vdb_intern.h"
+
+CCL_NAMESPACE_BEGIN
+
+typedef openvdb::math::Ray<float> vdb_ray_t;
+typedef openvdb::math::Transform vdb_transform_t;
+
+struct OpenVDBGlobals {
+	typedef openvdb::FloatGrid scalar_grid_t;
+	typedef openvdb::Vec3SGrid vector_grid_t;
+	typedef openvdb::tools::VolumeRayIntersector<scalar_grid_t, scalar_grid_t::TreeType::RootNodeType::ChildNodeType::LEVEL, vdb_ray_t> scalar_isector_t;
+	typedef openvdb::tools::VolumeRayIntersector<vector_grid_t, vector_grid_t::TreeType::RootNodeType::ChildNodeType::LEVEL, vdb_ray_t> vector_isector_t;
+	
+	vector<const scalar_grid_t *> scalar_grids;
+	vector<const vector_grid_t *> vector_grids;
+	/* Main intersectors, which initialize the voxels' bounding box
+	 * so the ones for the various threads do not do this,
+	 * rather they are generated from a copy of these
+	 */
+	vector<scalar_isector_t *> scalar_main_isectors;
+	vector<vector_isector_t *> vector_main_isectors;
+};
+
+CCL_NAMESPACE_END
+
+#endif /* __VDB_GLOBALS_H__ */

intern/cycles/kernel/openvdb/vdb_intern.h

@@ -0,0 +1,48 @@
+/*
+ * Copyright 2016 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 __VDB_INTERN_H__
+#define __VDB_INTERN_H__
+
+/* They are too many implicit float conversions happening in OpenVDB, disabling
+ * errors for now (kevin) */
+#ifdef __GNUC__
+#	pragma GCC diagnostic push
+#	pragma GCC diagnostic ignored "-Wfloat-conversion"
+#	pragma GCC diagnostic ignored "-Wdouble-promotion"
+#endif
+
+#include <openvdb/openvdb.h>
+#include <openvdb/tools/Interpolation.h>
+#include <openvdb/tools/RayIntersector.h>
+
+#ifdef __GNUC__
+#	pragma GCC diagnostic pop
+#endif
+
+#include "util/util_vector.h"
+
+CCL_NAMESPACE_BEGIN
+
+#if defined(HAS_CPP11_FEATURES)
+using std::isfinite;
+#else
+using boost::math::isfinite;
+#endif
+
+CCL_NAMESPACE_END
+
+#endif /* __VDB_INTERN_H__ */

intern/cycles/kernel/openvdb/vdb_thread.cpp

@@ -0,0 +1,228 @@
+/*
+ * Copyright 2016 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 "kernel_compat_cpu.h"
+#include "kernel_types.h"
+#include "kernel_globals.h"
+
+#include "vdb_globals.h"
+#include "vdb_intern.h"
+#include "vdb_thread.h"
+
+#include "util_vector.h"
+
+CCL_NAMESPACE_BEGIN
+
+/* Manage thread-local data associated with volumes */
+
+struct OpenVDBScalarThreadData {
+	typedef openvdb::FloatGrid grid_t;
+	typedef openvdb::FloatGrid::ConstAccessor accessor_t;
+	typedef openvdb::tools::GridSampler<accessor_t, openvdb::tools::PointSampler> point_sampler_t;
+	typedef openvdb::tools::GridSampler<accessor_t, openvdb::tools::BoxSampler> box_sampler_t;
+	typedef openvdb::tools::VolumeRayIntersector<grid_t, grid_t::TreeType::RootNodeType::ChildNodeType::LEVEL, vdb_ray_t> isector_t;
+	
+	void init(const grid_t &grid, const isector_t &main_isector)
+	{
+		accessor = new accessor_t(grid.getConstAccessor());
+		point_sampler = new point_sampler_t(*accessor, grid.transform());
+		box_sampler = new box_sampler_t(*accessor, grid.transform());
+		isector = new isector_t(main_isector);
+	}
+	
+	void free()
+	{
+		delete accessor;
+		delete point_sampler;
+		delete box_sampler;
+		delete isector;
+	}
+	
+	accessor_t *accessor;
+	point_sampler_t *point_sampler;
+	box_sampler_t *box_sampler;
+	isector_t *isector;
+};
+
+struct OpenVDBVectorThreadData {
+	typedef openvdb::Vec3SGrid grid_t;
+	typedef openvdb::Vec3SGrid::ConstAccessor accessor_t;
+	typedef openvdb::tools::GridSampler<accessor_t, openvdb::tools::PointSampler> point_sampler_t;
+	typedef openvdb::tools::GridSampler<accessor_t, openvdb::tools::BoxSampler> box_sampler_t;
+	typedef openvdb::tools::GridSampler<accessor_t, openvdb::tools::StaggeredPointSampler> stag_point_sampler_t;
+	typedef openvdb::tools::GridSampler<accessor_t, openvdb::tools::StaggeredBoxSampler> stag_box_sampler_t;
+	
+	typedef openvdb::tools::VolumeRayIntersector<grid_t, grid_t::TreeType::RootNodeType::ChildNodeType::LEVEL, vdb_ray_t> isector_t;
+	
+	void init (const grid_t &grid, const isector_t &main_isector)
+	{
+		accessor = new accessor_t(grid.getConstAccessor());
+		point_sampler = new point_sampler_t(*accessor, grid.transform());
+		box_sampler = new box_sampler_t(*accessor, grid.transform());
+		stag_point_sampler = new stag_point_sampler_t(*accessor, grid.transform());
+		stag_box_sampler = new stag_box_sampler_t(*accessor, grid.transform());
+		isector = new isector_t(main_isector);
+	}
+	
+	void free()
+	{
+		delete accessor;
+		delete point_sampler;
+		delete box_sampler;
+		delete stag_point_sampler;
+		delete stag_box_sampler;
+		delete isector;
+	}
+	
+	accessor_t *accessor;
+	point_sampler_t *point_sampler;
+	box_sampler_t *box_sampler;
+	stag_point_sampler_t *stag_point_sampler;
+	stag_box_sampler_t *stag_box_sampler;
+	isector_t *isector;
+};
+
+struct OpenVDBThreadData {
+	std::vector<OpenVDBScalarThreadData> scalar_data;
+	std::vector<OpenVDBVectorThreadData> vector_data;
+};
+
+void vdb_thread_init(KernelGlobals *kg, const KernelGlobals *kernel_globals, OpenVDBGlobals *vdb_globals)
+{
+	kg->vdb = vdb_globals;
+	
+	OpenVDBThreadData *tdata = new OpenVDBThreadData;
+	
+	tdata->scalar_data.resize(vdb_globals->scalar_grids.size());
+	tdata->vector_data.resize(vdb_globals->vector_grids.size());
+	for (size_t i = 0; i < vdb_globals->scalar_grids.size(); ++i) {
+		tdata->scalar_data[i].init(*vdb_globals->scalar_grids[i], *vdb_globals->scalar_main_isectors[i]);
+	}
+	for (size_t i = 0; i < vdb_globals->vector_grids.size(); ++i) {
+		tdata->vector_data[i].init(*vdb_globals->vector_grids[i], *vdb_globals->vector_main_isectors[i]);
+	}
+	kg->vdb_tdata = tdata;
+}
+
+void vdb_thread_free(KernelGlobals *kg)
+{
+	OpenVDBThreadData *tdata = kg->vdb_tdata;
+	kg->vdb_tdata = NULL;
+	
+	for (size_t i = 0; i < tdata->scalar_data.size(); ++i) {
+		tdata->scalar_data[i].free();
+	}
+	for (size_t i = 0; i < tdata->vector_data.size(); ++i) {
+		tdata->vector_data[i].free();
+	}
+	delete tdata;
+}
+
+bool vdb_volume_scalar_has_uniform_voxels(OpenVDBGlobals *vdb, int vdb_index)
+{
+	return vdb->scalar_grids[vdb_index]->hasUniformVoxels();
+}
+
+bool vdb_volume_vector_has_uniform_voxels(OpenVDBGlobals *vdb, int vdb_index)
+{
+	return vdb->vector_grids[vdb_index]->hasUniformVoxels();
+}
+
+float vdb_volume_sample_scalar(OpenVDBGlobals */*vdb*/, OpenVDBThreadData *vdb_thread, int vdb_index,
+                               float x, float y, float z, int sampling)
+{
+	OpenVDBScalarThreadData &data = vdb_thread->scalar_data[vdb_index];
+	
+	switch (sampling) {
+		case OPENVDB_SAMPLE_POINT:
+			return data.point_sampler->wsSample(openvdb::Vec3d(x, y, z));
+		case OPENVDB_SAMPLE_BOX:
+			return data.box_sampler->wsSample(openvdb::Vec3d(x, y, z));
+	}
+	
+	return 0.0f;
+}
+
+float3 vdb_volume_sample_vector(OpenVDBGlobals *vdb, OpenVDBThreadData *vdb_thread, int vdb_index,
+                                float x, float y, float z, int sampling)
+{
+	bool staggered = (vdb->vector_grids[vdb_index]->getGridClass() == openvdb::GRID_STAGGERED);
+	OpenVDBVectorThreadData &data = vdb_thread->vector_data[vdb_index];
+	openvdb::Vec3s r;
+	
+	if (staggered) {
+		switch (sampling) {
+			case OPENVDB_SAMPLE_POINT:
+				r = data.stag_point_sampler->wsSample(openvdb::Vec3d(x, y, z));
+				break;
+			case OPENVDB_SAMPLE_BOX:
+				r = data.stag_box_sampler->wsSample(openvdb::Vec3d(x, y, z));
+				break;
+		}
+	}
+	else {
+		switch (sampling) {
+			case OPENVDB_SAMPLE_POINT:
+				r = data.point_sampler->wsSample(openvdb::Vec3d(x, y, z));
+				break;
+			case OPENVDB_SAMPLE_BOX:
+				r = data.box_sampler->wsSample(openvdb::Vec3d(x, y, z));
+				break;
+		}
+	}
+	
+	return make_float3(r.x(), r.y(), r.z());
+}
+
+bool vdb_volume_intersect(OpenVDBThreadData *vdb_thread, int vdb_index,
+                          const Ray *ray, float *isect)
+{
+	OpenVDBScalarThreadData &data = vdb_thread->scalar_data[vdb_index];
+	
+	vdb_ray_t::Vec3Type P(ray->P.x, ray->P.y, ray->P.z);
+	vdb_ray_t::Vec3Type D(ray->D.x, ray->D.y, ray->D.z);
+	D.normalize();
+	
+	vdb_ray_t vdb_ray(P, D, 1e-5f, ray->t);
+	
+	if(data.isector->setWorldRay(vdb_ray)) {
+		// TODO(kevin): is this correct?
+		*isect = static_cast<float>(vdb_ray.t1());
+		
+		return true;
+	}
+	
+	return false;
+}
+
+bool vdb_volume_march(OpenVDBThreadData *vdb_thread, int vdb_index,
+                      float *t0, float *t1)
+{
+	OpenVDBScalarThreadData &data = vdb_thread->scalar_data[vdb_index];
+	
+	float vdb_t0(*t0), vdb_t1(*t1);
+
+	if(data.isector->march(vdb_t0, vdb_t1)) {
+		*t0 = data.isector->getWorldTime(vdb_t0);
+		*t1 = data.isector->getWorldTime(vdb_t1);
+
+		return true;
+	}
+
+	return false;
+}
+
+CCL_NAMESPACE_END

intern/cycles/kernel/openvdb/vdb_thread.h

@@ -0,0 +1,49 @@
+/*
+ * Copyright 2016 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 __VDB_THREAD_H__
+#define __VDB_THREAD_H__
+
+CCL_NAMESPACE_BEGIN
+
+struct Intersection;
+struct KernelGlobals;
+struct OpenVDBGlobals;
+struct OpenVDBThreadData;
+struct Ray;
+
+void vdb_thread_init(KernelGlobals *kg, const KernelGlobals *kernel_globals, OpenVDBGlobals *vdb_globals);
+void vdb_thread_free(KernelGlobals *kg);
+
+enum OpenVDB_SampleType {
+	OPENVDB_SAMPLE_POINT = 0,
+	OPENVDB_SAMPLE_BOX   = 1,
+};
+
+bool vdb_volume_scalar_has_uniform_voxels(OpenVDBGlobals *vdb, int vdb_index);
+bool vdb_volume_vector_has_uniform_voxels(OpenVDBGlobals *vdb, int vdb_index);
+float vdb_volume_sample_scalar(OpenVDBGlobals *vdb, OpenVDBThreadData *vdb_thread, int vdb_index,
+                               float x, float y, float z, int sampling);
+float3 vdb_volume_sample_vector(OpenVDBGlobals *vdb, OpenVDBThreadData *vdb_thread, int vdb_index,
+                                float x, float y, float z, int sampling);
+bool vdb_volume_intersect(OpenVDBThreadData *vdb_thread, int vdb_index,
+                          const Ray *ray, float *isect);
+bool vdb_volume_march(OpenVDBThreadData *vdb_thread, int vdb_index,
+                      float *t0, float *t1);
+
+CCL_NAMESPACE_END
+
+#endif /* __VDB_THREAD_H__ */

intern/cycles/kernel/osl/osl_services.cpp

@@ -34,10 +34,6 @@
 #include "kernel/osl/osl_services.h"
 #include "kernel/osl/osl_shader.h"
 
-#include "util/util_foreach.h"
-#include "util/util_logging.h"
-#include "util/util_string.h"
-
 #include "kernel/kernel_compat_cpu.h"
 #include "kernel/split/kernel_split_data_types.h"
 #include "kernel/kernel_globals.h"
@@ -50,6 +46,14 @@
 #include "kernel/geom/geom.h"
 #include "kernel/bvh/bvh.h"
 
+/* Note: "util_foreach.h" needs to be included after "kernel_compat_cpu.h", as
+ * for some reason ccl::foreach conflicts with openvdb::tools::foreach, which is
+ * indirectly included through "kernel_compat_cpu.h".
+ */
+#include "util/util_foreach.h"
+#include "util/util_logging.h"
+#include "util/util_string.h"
+
 #include "kernel/kernel_projection.h"
 #include "kernel/kernel_accumulate.h"
 #include "kernel/kernel_shader.h"

intern/cycles/render/CMakeLists.txt

@@ -8,6 +8,12 @@ set(INC_SYS
 	${GLEW_INCLUDE_DIR}
 )
 
+if(WITH_OPENVDB)
+	list(APPEND INC_SYS
+		${OPENVDB_INCLUDE_DIRS}
+	)
+endif()
+
 set(SRC
 	attribute.cpp
 	background.cpp
@@ -35,6 +41,7 @@ set(SRC
 	svm.cpp
 	tables.cpp
 	tile.cpp
+	volume.cpp
 )
 
 set(SRC_HEADERS
@@ -62,6 +69,7 @@ set(SRC_HEADERS
 	svm.h
 	tables.h
 	tile.h
+	volume.h
 )
 
 set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${RTTI_DISABLE_FLAGS}")

intern/cycles/render/attribute.cpp

@@ -33,7 +33,8 @@ Attribute::~Attribute()
 		VoxelAttribute *voxel_data = data_voxel();
 
 		if(voxel_data && voxel_data->slot != -1) {
-			voxel_data->manager->remove_image(voxel_data->slot);
+			if (voxel_data->manager)
+				voxel_data->manager->remove_image(voxel_data->slot);
 		}
 	}
 }
@@ -344,10 +345,13 @@ Attribute *AttributeSet::add(ustring name, TypeDesc type, AttributeElement eleme
 	/* this is weak .. */
 	if(triangle_mesh)
 		attr->resize(triangle_mesh, ATTR_PRIM_TRIANGLE, false);
-	if(curve_mesh)
+	else if(curve_mesh)
 		attr->resize(curve_mesh, ATTR_PRIM_CURVE, false);
-	if(subd_mesh)
+	else if(subd_mesh)
 		attr->resize(subd_mesh, ATTR_PRIM_SUBD, false);
+	else if(element == ATTR_ELEMENT_VOXEL) {
+		attr->resize(1);
+	}
 
 	return attr;
 }

intern/cycles/render/scene.cpp

@@ -32,6 +32,7 @@
 #include "render/shader.h"
 #include "render/svm.h"
 #include "render/tables.h"
+#include "volume.h"
 
 #include "util/util_foreach.h"
 #include "util/util_guarded_allocator.h"
@@ -58,6 +59,7 @@ Scene::Scene(const SceneParams& params_, const DeviceInfo& device_info_)
 	particle_system_manager = new ParticleSystemManager();
 	curve_system_manager = new CurveSystemManager();
 	bake_manager = new BakeManager();
+	volume_manager = new VolumeManager();
 
 	/* OSL only works on the CPU */
 	if(device_info_.type == DEVICE_CPU)
@@ -83,12 +85,15 @@ void Scene::free_memory(bool final)
 		delete l;
 	foreach(ParticleSystem *p, particle_systems)
 		delete p;
+	foreach(Volume *v, volumes)
+		delete v;
 
 	shaders.clear();
 	meshes.clear();
 	objects.clear();
 	lights.clear();
 	particle_systems.clear();
+	volumes.clear();
 
 	if(device) {
 		camera->device_free(device, &dscene, this);
@@ -112,6 +117,7 @@ void Scene::free_memory(bool final)
 			image_manager->device_free_builtin(device, &dscene);
 
 		lookup_tables->device_free(device, &dscene);
+		volume_manager->device_free(device, &dscene);
 	}
 
 	if(final) {
@@ -128,6 +134,7 @@ void Scene::free_memory(bool final)
 		delete curve_system_manager;
 		delete image_manager;
 		delete bake_manager;
+		delete volume_manager;
 	}
 }
 
@@ -234,6 +241,11 @@ void Scene::device_update(Device *device_, Progress& progress)
 
 	if(progress.get_cancel() || device->have_error()) return;
 
+	progress.set_status("Updating OpenVDB Volumes");
+	volume_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
 	if(device->have_error() == false) {
 		progress.set_status("Updating Device", "Writing constant memory");
 		device->const_copy_to("__data", &dscene.data, sizeof(dscene.data));
@@ -306,6 +318,7 @@ bool Scene::need_data_update()
 		|| particle_system_manager->need_update
 		|| curve_system_manager->need_update
 		|| bake_manager->need_update
+		|| volume_manager->need_update
 		|| film->need_update);
 }
 

intern/cycles/render/scene.cpp.orig

@@ -0,0 +1,377 @@
+/*
+ * Copyright 2011-2013 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 <stdlib.h>
+
+<<<<<<< HEAD
+#include "background.h"
+#include "bake.h"
+#include "camera.h"
+#include "curves.h"
+#include "device.h"
+#include "film.h"
+#include "integrator.h"
+#include "light.h"
+#include "mesh.h"
+#include "object.h"
+#include "osl.h"
+#include "particles.h"
+#include "scene.h"
+#include "shader.h"
+#include "svm.h"
+#include "tables.h"
+#include "volume.h"
+
+#include "util_foreach.h"
+#include "util_guarded_allocator.h"
+#include "util_logging.h"
+#include "util_progress.h"
+=======
+#include "render/background.h"
+#include "render/bake.h"
+#include "render/camera.h"
+#include "render/curves.h"
+#include "device/device.h"
+#include "render/film.h"
+#include "render/integrator.h"
+#include "render/light.h"
+#include "render/mesh.h"
+#include "render/object.h"
+#include "render/osl.h"
+#include "render/particles.h"
+#include "render/scene.h"
+#include "render/shader.h"
+#include "render/svm.h"
+#include "render/tables.h"
+
+#include "util/util_foreach.h"
+#include "util/util_guarded_allocator.h"
+#include "util/util_logging.h"
+#include "util/util_progress.h"
+>>>>>>> origin/master
+
+CCL_NAMESPACE_BEGIN
+
+Scene::Scene(const SceneParams& params_, const DeviceInfo& device_info_)
+: params(params_)
+{
+	device = NULL;
+	memset(&dscene.data, 0, sizeof(dscene.data));
+
+	camera = new Camera();
+	lookup_tables = new LookupTables();
+	film = new Film();
+	background = new Background();
+	light_manager = new LightManager();
+	mesh_manager = new MeshManager();
+	object_manager = new ObjectManager();
+	integrator = new Integrator();
+	image_manager = new ImageManager(device_info_);
+	particle_system_manager = new ParticleSystemManager();
+	curve_system_manager = new CurveSystemManager();
+	bake_manager = new BakeManager();
+	volume_manager = new VolumeManager();
+
+	/* OSL only works on the CPU */
+	if(device_info_.type == DEVICE_CPU)
+		shader_manager = ShaderManager::create(this, params.shadingsystem);
+	else
+		shader_manager = ShaderManager::create(this, SHADINGSYSTEM_SVM);
+}
+
+Scene::~Scene()
+{
+	free_memory(true);
+}
+
+void Scene::free_memory(bool final)
+{
+	foreach(Shader *s, shaders)
+		delete s;
+	foreach(Mesh *m, meshes)
+		delete m;
+	foreach(Object *o, objects)
+		delete o;
+	foreach(Light *l, lights)
+		delete l;
+	foreach(ParticleSystem *p, particle_systems)
+		delete p;
+	foreach(Volume *v, volumes)
+		delete v;
+
+	shaders.clear();
+	meshes.clear();
+	objects.clear();
+	lights.clear();
+	particle_systems.clear();
+	volumes.clear();
+
+	if(device) {
+		camera->device_free(device, &dscene, this);
+		film->device_free(device, &dscene, this);
+		background->device_free(device, &dscene);
+		integrator->device_free(device, &dscene);
+
+		object_manager->device_free(device, &dscene);
+		mesh_manager->device_free(device, &dscene);
+		shader_manager->device_free(device, &dscene, this);
+		light_manager->device_free(device, &dscene);
+
+		particle_system_manager->device_free(device, &dscene);
+		curve_system_manager->device_free(device, &dscene);
+
+		bake_manager->device_free(device, &dscene);
+
+		if(!params.persistent_data || final)
+			image_manager->device_free(device, &dscene);
+		else
+			image_manager->device_free_builtin(device, &dscene);
+
+		lookup_tables->device_free(device, &dscene);
+		volume_manager->device_free(device, &dscene);
+	}
+
+	if(final) {
+		delete lookup_tables;
+		delete camera;
+		delete film;
+		delete background;
+		delete integrator;
+		delete object_manager;
+		delete mesh_manager;
+		delete shader_manager;
+		delete light_manager;
+		delete particle_system_manager;
+		delete curve_system_manager;
+		delete image_manager;
+		delete bake_manager;
+		delete volume_manager;
+	}
+}
+
+void Scene::device_update(Device *device_, Progress& progress)
+{
+	if(!device)
+		device = device_;
+
+	bool print_stats = need_data_update();
+
+	/* The order of updates is important, because there's dependencies between
+	 * the different managers, using data computed by previous managers.
+	 *
+	 * - Image manager uploads images used by shaders.
+	 * - Camera may be used for adaptive subdivision.
+	 * - Displacement shader must have all shader data available.
+	 * - Light manager needs lookup tables and final mesh data to compute emission CDF.
+	 * - Film needs light manager to run for use_light_visibility
+	 * - Lookup tables are done a second time to handle film tables
+	 */
+
+	progress.set_status("Updating Shaders");
+	shader_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Background");
+	background->device_update(device, &dscene, this);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Camera");
+	camera->device_update(device, &dscene, this);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Meshes Flags");
+	mesh_manager->device_update_flags(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Objects");
+	object_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Meshes");
+	mesh_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Objects Flags");
+	object_manager->device_update_flags(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Images");
+	image_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Camera Volume");
+	camera->device_update_volume(device, &dscene, this);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Hair Systems");
+	curve_system_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Lookup Tables");
+	lookup_tables->device_update(device, &dscene);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Lights");
+	light_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Particle Systems");
+	particle_system_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Integrator");
+	integrator->device_update(device, &dscene, this);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Film");
+	film->device_update(device, &dscene, this);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Lookup Tables");
+	lookup_tables->device_update(device, &dscene);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating Baking");
+	bake_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	progress.set_status("Updating OpenVDB Volumes");
+	volume_manager->device_update(device, &dscene, this, progress);
+
+	if(progress.get_cancel() || device->have_error()) return;
+
+	if(device->have_error() == false) {
+		progress.set_status("Updating Device", "Writing constant memory");
+		device->const_copy_to("__data", &dscene.data, sizeof(dscene.data));
+	}
+
+	if(print_stats) {
+		size_t mem_used = util_guarded_get_mem_used();
+		size_t mem_peak = util_guarded_get_mem_peak();
+
+		VLOG(1) << "System memory statistics after full device sync:\n"
+		        << "  Usage: " << string_human_readable_number(mem_used)
+		        << " (" << string_human_readable_size(mem_used) << ")\n"
+		        << "  Peak: " << string_human_readable_number(mem_peak)
+		        << " (" << string_human_readable_size(mem_peak) << ")";
+	}
+}
+
+Scene::MotionType Scene::need_motion(bool advanced_shading)
+{
+	if(integrator->motion_blur)
+		return (advanced_shading)? MOTION_BLUR: MOTION_NONE;
+	else if(Pass::contains(film->passes, PASS_MOTION))
+		return MOTION_PASS;
+	else
+		return MOTION_NONE;
+}
+
+float Scene::motion_shutter_time()
+{
+	if(need_motion() == Scene::MOTION_PASS)
+		return 2.0f;
+	else
+		return camera->shuttertime;
+}
+
+bool Scene::need_global_attribute(AttributeStandard std)
+{
+	if(std == ATTR_STD_UV)
+		return Pass::contains(film->passes, PASS_UV);
+	else if(std == ATTR_STD_MOTION_VERTEX_POSITION)
+		return need_motion() != MOTION_NONE;
+	else if(std == ATTR_STD_MOTION_VERTEX_NORMAL)
+		return need_motion() == MOTION_BLUR;
+	
+	return false;
+}
+
+void Scene::need_global_attributes(AttributeRequestSet& attributes)
+{
+	for(int std = ATTR_STD_NONE; std < ATTR_STD_NUM; std++)
+		if(need_global_attribute((AttributeStandard)std))
+			attributes.add((AttributeStandard)std);
+}
+
+bool Scene::need_update()
+{
+	return (need_reset() || film->need_update);
+}
+
+bool Scene::need_data_update()
+{
+	return (background->need_update
+		|| image_manager->need_update
+		|| object_manager->need_update
+		|| mesh_manager->need_update
+		|| light_manager->need_update
+		|| lookup_tables->need_update
+		|| integrator->need_update
+		|| shader_manager->need_update
+		|| particle_system_manager->need_update
+		|| curve_system_manager->need_update
+		|| bake_manager->need_update
+		|| volume_manager->need_update
+		|| film->need_update);
+}
+
+bool Scene::need_reset()
+{
+	return need_data_update() || camera->need_update;
+}
+
+void Scene::reset()
+{
+	shader_manager->reset(this);
+	shader_manager->add_default(this);
+
+	/* ensure all objects are updated */
+	camera->tag_update();
+	film->tag_update(this);
+	background->tag_update(this);
+	integrator->tag_update(this);
+	object_manager->tag_update(this);
+	mesh_manager->tag_update(this);
+	light_manager->tag_update(this);
+	particle_system_manager->tag_update(this);
+	curve_system_manager->tag_update(this);
+}
+
+void Scene::device_free()
+{
+	free_memory(false);
+}
+
+CCL_NAMESPACE_END
+

intern/cycles/render/scene.h

@@ -54,6 +54,8 @@ class ShaderManager;
 class Progress;
 class BakeManager;
 class BakeData;
+class Volume;
+class VolumeManager;
 
 /* Scene Device Data */
 
@@ -121,6 +123,9 @@ public:
 	vector<device_vector<uchar>* > tex_byte_image;
 	vector<device_vector<half>* > tex_half_image;
 
+	/* volume */
+	device_vector<uint> vol_shader;
+
 	KernelData data;
 };
 
@@ -182,6 +187,7 @@ public:
 	vector<Shader*> shaders;
 	vector<Light*> lights;
 	vector<ParticleSystem*> particle_systems;
+	vector<Volume*> volumes;
 
 	/* data managers */
 	ImageManager *image_manager;
@@ -192,12 +198,14 @@ public:
 	ParticleSystemManager *particle_system_manager;
 	CurveSystemManager *curve_system_manager;
 	BakeManager *bake_manager;
+	VolumeManager *volume_manager;
 
 	/* default shaders */
 	Shader *default_surface;
 	Shader *default_light;
 	Shader *default_background;
 	Shader *default_empty;
+	Shader *default_volume;
 
 	/* device */
 	Device *device;

intern/cycles/render/session.cpp

@@ -701,11 +701,12 @@ DeviceRequestedFeatures Session::get_requested_device_features()
 	requested_features.use_camera_motion = scene->camera->use_motion;
 	foreach(Object *object, scene->objects) {
 		Mesh *mesh = object->mesh;
-		if(mesh->num_curves()) {
-			requested_features.use_hair = true;
+		if(mesh) {
+			if(mesh->num_curves()) {
+				requested_features.use_hair = true;
+			}
 		}
-		requested_features.use_object_motion |= object->use_motion | mesh->use_motion_blur;
-		requested_features.use_camera_motion |= mesh->use_motion_blur;
+		requested_features.use_object_motion |= object->use_motion;
 #ifdef WITH_OPENSUBDIV
 		if(mesh->subdivision_type != Mesh::SUBDIVISION_NONE) {
 			requested_features.use_patch_evaluation = true;

intern/cycles/render/shader.cpp

@@ -398,6 +398,7 @@ void ShaderManager::device_update_shaders_used(Scene *scene)
 	scene->default_light->used = true;
 	scene->default_background->used = true;
 	scene->default_empty->used = true;
+	scene->default_volume->used = true;
 
 	if(scene->background->shader)
 		scene->background->shader->used = true;
@@ -570,6 +571,23 @@ void ShaderManager::add_default(Scene *scene)
 		scene->shaders.push_back(shader);
 		scene->default_empty = shader;
 	}
+
+	/* default empty */
+	{
+		ShaderGraph *graph = new ShaderGraph();
+
+		ScatterVolumeNode *scatter = new ScatterVolumeNode();
+		scatter->input("Density")->set(1.0);
+		graph->add(scatter);
+
+		graph->connect(scatter->output("Volume"), graph->output()->input("Volume"));
+
+		Shader *shader = new Shader();
+		shader->name = "default_volume";
+		shader->graph = graph;
+		scene->shaders.push_back(shader);
+		scene->default_volume = shader;
+	}
 }
 
 void ShaderManager::get_requested_graph_features(ShaderGraph *graph,

intern/cycles/render/svm.cpp

@@ -59,7 +59,7 @@ void SVMShaderManager::device_update_shader(Scene *scene,
 	svm_nodes.push_back(make_int4(NODE_SHADER_JUMP, 0, 0, 0));
 
 	SVMCompiler::Summary summary;
-	SVMCompiler compiler(scene->shader_manager, scene->image_manager);
+	SVMCompiler compiler(scene->shader_manager, scene->image_manager, scene->volume_manager);
 	compiler.background = (shader == scene->default_background);
 	compiler.compile(scene, shader, svm_nodes, 0, &summary);
 
@@ -156,10 +156,11 @@ void SVMShaderManager::device_free(Device *device, DeviceScene *dscene, Scene *s
 
 /* Graph Compiler */
 
-SVMCompiler::SVMCompiler(ShaderManager *shader_manager_, ImageManager *image_manager_)
+SVMCompiler::SVMCompiler(ShaderManager *shader_manager_, ImageManager *image_manager_, VolumeManager *volume_manager_)
 {
 	shader_manager = shader_manager_;
 	image_manager = image_manager_;
+	volume_manager = volume_manager_;
 	max_stack_use = 0;
 	current_type = SHADER_TYPE_SURFACE;
 	current_shader = NULL;

intern/cycles/render/svm.h

@@ -35,6 +35,7 @@ class ShaderGraph;
 class ShaderInput;
 class ShaderNode;
 class ShaderOutput;
+class VolumeManager;
 
 /* Shader Manager */
 
@@ -95,7 +96,8 @@ public:
 		string full_report() const;
 	};
 
-	SVMCompiler(ShaderManager *shader_manager, ImageManager *image_manager);
+	SVMCompiler(ShaderManager *shader_manager, ImageManager *image_manager, VolumeManager *volume_manager_);
+
 	void compile(Scene *scene,
 	             Shader *shader,
 	             vector<int4>& svm_nodes,
@@ -123,6 +125,7 @@ public:
 	ShaderType output_type() { return current_type; }
 
 	ImageManager *image_manager;
+	VolumeManager *volume_manager;
 	ShaderManager *shader_manager;
 	bool background;
 

intern/cycles/render/volume.cpp

@@ -0,0 +1,426 @@
+/*
+ * Copyright 2015 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 "scene.h"
+#include "volume.h"
+
+#include "util/util_foreach.h"
+#include "util/util_logging.h"
+#include "util/util_progress.h"
+#include "util/util_task.h"
+
+#include "../kernel/openvdb/vdb_globals.h"
+
+CCL_NAMESPACE_BEGIN
+
+#define MAX_VOLUME 1024
+
+void Volume::tag_update(Scene *scene, bool /*rebuild*/)
+{
+	scene->volume_manager->need_update = true;
+}
+
+/* ------------------------------------------------------------------------- */
+
+VolumeManager::VolumeManager()
+{
+#ifdef WITH_OPENVDB
+	openvdb::initialize();
+
+	current_grids.reserve(64);
+#endif
+
+	need_update = true;
+	num_float_volume = 0;
+	num_float3_volume = 0;
+}
+
+VolumeManager::~VolumeManager()
+{
+	current_grids.clear();
+}
+
+static inline void catch_exceptions()
+{
+#ifdef WITH_OPENVDB
+	try {
+		throw;
+	}
+	catch(const openvdb::IoError& e) {
+		std::cerr << e.what() << "\n";
+	}
+#endif
+}
+
+int VolumeManager::add_volume(Volume *volume, const std::string &filename, const std::string &name)
+{
+	size_t slot = -1;
+
+	if((slot = find_existing_slot(volume, filename, name)) != -1) {
+		return slot;
+	}
+
+	if((num_float_volume + num_float3_volume + 1) > MAX_VOLUME) {
+		printf("VolumeManager::add_volume: volume limit reached %d!\n", MAX_VOLUME);
+		return -1;
+	}
+
+	try {
+		if(is_openvdb_file(filename)) {
+			slot = add_openvdb_volume(volume, filename, name);
+		}
+
+		add_grid_description(volume, filename, name, slot);
+
+		volumes.push_back(volume);
+	}
+	catch(...) {
+		catch_exceptions();
+		slot = -1;
+	}
+
+	return slot;
+}
+
+int VolumeManager::find_existing_slot(Volume *volume, const std::string &filename, const std::string &name)
+{
+	for(size_t i = 0; i < current_grids.size(); ++i) {
+		GridDescription grid = current_grids[i];
+
+		if(grid.volume == volume) {
+			if(grid.filename == filename && grid.name == name) {
+				return grid.slot;
+			}
+		}
+	}
+
+	return -1;
+}
+
+int VolumeManager::find_density_slot()
+{
+	/* first try finding a matching grid name */
+	for(size_t i = 0; i < current_grids.size(); ++i) {
+		GridDescription grid = current_grids[i];
+		
+		if(string_iequals(grid.name, "density") || string_iequals(grid.name, "density high"))
+			return grid.slot;
+	}
+	
+	/* try using the first scalar float grid instead */
+	for (size_t i = 0; i < volumes.size(); ++i) {
+		Volume *volume = volumes[i];
+		
+		if (!volume->scalar_grids.empty()) {
+			return 0;
+		}
+	}
+	
+	return -1;
+}
+
+bool VolumeManager::is_openvdb_file(const string& filename) const
+{
+	return string_endswith(filename, ".vdb");
+}
+
+template <typename Container>
+size_t find_empty_slot(Container container)
+{
+	size_t slot = 0;
+
+	for(; slot < container.size(); ++slot) {
+		if(!container[slot]) {
+			break;
+		}
+	}
+
+	if(slot == container.size()) {
+		if(slot == MAX_VOLUME) {
+			printf("VolumeManager::add_volume: volume limit reached %d!\n",
+			       MAX_VOLUME);
+			return -1;
+		}
+
+		container.resize(slot + 1);
+	}
+
+	return slot;
+}
+
+size_t VolumeManager::add_openvdb_volume(Volume *volume, const std::string &filename, const std::string &name)
+{
+	size_t slot = -1;
+
+#ifdef WITH_OPENVDB
+	openvdb::io::File file(filename);
+	file.open();
+
+	if(!file.hasGrid(name)) return -1;
+
+	openvdb::GridBase::Ptr grid = file.readGrid(name);
+	if(grid->getGridClass() == openvdb::GRID_LEVEL_SET) return -1;
+
+	if(grid->isType<openvdb::FloatGrid>()) {
+		openvdb::FloatGrid::Ptr fgrid = openvdb::gridPtrCast<openvdb::FloatGrid>(grid);
+		volume->scalar_grids.push_back(fgrid);
+
+		/* XXX Ray intersectors only support uniform grids.
+		 * Can we make this transparent somehow? - lukas
+		 */
+		assert(fgrid->hasUniformVoxels());
+
+		slot = num_float_volume++;
+	}
+	else if(grid->isType<openvdb::Vec3SGrid>()) {
+		openvdb::Vec3SGrid::Ptr vgrid = openvdb::gridPtrCast<openvdb::Vec3SGrid>(grid);
+		volume->vector_grids.push_back(vgrid);
+
+		slot = num_float3_volume++;
+	}
+#else
+	(void)volume;
+	(void)filename;
+	(void)name;
+#endif
+
+	return slot;
+}
+
+void VolumeManager::add_grid_description(Volume *volume, const std::string &filename, const std::string &name, int slot)
+{
+	GridDescription descr;
+	descr.filename = filename;
+	descr.name = name;
+	descr.volume = volume;
+	descr.slot = slot;
+
+	current_grids.push_back(descr);
+}
+
+static void update_attribute_element_offset(Attribute *vattr,
+                                            TypeDesc& type,
+                                            int& offset,
+                                            AttributeElement& element)
+{
+	if(vattr) {
+		/* store element and type */
+		element = vattr->element;
+		type = vattr->type;
+
+		/* store slot in offset value */
+		VoxelAttribute *voxel_data = vattr->data_voxel();
+		offset = voxel_data->slot;
+	}
+	else {
+		/* attribute not found */
+		element = ATTR_ELEMENT_NONE;
+		offset = 0;
+	}
+}
+
+void VolumeManager::device_update_attributes(Device *device, DeviceScene *dscene, Scene *scene, Progress& progress)
+{
+	progress.set_status("Updating Volume", "Computing attributes");
+
+	/* gather per volume requested attributes. as volumes may have multiple
+	 * shaders assigned, this merges the requested attributes that have
+	 * been set per shader by the shader manager */
+	vector<AttributeRequestSet> volume_attributes(volumes.size());
+
+	for(size_t i = 0; i < volumes.size(); i++) {
+		Volume *volume = volumes[i];
+
+		foreach(Shader *shader, volume->used_shaders) {
+			volume_attributes[i].add(shader->attributes);
+		}
+	}
+
+	for(size_t i = 0; i < volumes.size(); i++) {
+		Volume *volume = volumes[i];
+		AttributeRequestSet& attributes = volume_attributes[i];
+
+		/* todo: we now store std and name attributes from requests even if
+		 * they actually refer to the same mesh attributes, optimize */
+		foreach(AttributeRequest& req, attributes.requests) {
+			Attribute *vattr = volume->attributes.find(req);
+
+			update_attribute_element_offset(vattr,
+			                                req.triangle_type,
+			                                req.triangle_desc.offset,
+			                                req.triangle_desc.element);
+
+			if(progress.get_cancel()) return;
+		}
+	}
+
+	update_svm_attributes(device, dscene, scene, volume_attributes);
+}
+
+void VolumeManager::update_svm_attributes(Device *device, DeviceScene *dscene, Scene *scene, vector<AttributeRequestSet>& mesh_attributes)
+{
+	/* compute array stride */
+	int attr_map_stride = 0;
+
+	for(size_t i = 0; i < volumes.size(); i++) {
+		attr_map_stride = max(attr_map_stride, (mesh_attributes[i].size() + 1));
+	}
+
+	if(attr_map_stride == 0) {
+		return;
+	}
+
+	/* create attribute map */
+	uint4 *attr_map = dscene->attributes_map.resize(attr_map_stride*volumes.size());
+	memset(attr_map, 0, dscene->attributes_map.size()*sizeof(uint));
+
+	for(size_t i = 0; i < volumes.size(); i++) {
+		AttributeRequestSet& attributes = mesh_attributes[i];
+
+		/* set object attributes */
+		int index = i*attr_map_stride;
+
+		foreach(AttributeRequest& req, attributes.requests) {
+			uint id = scene->shader_manager->get_attribute_id(req.name);
+
+			attr_map[index].x = id;
+			attr_map[index].y = req.triangle_desc.element;
+			attr_map[index].z = as_uint(req.triangle_desc.offset);
+
+			if(req.triangle_type == TypeDesc::TypeFloat)
+				attr_map[index].w = NODE_ATTR_FLOAT;
+			else
+				attr_map[index].w = NODE_ATTR_FLOAT3;
+
+			index++;
+		}
+
+		/* terminator */
+		attr_map[index].x = ATTR_STD_NONE;
+		attr_map[index].y = 0;
+		attr_map[index].z = 0;
+		attr_map[index].w = 0;
+
+		index++;
+	}
+
+	device->tex_alloc("__attributes_map", dscene->attributes_map);
+}
+
+void VolumeManager::device_update(Device *device, DeviceScene *dscene, Scene *scene, Progress& progress)
+{
+	if(!need_update) {
+		return;
+	}
+
+	device_free(device, dscene);
+	progress.set_status("Updating OpenVDB volumes", "Sending volumes to device.");
+
+	uint *vol_shader = dscene->vol_shader.resize(num_float_volume + num_float3_volume);
+	int s = 0;
+
+	for (size_t i = 0; i < volumes.size(); ++i) {
+		Volume *volume = volumes[i];
+
+		for(size_t i = 0; i < volume->scalar_grids.size(); ++i) {
+			if(!volume->scalar_grids[i]) {
+				continue;
+			}
+
+			vol_shader[s++] = scene->shader_manager->get_shader_id(volume->used_shaders[0], false);
+		}
+
+		for(size_t i = 0; i < volume->vector_grids.size(); ++i) {
+			if(!volume->vector_grids[i]) {
+				continue;
+			}
+
+			vol_shader[s++] = scene->shader_manager->get_shader_id(volume->used_shaders[0], false);
+		}
+
+		if(progress.get_cancel()) {
+			return;
+		}
+	}
+
+	device->tex_alloc("__vol_shader", dscene->vol_shader);
+
+#ifdef WITH_OPENVDB
+	typedef typename OpenVDBGlobals::scalar_grid_t scalar_grid_t;
+	typedef typename OpenVDBGlobals::vector_grid_t vector_grid_t;
+	typedef typename OpenVDBGlobals::scalar_isector_t scalar_isector_t;
+	typedef typename OpenVDBGlobals::vector_isector_t vector_isector_t;
+	OpenVDBGlobals *vdb = device->vdb_memory();
+	
+	vdb->scalar_grids.reserve(num_float_volume);
+	vdb->vector_grids.reserve(num_float3_volume);
+	vdb->scalar_main_isectors.reserve(num_float_volume);
+	vdb->vector_main_isectors.reserve(num_float3_volume);
+	for (size_t i = 0; i < volumes.size(); ++i) {
+		Volume *volume = volumes[i];
+		
+		for (size_t k = 0; k < volume->scalar_grids.size(); ++k) {
+			scalar_grid_t *grid = volume->scalar_grids[k].get();
+			vdb->scalar_grids.push_back(grid);
+			vdb->scalar_main_isectors.push_back(new scalar_isector_t(*grid));
+			VLOG(1) << grid->getName().c_str() << " memory usage: " << grid->memUsage() / 1024.0f << " kilobytes.\n";
+		}
+		for (size_t k = 0; k < volume->vector_grids.size(); ++k) {
+			vector_grid_t *grid = volume->vector_grids[k].get();
+			vdb->vector_grids.push_back(grid);
+			vdb->vector_main_isectors.push_back(new vector_isector_t(*grid));
+			VLOG(1) << grid->getName().c_str() << " memory usage: " << grid->memUsage() / 1024.0f << " kilobytes.\n";
+		}
+	}
+#endif
+
+	if(progress.get_cancel()) {
+		return;
+	}
+
+	dscene->data.tables.num_volumes = num_float_volume/* + float3_volumes.size()*/;
+	dscene->data.tables.density_index = 0;
+
+	need_update = false;
+}
+
+void VolumeManager::device_free(Device *device, DeviceScene *dscene)
+{
+#ifdef WITH_OPENVDB
+	OpenVDBGlobals *vdb = device->vdb_memory();
+	for (size_t i = 0; i < vdb->scalar_main_isectors.size(); ++i) {
+		delete vdb->scalar_main_isectors[i];
+	}
+	vdb->scalar_grids.clear();
+	vdb->scalar_main_isectors.clear();
+	
+	for (size_t i = 0; i < vdb->vector_main_isectors.size(); ++i) {
+		delete vdb->vector_main_isectors[i];
+	}
+	vdb->vector_grids.clear();
+	vdb->vector_main_isectors.clear();
+#endif
+	
+	device->tex_free(dscene->vol_shader);
+	dscene->vol_shader.clear();
+}
+
+void VolumeManager::tag_update(Scene */*scene*/)
+{
+	need_update = true;
+}
+
+CCL_NAMESPACE_END

intern/cycles/render/volume.h

@@ -0,0 +1,95 @@
+/*
+ * Copyright 2015 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 __VOLUMEMANAGER_H__
+#define __VOLUMEMANAGER_H__
+
+#include "attribute.h"
+
+#include "util/util_string.h"
+#include "util/util_types.h"
+
+#ifdef WITH_OPENVDB
+#include <openvdb/openvdb.h>
+#include <openvdb/tools/Interpolation.h>
+#include <openvdb/tools/RayIntersector.h>
+#endif
+
+CCL_NAMESPACE_BEGIN
+
+class Device;
+class DeviceScene;
+class Progress;
+class Scene;
+class Shader;
+
+class Volume {
+public:
+	vector<Shader*> used_shaders;
+	AttributeSet attributes;
+	string name;
+
+#ifdef WITH_OPENVDB
+	vector<openvdb::FloatGrid::Ptr> scalar_grids;
+	vector<openvdb::Vec3SGrid::Ptr> vector_grids;
+#endif
+	
+	void tag_update(Scene *scene, bool rebuild);
+};
+
+class VolumeManager {
+	struct GridDescription {
+		Volume *volume;
+		string filename;
+		string name;
+		int sampling;
+		int slot;
+	};
+
+	vector<GridDescription> current_grids;
+	int num_float_volume;
+	int num_float3_volume;
+
+	void delete_volume(int grid_type, int sampling, size_t slot);
+
+	void add_grid_description(Volume *volume, const string& filename, const string& name, int slot);
+	int find_existing_slot(Volume *volume, const string& filename, const string& name);
+
+	bool is_openvdb_file(const string& filename) const;
+	size_t add_openvdb_volume(Volume *volume, const string& filename, const string& name);
+
+public:
+	VolumeManager();
+	~VolumeManager();
+
+	int add_volume(Volume *volume, const string& filename, const string& name);
+	int find_density_slot();
+
+	void device_update(Device *device, DeviceScene *dscene, Scene *scene, Progress& progress);
+	void device_update_attributes(Device *device, DeviceScene *dscene, Scene *scene, Progress& progress);
+	void update_svm_attributes(Device *device, DeviceScene *dscene, Scene *scene, vector<AttributeRequestSet>& mesh_attributes);
+	void device_free(Device *device, DeviceScene *dscene);
+
+	void tag_update(Scene *scene);
+
+	bool need_update;
+
+	vector<Volume*> volumes;
+};
+
+CCL_NAMESPACE_END
+
+#endif /* __VOLUMEMANAGER_H__ */

intern/openvdb/openvdb_capi.cc

@@ -36,6 +36,33 @@ int OpenVDB_getVersionHex()
 	return openvdb::OPENVDB_LIBRARY_VERSION;
 }
 
+void OpenVDB_get_grid_info(const char *filename, OpenVDBGridInfoCallback cb, void *userdata)
+{
+	Timer(__func__);
+
+	using namespace openvdb;
+
+	initialize();
+
+	io::File file(filename);
+	file.open();
+
+	GridPtrVecPtr grids = file.getGrids();
+	int grid_num = grids->size();
+
+	for (size_t i = 0; i < grid_num; ++i) {
+		GridBase::ConstPtr grid = (*grids)[i];
+
+		Name name = grid->getName();
+		Name value_type = grid->valueType();
+		bool is_color = false;
+		if (grid->getMetadata< TypedMetadata<bool> >("is_color"))
+			is_color = grid->metaValue<bool>("is_color");
+
+		cb(userdata, name.c_str(), value_type.c_str(), is_color);
+	}
+}
+
 OpenVDBFloatGrid *OpenVDB_export_grid_fl(
         OpenVDBWriter *writer,
         const char *name, float *data,

intern/openvdb/openvdb_capi.h

@@ -38,6 +38,9 @@ struct OpenVDBVectorGrid;
 
 int OpenVDB_getVersionHex(void);
 
+typedef void (*OpenVDBGridInfoCallback)(void *userdata, const char *name, const char *value_type, bool is_color);
+void OpenVDB_get_grid_info(const char *filename, OpenVDBGridInfoCallback cb, void *userdata);
+
 enum {
 	VEC_INVARIANT = 0,
 	VEC_COVARIANT = 1,

source/blender/blenkernel/BKE_pointcache.h

@@ -280,6 +280,7 @@ void BKE_ptcache_ids_from_object(struct ListBase *lb, struct Object *ob, struct
 
 /***************** Global funcs ****************************/
 void BKE_ptcache_remove(void);
+int ptcache_filename(PTCacheID *pid, char *filename, int cfra, short do_path, short do_ext);
 
 /************ ID specific functions ************************/
 void    BKE_ptcache_id_clear(PTCacheID *id, int mode, unsigned int cfra);

source/blender/blenkernel/intern/pointcache.c

@@ -1814,7 +1814,7 @@ static int ptcache_path(PTCacheID *pid, char *filename)
 	return BLI_add_slash(filename); /* new strlen() */
 }
 
-static int ptcache_filename(PTCacheID *pid, char *filename, int cfra, short do_path, short do_ext)
+int ptcache_filename(PTCacheID *pid, char *filename, int cfra, short do_path, short do_ext)
 {
 	int len=0;
 	char *idname;

source/blender/makesdna/DNA_smoke_types.h

@@ -217,6 +217,8 @@ typedef struct SmokeDomainSettings {
 	char use_coba;
 	char coba_field;  /* simulation field used for the color mapping */
 	char pad2;
+
+	char cache_filename[1024];
 } SmokeDomainSettings;
 
 

source/blender/makesrna/intern/makesrna.c

@@ -3523,6 +3523,7 @@ static void rna_generate_header(BlenderRNA *UNUSED(brna), FILE *f)
 
 static const char *cpp_classes = ""
 "\n"
+"#include <cstdlib>\n"
 "#include <string>\n"
 "#include <string.h> /* for memcpy */\n"
 "\n"

source/blender/makesrna/intern/rna_smoke.c

@@ -396,6 +396,15 @@ static void rna_Smoke_use_color_ramp_set(PointerRNA *ptr, int value)
 	}
 }
 
+static void rna_SmokeModifier_cache_filename_get(PointerRNA *ptr, char *filename)
+{
+	SmokeDomainSettings *sds = (SmokeDomainSettings *)ptr->data;
+	PTCacheID pid;
+
+	BKE_ptcache_id_from_smoke(&pid, ptr->id.data, sds->smd);
+	ptcache_filename(&pid, filename, sds->smd->time, 1, 1);
+}
+
 #else
 
 static void rna_def_smoke_domain_settings(BlenderRNA *brna)
@@ -853,6 +862,12 @@ static void rna_def_smoke_domain_settings(BlenderRNA *brna)
 	RNA_def_property_struct_type(prop, "ColorRamp");
 	RNA_def_property_ui_text(prop, "Color Ramp", "");
 	RNA_def_property_update(prop, NC_OBJECT | ND_DRAW, NULL);
+
+	prop = RNA_def_property(srna, "cache_filename", PROP_STRING, PROP_NONE);
+	RNA_def_property_string_maxlength(prop, 1024);
+	RNA_def_property_string_funcs(prop, "rna_SmokeModifier_cache_filename_get", NULL, NULL);
+	RNA_def_property_clear_flag(prop, PROP_EDITABLE);
+	RNA_def_property_ui_text(prop, "Filename", "Path to the .blend file");
 }
 
 static void rna_def_smoke_flow_settings(BlenderRNA *brna)