diff --git a/.github/workflows/build-and-test.yml b/.github/workflows/build-and-test.yml
index 02946545..98d1ad94 100644
--- a/.github/workflows/build-and-test.yml
+++ b/.github/workflows/build-and-test.yml
@@ -12,6 +12,7 @@ jobs:
     strategy:
       matrix:
         build-type: [Debug, Release]
+        c-api: [ON, OFF]
         os: [ubuntu-latest, windows-latest, macos-latest]
 
     steps:
@@ -23,7 +24,7 @@ jobs:
     - name: Configure CMake
       working-directory: ${{runner.workspace}}/build
       shell: bash # Necessary because of the $GITHUB_WORKSPACE variable
-      run: cmake -DCMAKE_BUILD_TYPE=${{matrix.build-type}} -DBUILD_TESTING=ON -S $GITHUB_WORKSPACE
+      run: cmake -DCMAKE_BUILD_TYPE=${{matrix.build-type}} -DBVH_BUILD_C_API=${{matrix.c-api}} -DBUILD_TESTING=ON -S $GITHUB_WORKSPACE
 
     - name: Build
       working-directory: ${{runner.workspace}}/build
diff --git a/CMakeLists.txt b/CMakeLists.txt
index c5a1a770..1ab7c2e5 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -1,5 +1,12 @@
 cmake_minimum_required(VERSION 3.21)
-project(bvh)
+project(bvh VERSION 2.0)
+
+option(BVH_BUILD_C_API "Builds the C API library wrapper" OFF)
+option(BVH_STATIC_LINK_STDLIB_C_API "Link the C API library statically against the standard C++ library" OFF)
+
+set(CMAKE_RUNTIME_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/bin)
+set(CMAKE_LIBRARY_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
+set(CMAKE_ARCHIVE_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
 
 add_subdirectory(src/bvh/v2)
 
@@ -8,4 +15,6 @@ if (PROJECT_IS_TOP_LEVEL)
     if (BUILD_TESTING)
         add_subdirectory(test)
     endif()
+
+    include(cmake/Install.cmake)
 endif()
diff --git a/README.md b/README.md
index 1992bf36..bf6ec72d 100644
--- a/README.md
+++ b/README.md
@@ -46,7 +46,8 @@ Here is a list of features supported by this library (changes from `v1` are indi
 - [NEW] Variable amount of dimensions (e.g. 2D, 3D, 4D BVHs are supported) and different scalar types
   (e.g. `float` or `double`),
 - [NEW] Only depends on the standard library (parallelization uses a custom thread pool based on
-  `std::thread`).
+  `std::thread`),
+- [NEW] C API for the high-level parts of the library are available.
 
 ## Building
 
@@ -64,6 +65,14 @@ If you want to build the examples, use:
     cmake .. -DCMAKE_BUILD_TYPE=<Debug|Release> -DENABLE_TESTING=ON
     cmake --build .
 
+## C API
+
+The library can be used via a small set of high-level C bindings. These bindings are not enabled by
+default, but can be built by configuring CMake with `-DBVH_BUILD_C_API=ON`. Additionally, if the
+intent is to use the library in a pure C environment which does not have the C++ standard library as
+a dependency, it might be a good idea to statically the C++ standard library. That can be done by
+adding the flag `-DBVH_STATIC_LINK_STDLIB_C_API=ON` to the CMake command line.
+
 ## Usage
 
 The library contains several examples that are kept up-to-date with the API:
@@ -71,3 +80,4 @@ The library contains several examples that are kept up-to-date with the API:
 - A [basic example](test/simple_example.cpp) that traces one ray on a scene made of a couple of triangles,
 - A [benchmarking utility](test/benchmark.cpp) that showcases what the library can do.
 - A [serialization test](test/serialize.cpp) that shows how to save and load a BVH from a file.
+- A [C API example](test/c_api_example.cpp) that shows how to use the C bindings to this library.
diff --git a/cmake/Install.cmake b/cmake/Install.cmake
new file mode 100644
index 00000000..08df5ca0
--- /dev/null
+++ b/cmake/Install.cmake
@@ -0,0 +1,34 @@
+set(bvh_targets bvh)
+if (BVH_BUILD_C_API)
+    list(APPEND bvh_targets bvh_c)
+endif()
+
+install(
+    TARGETS ${bvh_targets}
+    EXPORT bvh_exports
+    RUNTIME DESTINATION bin/
+    LIBRARY DESTINATION lib/
+    ARCHIVE DESTINATION lib/
+    INCLUDES DESTINATION include/)
+install(
+    EXPORT bvh_exports
+    FILE bvh-targets.cmake
+    NAMESPACE bvh::v2::
+    DESTINATION lib/cmake/bvh/v2/)
+
+include(CMakePackageConfigHelpers)
+
+configure_package_config_file(
+    "${PROJECT_SOURCE_DIR}/cmake/bvh-config.cmake.in"
+    "${CMAKE_CURRENT_BINARY_DIR}/bvh-config.cmake"
+    INSTALL_DESTINATION lib/cmake/bvh/v2/)
+
+write_basic_package_version_file(
+    "${CMAKE_CURRENT_BINARY_DIR}/bvh-config-version.cmake"
+    COMPATIBILITY AnyNewerVersion)
+
+install(
+    FILES
+        "${CMAKE_CURRENT_BINARY_DIR}/bvh-config.cmake"
+        "${CMAKE_CURRENT_BINARY_DIR}/bvh-config-version.cmake"
+    DESTINATION lib/cmake/bvh/v2/)
diff --git a/cmake/bvh-config.cmake.in b/cmake/bvh-config.cmake.in
new file mode 100644
index 00000000..7aa073b0
--- /dev/null
+++ b/cmake/bvh-config.cmake.in
@@ -0,0 +1,5 @@
+@PACKAGE_INIT@
+
+include("${CMAKE_CURRENT_LIST_DIR}/bvh-targets.cmake")
+
+check_required_components(bvh_c)
diff --git a/src/bvh/v2/CMakeLists.txt b/src/bvh/v2/CMakeLists.txt
index 1cbb27b0..8bbfd345 100644
--- a/src/bvh/v2/CMakeLists.txt
+++ b/src/bvh/v2/CMakeLists.txt
@@ -6,5 +6,18 @@ if (Threads_FOUND)
     target_link_libraries(bvh INTERFACE Threads::Threads)
 endif()
 
-target_include_directories(bvh INTERFACE ../..)
+target_include_directories(bvh INTERFACE
+    $<BUILD_INTERFACE:${PROJECT_SOURCE_DIR}/src>
+    $<INSTALL_INTERFACE:include>)
+
 set_target_properties(bvh PROPERTIES CXX_STANDARD 20)
+
+install(
+    DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR}/
+    DESTINATION include/bvh/v2
+    FILES_MATCHING PATTERN "*.h"
+    PATTERN "c_api" EXCLUDE)
+
+if (BVH_BUILD_C_API)
+    add_subdirectory(c_api)
+endif()
diff --git a/src/bvh/v2/c_api/CMakeLists.txt b/src/bvh/v2/c_api/CMakeLists.txt
new file mode 100644
index 00000000..b2c0b309
--- /dev/null
+++ b/src/bvh/v2/c_api/CMakeLists.txt
@@ -0,0 +1,17 @@
+add_library(bvh_c SHARED bvh.cpp)
+target_link_libraries(bvh_c PRIVATE bvh)
+target_compile_definitions(bvh_c PRIVATE -DBVH_BUILD_API)
+target_include_directories(bvh_c INTERFACE
+    $<BUILD_INTERFACE:${PROJECT_SOURCE_DIR}/src>
+    $<INSTALL_INTERFACE:include>)
+set_target_properties(bvh_c PROPERTIES
+    CXX_STANDARD 20
+    CXX_VISIBILITY_PRESET hidden
+    INTERPROCEDURAL_OPTIMIZATION_RELEASE ON)
+
+if (BVH_STATIC_LINK_STDLIB_C_API)
+    # Link statically against standard C++ library
+    target_link_options(bvh_c PRIVATE $<$<C_COMPILER_ID:GNU,Clang>:-static-libstdc++>)
+endif()
+
+install(FILES bvh.h DESTINATION include/bvh/v2/c_api/)
diff --git a/src/bvh/v2/c_api/bvh.cpp b/src/bvh/v2/c_api/bvh.cpp
new file mode 100644
index 00000000..df48f93e
--- /dev/null
+++ b/src/bvh/v2/c_api/bvh.cpp
@@ -0,0 +1,28 @@
+#include <bvh/v2/c_api/bvh.h>
+#include <bvh/v2/c_api/bvh_impl.h>
+#include <bvh/v2/thread_pool.h>
+
+namespace bvh::v2::c_api {
+
+BVH_TYPES(float, 2, 2f)
+BVH_TYPES(float, 3, 3f)
+BVH_TYPES(double, 2, 2d)
+BVH_TYPES(double, 3, 3d)
+
+extern "C" {
+
+BVH_EXPORT struct bvh_thread_pool* bvh_thread_pool_create(size_t thread_count) {
+    return reinterpret_cast<bvh_thread_pool*>(new bvh::v2::ThreadPool(thread_count));
+}
+
+BVH_EXPORT void bvh_thread_pool_destroy(bvh_thread_pool* thread_pool) {
+    return delete reinterpret_cast<bvh::v2::ThreadPool*>(thread_pool);
+}
+
+BVH_IMPL(float, 2, 2f)
+BVH_IMPL(float, 3, 3f)
+BVH_IMPL(double, 2, 2d)
+BVH_IMPL(double, 3, 3d)
+
+} // extern "C"
+} // namespace bvh::v2::c_api
diff --git a/src/bvh/v2/c_api/bvh.h b/src/bvh/v2/c_api/bvh.h
new file mode 100644
index 00000000..819b0f9d
--- /dev/null
+++ b/src/bvh/v2/c_api/bvh.h
@@ -0,0 +1,301 @@
+#ifndef BVH_V2_C_API_BVH_H
+#define BVH_V2_C_API_BVH_H
+
+// C API to the BVH library, providing high-level access to BVH construction and traversal. This API
+// may have an overhead compared to using the C++ API directly, as it may have to translate types
+// between the C interface and the C++ one, and because callbacks are not going to be inlined at the
+// API boundary. To mitigate that problem, specialized functions are provided which provide faster
+// operation for a small set of situations.
+
+#include <stddef.h>
+#include <stdio.h>
+#include <stdbool.h>
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#ifdef _MSC_VER
+#define BVH_EXPORT __declspec(dllexport)
+#define BVH_IMPORT __declspec(dllimport)
+#else
+#define BVH_EXPORT __attribute__((visibility("default")))
+#define BVH_IMPORT BVH_EXPORT
+#endif
+
+#ifdef BVH_BUILD_API
+#define BVH_API BVH_EXPORT
+#else
+#define BVH_API BVH_IMPORT
+#endif
+
+#define BVH_ROOT_INDEX 0
+#define BVH_INVALID_PRIM_ID SIZE_MAX
+
+struct bvh2f;
+struct bvh3f;
+struct bvh_node2f;
+struct bvh_node3f;
+
+struct bvh2d;
+struct bvh3d;
+struct bvh_node2d;
+struct bvh_node3d;
+
+struct bvh_thread_pool;
+
+enum bvh_build_quality {
+    BVH_BUILD_QUALITY_LOW,
+    BVH_BUILD_QUALITY_MEDIUM,
+    BVH_BUILD_QUALITY_HIGH
+};
+
+struct bvh_build_config {
+    enum bvh_build_quality quality;
+    size_t min_leaf_size;
+    size_t max_leaf_size;
+    size_t parallel_threshold;
+};
+
+struct bvh_vec2f { float x, y; };
+struct bvh_vec3f { float x, y, z; };
+struct bvh_vec2d { double x, y; };
+struct bvh_vec3d { double x, y, z; };
+
+struct bvh_bbox2f { struct bvh_vec2f min, max; };
+struct bvh_bbox3f { struct bvh_vec3f min, max; };
+struct bvh_bbox2d { struct bvh_vec2d min, max; };
+struct bvh_bbox3d { struct bvh_vec3d min, max; };
+
+struct bvh_ray2f { struct bvh_vec2f org, dir; float tmin, tmax; };
+struct bvh_ray3f { struct bvh_vec3f org, dir; float tmin, tmax; };
+struct bvh_ray2d { struct bvh_vec2d org, dir; double tmin, tmax; };
+struct bvh_ray3d { struct bvh_vec3d org, dir; double tmin, tmax; };
+
+struct bvh_intersect_callbackf {
+    void* user_data;
+    bool (*user_fn)(void*, float*, size_t begin, size_t end);
+};
+
+struct bvh_intersect_callbackd {
+    void* user_data;
+    bool (*user_fn)(void*, double*, size_t begin, size_t end);
+};
+
+// Thread Pool ------------------------------------------------------------------------------------
+
+// A thread count of zero instructs the thread pool to detect the number of threads available on the
+// machine.
+
+BVH_API struct bvh_thread_pool* bvh_thread_pool_create(size_t thread_count);
+BVH_API void bvh_thread_pool_destroy(struct bvh_thread_pool*);
+
+// BVH Construction -------------------------------------------------------------------------------
+
+// These construction functions can be called with a `NULL` thread pool, in which case the BVH is
+// constructed serially. The configuration can also be `NULL`, in which case the default
+// configuration is used.
+
+BVH_API struct bvh2f* bvh2f_build(
+    struct bvh_thread_pool* thread_pool,
+    const struct bvh_bbox2f* bboxes,
+    const struct bvh_vec2f* centers,
+    size_t prim_count,
+    const struct bvh_build_config* config);
+
+BVH_API struct bvh3f* bvh3f_build(
+    struct bvh_thread_pool* thread_pool,
+    const struct bvh_bbox3f* bboxes,
+    const struct bvh_vec3f* centers,
+    size_t prim_count,
+    const struct bvh_build_config* config);
+
+BVH_API struct bvh2d* bvh2d_build(
+    struct bvh_thread_pool* thread_pool,
+    const struct bvh_bbox2d* bboxes,
+    const struct bvh_vec2d* centers,
+    size_t prim_count,
+    const struct bvh_build_config* config);
+
+BVH_API struct bvh3d* bvh3d_build(
+    struct bvh_thread_pool* thread_pool,
+    const struct bvh_bbox3d* bboxes,
+    const struct bvh_vec3d* centers,
+    size_t prim_count,
+    const struct bvh_build_config* config);
+
+// BVH Destruction --------------------------------------------------------------------------------
+
+BVH_API void bvh2f_destroy(struct bvh2f*);
+BVH_API void bvh3f_destroy(struct bvh3f*);
+BVH_API void bvh2d_destroy(struct bvh2d*);
+BVH_API void bvh3d_destroy(struct bvh3d*);
+
+// Serialization/Deserialization to Files ---------------------------------------------------------
+
+BVH_API void bvh2f_save(const struct bvh2f*, FILE*);
+BVH_API void bvh3f_save(const struct bvh3f*, FILE*);
+BVH_API void bvh2d_save(const struct bvh2d*, FILE*);
+BVH_API void bvh3d_save(const struct bvh3d*, FILE*);
+
+BVH_API struct bvh2f* bvh2f_load(FILE*);
+BVH_API struct bvh3f* bvh3f_load(FILE*);
+BVH_API struct bvh2d* bvh2d_load(FILE*);
+BVH_API struct bvh3d* bvh3d_load(FILE*);
+
+// Accessing Nodes and Primitive Indices ----------------------------------------------------------
+
+BVH_API struct bvh_node2f* bvh2f_get_node(struct bvh2f*, size_t);
+BVH_API struct bvh_node3f* bvh3f_get_node(struct bvh3f*, size_t);
+BVH_API struct bvh_node2d* bvh2d_get_node(struct bvh2d*, size_t);
+BVH_API struct bvh_node3d* bvh3d_get_node(struct bvh3d*, size_t);
+
+BVH_API size_t bvh2f_get_prim_id(const struct bvh2f*, size_t);
+BVH_API size_t bvh3f_get_prim_id(const struct bvh3f*, size_t);
+BVH_API size_t bvh2d_get_prim_id(const struct bvh2d*, size_t);
+BVH_API size_t bvh3d_get_prim_id(const struct bvh3d*, size_t);
+
+BVH_API size_t bvh2f_get_prim_count(const struct bvh2f*);
+BVH_API size_t bvh3f_get_prim_count(const struct bvh3f*);
+BVH_API size_t bvh2d_get_prim_count(const struct bvh2d*);
+BVH_API size_t bvh3d_get_prim_count(const struct bvh3d*);
+
+BVH_API size_t bvh2f_get_node_count(const struct bvh2f*);
+BVH_API size_t bvh3f_get_node_count(const struct bvh3f*);
+BVH_API size_t bvh2d_get_node_count(const struct bvh2d*);
+BVH_API size_t bvh3d_get_node_count(const struct bvh3d*);
+
+// Accessing and Modifying Node properties --------------------------------------------------------
+
+BVH_API bool bvh_node2f_is_leaf(const struct bvh_node2f*);
+BVH_API bool bvh_node3f_is_leaf(const struct bvh_node3f*);
+BVH_API bool bvh_node2d_is_leaf(const struct bvh_node2d*);
+BVH_API bool bvh_node3d_is_leaf(const struct bvh_node3d*);
+
+BVH_API size_t bvh_node2f_get_prim_count(const struct bvh_node2f*);
+BVH_API size_t bvh_node3f_get_prim_count(const struct bvh_node3f*);
+BVH_API size_t bvh_node2d_get_prim_count(const struct bvh_node2d*);
+BVH_API size_t bvh_node3d_get_prim_count(const struct bvh_node3d*);
+
+BVH_API void bvh_node2f_set_prim_count(struct bvh_node2f*, size_t);
+BVH_API void bvh_node3f_set_prim_count(struct bvh_node3f*, size_t);
+BVH_API void bvh_node2d_set_prim_count(struct bvh_node2d*, size_t);
+BVH_API void bvh_node3d_set_prim_count(struct bvh_node3d*, size_t);
+
+BVH_API size_t bvh_node2f_get_first_id(const struct bvh_node2f*);
+BVH_API size_t bvh_node3f_get_first_id(const struct bvh_node3f*);
+BVH_API size_t bvh_node2d_get_first_id(const struct bvh_node2d*);
+BVH_API size_t bvh_node3d_get_first_id(const struct bvh_node3d*);
+
+BVH_API void bvh_node2f_set_first_id(struct bvh_node2f*, size_t);
+BVH_API void bvh_node3f_set_first_id(struct bvh_node3f*, size_t);
+BVH_API void bvh_node2d_set_first_id(struct bvh_node2d*, size_t);
+BVH_API void bvh_node3d_set_first_id(struct bvh_node3d*, size_t);
+
+BVH_API struct bvh_bbox2f bvh_node2f_get_bbox(const struct bvh_node2f*);
+BVH_API struct bvh_bbox3f bvh_node3f_get_bbox(const struct bvh_node3f*);
+BVH_API struct bvh_bbox2d bvh_node2d_get_bbox(const struct bvh_node2d*);
+BVH_API struct bvh_bbox3d bvh_node3d_get_bbox(const struct bvh_node3d*);
+
+BVH_API void bvh_node2f_set_bbox(struct bvh_node2f*, const struct bvh_bbox2f*);
+BVH_API void bvh_node3f_set_bbox(struct bvh_node3f*, const struct bvh_bbox3f*);
+BVH_API void bvh_node2d_set_bbox(struct bvh_node2d*, const struct bvh_bbox2d*);
+BVH_API void bvh_node3d_set_bbox(struct bvh_node3d*, const struct bvh_bbox3d*);
+
+// Refitting and BVH Modification -----------------------------------------------------------------
+
+// Refitting functions resize the bounding boxes of the BVH in a bottom-up fashion. This can be
+// combined with the optimization function below to allow updating the BVH in an incremental manner.
+// IMPORTANT: Appending a node to the BVH invalidates all the node pointers.
+
+BVH_API void bvh2f_append_node(struct bvh2f*);
+BVH_API void bvh3f_append_node(struct bvh3f*);
+BVH_API void bvh2d_append_node(struct bvh2d*);
+BVH_API void bvh3d_append_node(struct bvh3d*);
+
+BVH_API void bvh2f_remove_last_node(struct bvh2f*);
+BVH_API void bvh3f_remove_last_node(struct bvh3f*);
+BVH_API void bvh2d_remove_last_node(struct bvh2d*);
+BVH_API void bvh3d_remove_last_node(struct bvh3d*);
+
+BVH_API void bvh2f_refit(struct bvh2f*);
+BVH_API void bvh3f_refit(struct bvh3f*);
+BVH_API void bvh2d_refit(struct bvh2d*);
+BVH_API void bvh3d_refit(struct bvh3d*);
+
+BVH_API void bvh2f_optimize(struct bvh_thread_pool*, struct bvh2f*);
+BVH_API void bvh3f_optimize(struct bvh_thread_pool*, struct bvh3f*);
+BVH_API void bvh2d_optimize(struct bvh_thread_pool*, struct bvh2d*);
+BVH_API void bvh3d_optimize(struct bvh_thread_pool*, struct bvh3d*);
+
+// BVH Intersection -------------------------------------------------------------------------------
+
+// Intersection routines: Intersects the BVH with a ray using a callback to intersect the primitives
+// contained in the leaves.
+//
+// The callback takes a pointer to user data, a pointer to the current distance along the ray, and a
+// range of primitives. It returns true when an intersection is found, in which case it writes the
+// intersection distance in the provided pointer, otherwise, if no intersection is found, it returns
+// false and leaves the intersection distance unchanged. The given range is expressed in terms of
+// *BVH primitives*, which means that it does not correspond to a range in the original set of
+// primitives. Instead, the user has two options: pre-permute the primitives according to the BVH
+// primitive indices, or perform an indirection everytime a primitive is accessed.
+//
+// Here is a basic example of an intersection callback:
+//
+//     struct my_user_data {
+//         struct bvh3f* bvh;
+//         struct bvh_ray3f* ray;
+//         struct my_prim* prims;
+//     };
+//
+//     struct my_prim_hit {
+//         float t;
+//     };
+//
+//     // Declared & defined somewhere else.
+//     bool intersect_prim(const struct my_prim*, const struct bvh_ray3f*, struct my_prim_hit*);
+//
+//     bool my_user_fn(void* user_data, float* t, size_t begin, size_t end) {
+//         struct my_user_data* my_user_data = (my_user_data)data;
+//         bool was_hit = false;
+//         for (size_t i = begin; i < end; ++i) {
+//             const size_t prim_id = bvh3f_get_prim_id(my_user_data->bvh, i);
+//             struct my_prim_hit hit;
+//             if (intersect_prim(&my_user_data->prims[i], my_user_data->ray, &hit)) {
+//                 // Note: It is important to remember the maximum distance so that the BVH
+//                 // traversal routine can cull nodes that are too far away, and so that the primitive
+//                 // intersection routine can exit earlier when that is possible.
+//                 *t = my_user_data->ray->tmax = hit.t;
+//                 was_hit = true;
+//             }
+//         }
+//         return was_hit;
+//     }
+//
+
+BVH_API void bvh2f_intersect_ray_any(const struct bvh2f*, const struct bvh_ray2f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh3f_intersect_ray_any(const struct bvh3f*, const struct bvh_ray3f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh2d_intersect_ray_any(const struct bvh2d*, const struct bvh_ray2d*, const struct bvh_intersect_callbackd*);
+BVH_API void bvh3d_intersect_ray_any(const struct bvh3d*, const struct bvh_ray3d*, const struct bvh_intersect_callbackd*);
+
+BVH_API void bvh2f_intersect_ray_any_robust(const struct bvh2f*, const struct bvh_ray2f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh3f_intersect_ray_any_robust(const struct bvh3f*, const struct bvh_ray3f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh2d_intersect_ray_any_robust(const struct bvh2d*, const struct bvh_ray2d*, const struct bvh_intersect_callbackd*);
+BVH_API void bvh3d_intersect_ray_any_robust(const struct bvh3d*, const struct bvh_ray3d*, const struct bvh_intersect_callbackd*);
+
+BVH_API void bvh2f_intersect_ray(const struct bvh2f*, const struct bvh_ray2f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh3f_intersect_ray(const struct bvh3f*, const struct bvh_ray3f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh2d_intersect_ray(const struct bvh2d*, const struct bvh_ray2d*, const struct bvh_intersect_callbackd*);
+BVH_API void bvh3d_intersect_ray(const struct bvh3d*, const struct bvh_ray3d*, const struct bvh_intersect_callbackd*);
+
+BVH_API void bvh2f_intersect_ray_robust(const struct bvh2f*, const struct bvh_ray2f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh3f_intersect_ray_robust(const struct bvh3f*, const struct bvh_ray3f*, const struct bvh_intersect_callbackf*);
+BVH_API void bvh2d_intersect_ray_robust(const struct bvh2d*, const struct bvh_ray2d*, const struct bvh_intersect_callbackd*);
+BVH_API void bvh3d_intersect_ray_robust(const struct bvh3d*, const struct bvh_ray3d*, const struct bvh_intersect_callbackd*);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif
diff --git a/src/bvh/v2/c_api/bvh_impl.h b/src/bvh/v2/c_api/bvh_impl.h
new file mode 100644
index 00000000..cae5b3c2
--- /dev/null
+++ b/src/bvh/v2/c_api/bvh_impl.h
@@ -0,0 +1,357 @@
+#ifndef BVH_V2_C_API_BVH_IMPL_H
+#define BVH_V2_C_API_BVH_IMPL_H
+
+#include <bvh/v2/stack.h>
+#include <bvh/v2/default_builder.h>
+#include <bvh/v2/reinsertion_optimizer.h>
+
+namespace bvh::v2::c_api {
+
+template <typename T, size_t Dim>
+struct BvhTypes {};
+
+template <typename T>
+struct BvhCallback {};
+
+template <> struct BvhCallback<float>  { using Type = bvh_intersect_callbackf; };
+template <> struct BvhCallback<double> { using Type = bvh_intersect_callbackd; };
+
+template <typename T, size_t Dim>
+static auto translate(enum bvh_build_quality quality) {
+    switch (quality) {
+#ifndef BVH_C_UNSAFE_CASTS
+        case BVH_BUILD_QUALITY_LOW:
+            return bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::Quality::Low;
+        case BVH_BUILD_QUALITY_MEDIUM:
+            return bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::Quality::Medium;
+        case BVH_BUILD_QUALITY_HIGH:
+            return bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::Quality::High;
+#endif
+        default:
+            return static_cast<bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::Quality>(quality);
+    }
+}
+
+template <typename T, size_t Dim>
+static auto translate(const bvh_build_config* config) {
+    auto translated_config = typename bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::Config {};
+    if (config) {
+        translated_config.quality = translate<T, Dim>(config->quality);
+        translated_config.min_leaf_size = config->min_leaf_size;
+        translated_config.max_leaf_size = config->max_leaf_size;
+        translated_config.parallel_threshold = config->parallel_threshold;
+    }
+    return translated_config;
+}
+
+template <typename T, size_t Dim>
+static auto translate(const typename BvhTypes<T, Dim>::Vec& vec) {
+    static_assert(Dim == 2 || Dim == 3);
+    if constexpr (Dim == 2) {
+        return bvh::v2::Vec<T, 2> { vec.x, vec.y };
+    } else {
+        return bvh::v2::Vec<T, Dim> { vec.x, vec.y, vec.z };
+    }
+}
+
+template <typename T, size_t Dim>
+static auto translate(const bvh::v2::Vec<T, Dim>& vec) {
+    static_assert(Dim == 2 || Dim == 3);
+    if constexpr (Dim == 2) {
+        return typename BvhTypes<T, Dim>::Vec { vec[0], vec[1] };
+    } else {
+        return typename BvhTypes<T, Dim>::Vec { vec[0], vec[1], vec[2] };
+    }
+}
+
+template <typename T, size_t Dim>
+static auto translate(const typename BvhTypes<T, Dim>::BBox& bbox) {
+    return bvh::v2::BBox<T, Dim> { translate<T, Dim>(bbox.min), translate<T, Dim>(bbox.max) };
+}
+
+template <typename T, size_t Dim>
+static auto translate(const bvh::v2::BBox<T, Dim>& bbox) {
+    return typename BvhTypes<T, Dim>::BBox { translate<T, Dim>(bbox.min), translate<T, Dim>(bbox.max) };
+}
+
+template <typename T, size_t Dim>
+static auto translate(const typename BvhTypes<T, Dim>::Ray& ray) {
+    return bvh::v2::Ray<T, Dim> { translate<T, Dim>(ray.org), translate<T, Dim>(ray.dir), ray.tmin, ray.tmax };
+}
+
+template <typename T, size_t Dim>
+static typename BvhTypes<T, Dim>::Bvh* bvh_build(
+    bvh_thread_pool* thread_pool,
+    const typename BvhTypes<T, Dim>::BBox* bboxes,
+    const typename BvhTypes<T, Dim>::Vec* centers,
+    size_t prim_count,
+    const bvh_build_config* config)
+{
+    bvh::v2::BBox<T, Dim>* translated_bboxes = nullptr;
+    bvh::v2::Vec<T, Dim>* translated_centers = nullptr;
+#ifdef BVH_C_UNSAFE_CASTS
+    translated_bboxes  = reinterpret_cast<decltype(translated_bboxes)>(bboxes);
+    translated_centers = reinterpret_cast<decltype(translated_centers)>(centers);
+#else
+    std::vector<bvh::v2::BBox<T, Dim>> bbox_buf(prim_count);
+    std::vector<bvh::v2::Vec<T, Dim>> center_buf(prim_count);
+    for (size_t i = 0; i < prim_count; ++i) {
+        bbox_buf[i]   = translate<T, Dim>(bboxes[i]);
+        center_buf[i] = translate<T, Dim>(centers[i]);
+    }
+    translated_bboxes = bbox_buf.data();
+    translated_centers = center_buf.data();
+#endif
+    auto bvh = thread_pool
+        ? bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::build(
+            *reinterpret_cast<bvh::v2::ThreadPool*>(thread_pool),
+            std::span { translated_bboxes,  translated_bboxes + prim_count },
+            std::span { translated_centers, translated_centers + prim_count },
+            translate<T, Dim>(config))
+        : bvh::v2::DefaultBuilder<bvh::v2::Node<T, Dim>>::build(
+            std::span { translated_bboxes,  translated_bboxes + prim_count },
+            std::span { translated_centers, translated_centers + prim_count },
+            translate<T, Dim>(config));
+    return reinterpret_cast<typename BvhTypes<T, Dim>::Bvh*>(new decltype(bvh) { std::move(bvh) });
+}
+
+template <typename T, size_t Dim>
+static void bvh_save(const typename BvhTypes<T, Dim>::Bvh* bvh, FILE* file) {
+    struct FileOutputStream : bvh::v2::OutputStream {
+        FILE* file;
+        FileOutputStream(FILE* file)
+            : file(file)
+        {}
+        bool write_raw(const void* data, size_t size) override {
+            return fwrite(data, 1, size, file) == size;
+        }
+    };
+    FileOutputStream stream { file };
+    reinterpret_cast<const bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->serialize(stream);
+}
+
+template <typename T, size_t Dim>
+static typename BvhTypes<T, Dim>::Bvh* bvh_load(FILE* file) {
+    struct FileInputStream : bvh::v2::InputStream {
+        FILE* file;
+        FileInputStream(FILE* file)
+            : file(file)
+        {}
+        size_t read_raw(void* data, size_t size) override {
+            return fread(data, 1, size, file);
+        }
+    };
+    FileInputStream stream { file };
+    return reinterpret_cast<typename BvhTypes<T, Dim>::Bvh*>(
+        new bvh::v2::Bvh<bvh::v2::Node<T, Dim>>(bvh::v2::Bvh<bvh::v2::Node<T, Dim>>::deserialize(stream)));
+}
+
+template <typename T, size_t Dim>
+static typename BvhTypes<T, Dim>::Node* bvh_get_node(typename BvhTypes<T, Dim>::Bvh* bvh, size_t node_id) {
+    auto& nodes = reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->nodes;
+    assert(node_id < nodes.size());
+    return reinterpret_cast<typename BvhTypes<T, Dim>::Node*>(&nodes[node_id]);
+}
+
+template <typename T, size_t Dim>
+static size_t bvh_get_prim_id(const typename BvhTypes<T, Dim>::Bvh* bvh, size_t i) {
+    auto& prim_ids = reinterpret_cast<const bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->prim_ids;
+    assert(i < prim_ids.size());
+    return prim_ids[i];
+}
+
+template <typename T, size_t Dim>
+static size_t bvh_get_prim_count(const typename BvhTypes<T, Dim>::Bvh* bvh) {
+    return reinterpret_cast<const bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->prim_ids.size();
+}
+
+template <typename T, size_t Dim>
+static size_t bvh_get_node_count(const typename BvhTypes<T, Dim>::Bvh* bvh) {
+    return reinterpret_cast<const bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->nodes.size();
+}
+
+template <typename T, size_t Dim>
+static bool bvh_node_is_leaf(const typename BvhTypes<T, Dim>::Node* node) {
+    return reinterpret_cast<const bvh::v2::Node<T, Dim>*>(node)->is_leaf();
+}
+
+template <typename T, size_t Dim>
+static size_t bvh_node_get_prim_count(const typename BvhTypes<T, Dim>::Node* node) {
+    return reinterpret_cast<const bvh::v2::Node<T, Dim>*>(node)->index.prim_count();
+}
+
+template <typename T, size_t Dim>
+static void bvh_node_set_prim_count(typename BvhTypes<T, Dim>::Node* node, size_t prim_count) {
+    return reinterpret_cast<bvh::v2::Node<T, Dim>*>(node)->index.set_prim_count(prim_count);
+}
+
+template <typename T, size_t Dim>
+static size_t bvh_node_get_first_id(const typename BvhTypes<T, Dim>::Node* node) {
+    return reinterpret_cast<const bvh::v2::Node<T, Dim>*>(node)->index.first_id();
+}
+
+template <typename T, size_t Dim>
+static void bvh_node_set_first_id(typename BvhTypes<T, Dim>::Node* node, size_t first_id) {
+    reinterpret_cast<bvh::v2::Node<T, Dim>*>(node)->index.set_first_id(first_id);
+}
+
+template <typename T, size_t Dim>
+static typename BvhTypes<T, Dim>::BBox bvh_node_get_bbox(const typename BvhTypes<T, Dim>::Node* node) {
+    return translate(reinterpret_cast<const bvh::v2::Node<T, Dim>*>(node)->get_bbox());
+}
+
+template <typename T, size_t Dim>
+static void bvh_node_set_bbox(typename BvhTypes<T, Dim>::Node* node, const typename BvhTypes<T, Dim>::BBox* bbox) {
+    reinterpret_cast<bvh::v2::Node<T, Dim>*>(node)->set_bbox(translate<T, Dim>(*bbox));
+}
+
+template <typename T, size_t Dim>
+static void bvh_append_node(typename BvhTypes<T, Dim>::Bvh* bvh) {
+    reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->nodes.emplace_back();
+}
+
+template <typename T, size_t Dim>
+static void bvh_remove_last_node(typename BvhTypes<T, Dim>::Bvh* bvh) {
+    reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->nodes.pop_back();
+}
+
+template <typename T, size_t Dim>
+static void bvh_refit(typename BvhTypes<T, Dim>::Bvh* bvh) {
+    reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh)->refit();
+}
+
+template <typename T, size_t Dim>
+static void bvh_optimize(bvh_thread_pool* thread_pool, typename BvhTypes<T, Dim>::Bvh* bvh) {
+    if (thread_pool) {
+        bvh::v2::ReinsertionOptimizer<bvh::v2::Node<T, Dim>>::optimize(
+            *reinterpret_cast<bvh::v2::ThreadPool*>(thread_pool),
+            *reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh));
+    } else {
+        bvh::v2::ReinsertionOptimizer<bvh::v2::Node<T, Dim>>::optimize(
+            *reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh));
+    }
+}
+
+template <typename T, size_t Dim, bool AnyHit, bool RobustTraversal>
+static void bvh_intersect_ray(
+    const typename BvhTypes<T, Dim>::Bvh* bvh,
+    const typename BvhTypes<T, Dim>::Ray* ray,
+    const typename BvhCallback<T>::Type* callback)
+{
+    static constexpr size_t stack_size = 64;
+    auto translated_ray = translate<T, Dim>(*ray);
+    auto& translated_bvh = *reinterpret_cast<const bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh);
+    bvh::v2::SmallStack<typename bvh::v2::Node<T, Dim>::Index, stack_size> stack;
+    translated_bvh.template intersect<AnyHit, RobustTraversal>(
+        translated_ray, translated_bvh.get_root().index, stack,
+        [&] (size_t begin, size_t end) {
+            return callback->user_fn(callback->user_data, &translated_ray.tmax, begin, end);
+        });
+}
+
+#define BVH_TYPES(T, Dim, suffix) \
+    template <> \
+    struct BvhTypes<T, Dim> { \
+        using Bvh  = bvh##suffix; \
+        using Vec  = bvh_vec##suffix; \
+        using BBox = bvh_bbox##suffix; \
+        using Node = bvh_node##suffix; \
+        using Ray  = bvh_ray##suffix; \
+    };
+
+#define BVH_IMPL(T, Dim, suffix) \
+    BVH_EXPORT typename BvhTypes<T, Dim>::Bvh* bvh##suffix##_build( \
+        bvh_thread_pool* thread_pool, \
+        const typename BvhTypes<T, Dim>::BBox* bboxes, \
+        const typename BvhTypes<T, Dim>::Vec* centers, \
+        size_t prim_count, \
+        const bvh_build_config* config) \
+    { \
+        return bvh_build<T, Dim>(thread_pool, bboxes, centers, prim_count, config); \
+    } \
+    BVH_EXPORT void bvh##suffix##_destroy(typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        delete reinterpret_cast<bvh::v2::Bvh<bvh::v2::Node<T, Dim>>*>(bvh); \
+    } \
+    BVH_EXPORT void bvh##suffix##_save(const typename BvhTypes<T, Dim>::Bvh* bvh, FILE* file) { \
+        bvh_save<T, Dim>(bvh, file); \
+    } \
+    BVH_EXPORT typename BvhTypes<T, Dim>::Bvh* bvh##suffix##_load(FILE* file) { \
+        return bvh_load<T, Dim>(file); \
+    } \
+    BVH_EXPORT typename BvhTypes<T, Dim>::Node* bvh##suffix##_get_node(typename BvhTypes<T, Dim>::Bvh* bvh, size_t node_id) { \
+        return bvh_get_node<T, Dim>(bvh, node_id); \
+    } \
+    BVH_EXPORT size_t bvh##suffix##_get_prim_id(const typename BvhTypes<T, Dim>::Bvh* bvh, size_t i) { \
+        return bvh_get_prim_id<T, Dim>(bvh, i); \
+    } \
+    BVH_EXPORT size_t bvh##suffix##_get_prim_count(const typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        return bvh_get_prim_count<T, Dim>(bvh); \
+    } \
+    BVH_EXPORT size_t bvh##suffix##_get_node_count(const typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        return bvh_get_node_count<T, Dim>(bvh); \
+    } \
+    BVH_EXPORT bool bvh_node##suffix##_is_leaf(const typename BvhTypes<T, Dim>::Node* node) { \
+        return bvh_node_is_leaf<T, Dim>(node); \
+    } \
+    BVH_EXPORT size_t bvh_node##suffix##_get_prim_count(const typename BvhTypes<T, Dim>::Node* node) { \
+        return bvh_node_get_prim_count<T, Dim>(node); \
+    } \
+    BVH_EXPORT void bvh_node##suffix##_set_prim_count(typename BvhTypes<T, Dim>::Node* node, size_t prim_count) { \
+        bvh_node_set_prim_count<T, Dim>(node, prim_count); \
+    } \
+    BVH_EXPORT size_t bvh_node##suffix##_get_first_id(const typename BvhTypes<T, Dim>::Node* node) { \
+        return bvh_node_get_first_id<T, Dim>(node); \
+    } \
+    BVH_EXPORT void bvh_node##suffix##_set_first_id(typename BvhTypes<T, Dim>::Node* node, size_t first_id) { \
+        bvh_node_set_first_id<T, Dim>(node, first_id); \
+    } \
+    BVH_EXPORT typename BvhTypes<T, Dim>::BBox bvh_node##suffix##_get_bbox(const typename BvhTypes<T, Dim>::Node* node) { \
+        return bvh_node_get_bbox<T, Dim>(node); \
+    } \
+    BVH_EXPORT void bvh_node##suffix##_set_bbox(typename BvhTypes<T, Dim>::Node* node, const typename BvhTypes<T, Dim>::BBox* bbox) { \
+        bvh_node_set_bbox<T, Dim>(node, bbox); \
+    } \
+    BVH_EXPORT void bvh##suffix##_append_node(typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        bvh_append_node<T, Dim>(bvh); \
+    } \
+    BVH_EXPORT void bvh##suffix##_remove_last_node(typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        bvh_remove_last_node<T, Dim>(bvh); \
+    } \
+    BVH_EXPORT void bvh##suffix##_refit(typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        bvh_refit<T, Dim>(bvh); \
+    } \
+    BVH_EXPORT void bvh##suffix##_optimize(bvh_thread_pool* thread_pool, typename BvhTypes<T, Dim>::Bvh* bvh) { \
+        bvh_optimize<T, Dim>(thread_pool, bvh); \
+    } \
+    BVH_EXPORT void bvh##suffix##_intersect_ray_any( \
+        const typename BvhTypes<T, Dim>::Bvh* bvh, \
+        const typename BvhTypes<T, Dim>::Ray* ray, \
+        const typename BvhCallback<T>::Type* callback) \
+    { \
+        bvh_intersect_ray<T, Dim, true, false>(bvh, ray, callback); \
+    } \
+    BVH_EXPORT void bvh##suffix##_intersect_ray_any_robust( \
+        const typename BvhTypes<T, Dim>::Bvh* bvh, \
+        const typename BvhTypes<T, Dim>::Ray* ray, \
+        const typename BvhCallback<T>::Type* callback) \
+    { \
+        bvh_intersect_ray<T, Dim, true, true>(bvh, ray, callback); \
+    } \
+    BVH_EXPORT void bvh##suffix##_intersect_ray( \
+        const typename BvhTypes<T, Dim>::Bvh* bvh, \
+        const typename BvhTypes<T, Dim>::Ray* ray, \
+        const typename BvhCallback<T>::Type* callback) \
+    { \
+        bvh_intersect_ray<T, Dim, false, false>(bvh, ray, callback); \
+    } \
+    BVH_EXPORT void bvh##suffix##_intersect_ray_robust( \
+        const typename BvhTypes<T, Dim>::Bvh* bvh, \
+        const typename BvhTypes<T, Dim>::Ray* ray, \
+        const typename BvhCallback<T>::Type* callback) \
+    { \
+        bvh_intersect_ray<T, Dim, false, true>(bvh, ray, callback); \
+    }
+
+} // namespace bvh::v2::c_api
+
+#endif
diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt
index bdf92254..61c58228 100644
--- a/test/CMakeLists.txt
+++ b/test/CMakeLists.txt
@@ -5,14 +5,35 @@ endif()
 add_executable(simple_example simple_example.cpp)
 target_link_libraries(simple_example PUBLIC bvh)
 set_target_properties(simple_example PROPERTIES CXX_STANDARD 20)
+add_test(NAME simple_example COMMAND simple_example)
 
 add_executable(serialize serialize.cpp)
 target_link_libraries(serialize PUBLIC bvh)
 set_target_properties(serialize PROPERTIES CXX_STANDARD 20)
+add_test(NAME serialize COMMAND serialize)
 
 add_executable(benchmark benchmark.cpp load_obj.cpp)
 target_link_libraries(benchmark PUBLIC bvh)
 set_target_properties(benchmark PROPERTIES CXX_STANDARD 20)
+add_test(
+    NAME benchmark
+    COMMAND
+        benchmark ${CMAKE_CURRENT_SOURCE_DIR}/scenes/cornell_box.obj
+        --eye 0 1 2
+        --dir 0 0 -1
+        --up 0 1 0)
 
-add_test(NAME simple_example COMMAND simple_example)
-add_test(NAME serialize COMMAND serialize)
+if (BVH_BUILD_C_API)
+    add_executable(c_api_example c_api_example.c load_obj.cpp)
+    target_link_libraries(c_api_example PUBLIC bvh_c)
+    set_target_properties(c_api_example PROPERTIES CXX_STANDARD 20 C_STANDARD 11)
+
+    add_test(
+        NAME c_api_example
+        COMMAND
+            c_api_example
+            ${CMAKE_CURRENT_SOURCE_DIR}/scenes/cornell_box.obj
+            --eye 0 1 2
+            --dir 0 0 -1
+            --up 0 1 0)
+endif()
diff --git a/test/benchmark.cpp b/test/benchmark.cpp
index 22b75109..0f46ead5 100644
--- a/test/benchmark.cpp
+++ b/test/benchmark.cpp
@@ -429,6 +429,6 @@ int main(int argc, char** argv) {
     }
 
     image.save(options.output_image);
-    std::cout << "Image saved as " << options.output_image << std::endl;
+    std::cout << "Image saved as '" << options.output_image << "'" << std::endl;
     return 0;
 }
diff --git a/test/c_api_example.c b/test/c_api_example.c
new file mode 100644
index 00000000..237bf1dd
--- /dev/null
+++ b/test/c_api_example.c
@@ -0,0 +1,366 @@
+#include <bvh/v2/c_api/bvh.h>
+
+#include "load_obj.h"
+
+#include <stddef.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <stdint.h>
+#include <string.h>
+#include <inttypes.h>
+#include <float.h>
+#include <time.h>
+#include <stdbool.h>
+#include <assert.h>
+#include <tgmath.h>
+
+#define INVALID_TRI_ID UINT32_MAX
+
+struct scene {
+    struct bvh3f* bvh;
+    struct tri* tris;
+    size_t tri_count;
+};
+
+struct options {
+    const char* input_scene;
+    const char* output_image;
+    uint32_t width, height;
+    struct bvh_vec3f eye, dir, up;
+};
+
+struct image {
+    uint8_t* pixels;
+    uint32_t width, height;
+};
+
+struct hit {
+    float t, u, v;
+    uint32_t tri_id;
+};
+
+struct intersect_user_data {
+    const struct bvh_ray3f* ray;
+    const struct scene* scene;
+    struct hit* hit;
+};
+
+static inline float dot(struct bvh_vec3f v, struct bvh_vec3f w) {
+    return v.x * w.x + v.y * w.y + v.z * w.z;
+}
+
+static inline float length(struct bvh_vec3f v) {
+    return dot(v, v);
+}
+
+static inline struct bvh_vec3f sub(struct bvh_vec3f v, struct bvh_vec3f w) {
+    return (struct bvh_vec3f) {
+        .x = v.x - w.x,
+        .y = v.y - w.y,
+        .z = v.z - w.z,
+    };
+}
+
+static inline struct bvh_vec3f cross(struct bvh_vec3f v, struct bvh_vec3f w) {
+    return (struct bvh_vec3f) {
+        .x = v.y * w.z - v.z * w.y,
+        .y = v.z * w.x - v.x * w.z,
+        .z = v.x * w.y - v.y * w.x,
+    };
+}
+
+static inline struct bvh_vec3f normalize(struct bvh_vec3f v) {
+    const float inv_len = 1.f / length(v);
+    return (struct bvh_vec3f) {
+        .x = v.x * inv_len,
+        .y = v.y * inv_len,
+        .z = v.z * inv_len,
+    };
+}
+
+static inline struct bvh_vec3f tri_center(const struct tri* tri) {
+    return (struct bvh_vec3f) {
+        .x = (tri->v[0].x + tri->v[1].x + tri->v[2].x) * (1.f / 3.f),
+        .y = (tri->v[0].y + tri->v[1].y + tri->v[2].y) * (1.f / 3.f),
+        .z = (tri->v[0].z + tri->v[1].z + tri->v[2].z) * (1.f / 3.f)
+    };
+}
+
+static inline struct bvh_bbox3f tri_bbox(const struct tri* tri) {
+    struct bvh_vec3f min = tri->v[0], max = min;
+    for (int i = 1; i < 3; ++i) {
+        min.x = min.x < tri->v[i].x ? min.x : tri->v[i].x;
+        min.y = min.y < tri->v[i].y ? min.y : tri->v[i].y;
+        min.z = min.z < tri->v[i].z ? min.z : tri->v[i].z;
+        max.x = max.x > tri->v[i].x ? max.x : tri->v[i].x;
+        max.y = max.y > tri->v[i].y ? max.y : tri->v[i].y;
+        max.z = max.z > tri->v[i].z ? max.z : tri->v[i].z;
+    }
+    return (struct bvh_bbox3f) { .min = min, .max = max };
+}
+
+static inline void build_bvh(struct scene* scene) {
+    assert(!scene->bvh);
+    struct bvh_bbox3f* bboxes = malloc(sizeof(struct bvh_bbox3f) * scene->tri_count);
+    struct bvh_vec3f* centers = malloc(sizeof(struct bvh_vec3f) * scene->tri_count);
+    for (size_t i = 0; i < scene->tri_count; ++i) {
+        bboxes[i]  = tri_bbox(&scene->tris[i]);
+        centers[i] = tri_center(&scene->tris[i]);
+    }
+
+    struct bvh_thread_pool* thread_pool = bvh_thread_pool_create(0);
+    struct timespec ts0, ts1;
+    timespec_get(&ts0, TIME_UTC);
+    scene->bvh = bvh3f_build(thread_pool, bboxes, centers, scene->tri_count, NULL);
+    timespec_get(&ts1, TIME_UTC);
+    bvh_thread_pool_destroy(thread_pool);
+    free(bboxes);
+    free(centers);
+    const uint64_t build_time = (ts1.tv_sec - ts0.tv_sec) * 1000 + (ts1.tv_nsec - ts0.tv_nsec) / 1000000;
+    printf("Built BVH with %zu node(s) in %"PRIu64"\n", bvh3f_get_node_count(scene->bvh), build_time);
+}
+
+static inline void destroy_scene(struct scene* scene) {
+    bvh3f_destroy(scene->bvh);
+    free(scene->tris);
+    memset(scene, 0, sizeof(struct scene));
+}
+
+struct image alloc_image(size_t width, size_t height) {
+    return (struct image) {
+        .pixels = malloc(sizeof(uint8_t) * width * height * 3),
+        .width = width,
+        .height = height
+    };
+}
+
+static void free_image(struct image* image) {
+    free(image->pixels);
+    memset(image, 0, sizeof(struct image));
+}
+
+static bool save_image(const struct image* image, const char* file_name) {
+    FILE* file = fopen(file_name, "wb");
+    if (!file)
+        return false;
+    fprintf(file, "P6 %"PRIu32" %"PRIu32" 255\n", image->width, image->height);
+    for (uint32_t y = image->height; y-- > 0;)
+        fwrite(image->pixels + y * 3 * image->width, sizeof(uint8_t), 3 * image->width, file);
+    fclose(file);
+    return true;
+}
+
+static void set_pixel(struct image* image, uint32_t x, uint32_t y, uint8_t r, uint8_t g, uint8_t b) {
+    uint8_t* p = &image->pixels[(y * image->width + x) * 3];
+    p[0] = r;
+    p[1] = g;
+    p[2] = b;
+}
+
+static void usage() {
+    printf(
+        "Usage: c_api_example [options] file.obj\n"
+        "\nOptions:\n"
+        "  -h  --help                  Shows this message.\n"
+        "  -e  --eye <x> <y> <z>       Sets the position of the camera.\n"
+        "  -d  --dir <x> <y> <z>       Sets the direction of the camera.\n"
+        "  -u  --up  <x> <y> <z>       Sets the up vector of the camera.\n"
+        "  -w  --width <pixels>        Sets the image width.\n"
+        "  -h  --height <pixels>       Sets the image height.\n"
+        "  -o <file.ppm>               Sets the output file name (defaults to 'render.ppm').\n");
+}
+
+static inline bool check_arg(int i, int n, int argc, char** argv) {
+    if (i + n > argc) {
+        fprintf(stderr, "Missing argument(s) for '%s'\n", argv[i]);
+        return false;
+    }
+    return true;
+}
+
+static inline bool parse_options(int argc, char** argv, struct options* options) {
+    for (int i = 1; i < argc; ++i) {
+        if (argv[i][0] == '-') {
+            if (!strcmp(argv[i], "-e") || !strcmp(argv[i], "--eye")) {
+                if (!check_arg(i, 3, argc, argv))
+                    return false;
+                options->eye.x = strtof(argv[++i], NULL);
+                options->eye.y = strtof(argv[++i], NULL);
+                options->eye.z = strtof(argv[++i], NULL);
+            } else if (!strcmp(argv[i], "-d") || !strcmp(argv[i], "--dir")) {
+                if (!check_arg(i, 3, argc, argv))
+                    return false;
+                options->dir.x = strtof(argv[++i], NULL);
+                options->dir.y = strtof(argv[++i], NULL);
+                options->dir.z = strtof(argv[++i], NULL);
+            } else if (!strcmp(argv[i], "-u") || !strcmp(argv[i], "--up")) {
+                if (!check_arg(i, 3, argc, argv))
+                    return false;
+                options->up.x = strtof(argv[++i], NULL);
+                options->up.y = strtof(argv[++i], NULL);
+                options->up.z = strtof(argv[++i], NULL);
+            } else if (!strcmp(argv[i], "-w") || !strcmp(argv[i], "--width")) {
+                if (!check_arg(i, 1, argc, argv))
+                    return false;
+                options->width = strtoul(argv[++i], NULL, 10);
+            } else if (!strcmp(argv[i], "-h") || !strcmp(argv[i], "--height")) {
+                if (!check_arg(i, 1, argc, argv))
+                    return false;
+                options->height = strtoul(argv[++i], NULL, 10);
+            } else if (!strcmp(argv[i], "-o")) {
+                if (!check_arg(i, 1, argc, argv))
+                    return false;
+                options->output_image = argv[++i];
+            } else if (!strcmp(argv[i], "-h") || !strcmp(argv[i], "--help")) {
+                usage();
+                return false;
+            }
+        } else {
+            if (options->input_scene) {
+                fprintf(stderr, "Input scene file specified twice\n");
+                return false;
+            }
+            options->input_scene = argv[i];
+        }
+    }
+    if (!options->input_scene) {
+        fprintf(stderr, "Missing input scene file\n");
+        return false;
+    }
+    return true;
+}
+
+static inline bool intersect_ray_tri(const struct bvh_ray3f* ray, const struct tri* tri, struct hit* hit) {
+    const struct bvh_vec3f e1 = sub(tri->v[0], tri->v[1]);
+    const struct bvh_vec3f e2 = sub(tri->v[2], tri->v[0]);
+    const struct bvh_vec3f n = cross(e1, e2);
+    const struct bvh_vec3f c = sub(tri->v[0], ray->org);
+    const struct bvh_vec3f r = cross(ray->dir, c);
+
+    const float inv_det = 1.f / dot(n, ray->dir);
+
+    const float u = dot(r, e2) * inv_det;
+    const float v = dot(r, e1) * inv_det;
+    const float w = 1.f - u - v;
+
+    static const float tolerance = -FLT_EPSILON;
+    if (u >= tolerance && v >= tolerance && w >= tolerance) {
+        const float t = dot(n, c) * inv_det;
+        if (t >= ray->tmin && t <= hit->t) {
+            hit->t = t;
+            hit->u = u;
+            hit->v = v;
+            return true;
+        }
+    }
+    return false;
+}
+
+static bool intersect_bvh_leaf(void* user_data, float* t, size_t begin, size_t end) {
+    const struct bvh_ray3f* ray = ((struct intersect_user_data*)user_data)->ray;
+    const struct scene* scene = ((struct intersect_user_data*)user_data)->scene;
+    struct hit* hit = ((struct intersect_user_data*)user_data)->hit;
+    bool was_hit = false;
+    for (size_t i = begin; i < end; ++i) {
+        const uint32_t tri_id = bvh3f_get_prim_id(scene->bvh, i);
+        if (intersect_ray_tri(ray, &scene->tris[tri_id], hit)) {
+            *t = hit->t;
+            hit->tri_id = tri_id;
+            was_hit = true;
+        }
+    }
+    return was_hit;
+}
+
+static inline void render_image(
+    const struct scene* scene,
+    const struct options* options,
+    struct image* image)
+{
+    const struct bvh_vec3f dir   = normalize(options->dir);
+    const struct bvh_vec3f right = normalize(cross(dir, options->up));
+    const struct bvh_vec3f up    = cross(right, dir);
+
+    struct intersect_user_data user_data = { .scene = scene };
+    const struct bvh_intersect_callbackf callback = {
+        .user_data = &user_data,
+        .user_fn = intersect_bvh_leaf 
+    };
+
+    size_t intr_count = 0;
+    struct timespec ts0, ts1;
+    timespec_get(&ts0, TIME_UTC);
+    for (uint32_t y = 0; y < options->height; ++y) {
+        for (uint32_t x = 0; x < options->width; ++x) {
+            const float u = 2.f * ((float)x) / ((float)options->width)  - 1.f;
+            const float v = 2.f * ((float)y) / ((float)options->height) - 1.f;
+
+            const struct bvh_ray3f ray = {
+                .org = options->eye,
+                .dir = {
+                    .x = options->dir.x + u * right.x + v * up.x,
+                    .y = options->dir.y + u * right.y + v * up.y,
+                    .z = options->dir.z + u * right.z + v * up.z
+                },
+                .tmin = 0.f,
+                .tmax = FLT_MAX
+            };
+
+            struct hit hit = { .tri_id = INVALID_TRI_ID, .t = FLT_MAX };
+
+            user_data.ray = &ray;
+            user_data.hit = &hit;
+            bvh3f_intersect_ray(scene->bvh, &ray, &callback);
+
+            uint8_t r = 0;
+            uint8_t g = 0;
+            uint8_t b = 0;
+            if (hit.tri_id != INVALID_TRI_ID) {
+                r = (hit.tri_id * 79) % 255 + 1;
+                g = (hit.tri_id * 43) % 255 + 1;
+                b = (hit.tri_id * 57) % 255 + 1;
+                intr_count++;
+            }
+            set_pixel(image, x, y, r, g, b);
+        }
+    }
+    timespec_get(&ts1, TIME_UTC);
+    const uint64_t render_time = (ts1.tv_sec - ts0.tv_sec) * 1000 + (ts1.tv_nsec - ts0.tv_nsec) / 1000000;
+    printf("%zu intersection(s) found in %"PRIu64"\n", intr_count, render_time);
+}
+
+int main(int argc, char** argv) {
+    struct options options = {
+        .eye = (struct bvh_vec3f) { 0, 0, 0 },
+        .dir = (struct bvh_vec3f) { 0, 0, 1 },
+        .up  = (struct bvh_vec3f) { 0, 1, 0 },
+        .output_image = "render.ppm",
+        .width  = 1024,
+        .height = 1024
+    };
+    if (!parse_options(argc, argv, &options))
+        return 1;
+
+    struct scene scene = { 0 };
+    scene.tris = load_obj(argv[1], &scene.tri_count);
+    if (!scene.tris) {
+        fprintf(stderr, "Invalid or empty OBJ file\n");
+        return 1;
+    }
+    printf("Loaded file with %zu triangle(s)\n", scene.tri_count);
+
+    build_bvh(&scene);
+
+    struct image image = alloc_image(options.width, options.height);
+    render_image(&scene, &options, &image);
+    if (!save_image(&image, options.output_image)) {
+        fprintf(stderr, "Could not save rendered image to '%s'\n", options.output_image);
+        return 1;
+    } else {
+        printf("Image saved as '%s'\n", options.output_image);
+    }
+
+    free_image(&image);
+    destroy_scene(&scene);
+    return 0;
+}
diff --git a/test/load_obj.cpp b/test/load_obj.cpp
index 13018b1d..5f77083f 100644
--- a/test/load_obj.cpp
+++ b/test/load_obj.cpp
@@ -112,3 +112,27 @@ std::vector<Tri<T, 3>> load_obj(const std::string& file) {
 
 template std::vector<Tri<float, 3>> load_obj(const std::string&);
 template std::vector<Tri<double, 3>> load_obj(const std::string&);
+
+extern "C" {
+
+tri* load_obj(const char* file_name, size_t* tri_count) {
+    auto tri_vec = load_obj<float>(file_name);
+    *tri_count = tri_vec.size();
+    if (tri_vec.size() == 0)
+        return NULL;
+    tri* tris = reinterpret_cast<tri*>(malloc(sizeof(tri) * tri_vec.size()));
+    for (size_t i = 0; i < tri_vec.size(); ++i) {
+        tris[i].v[0].x = tri_vec[i].p0[0];
+        tris[i].v[0].y = tri_vec[i].p0[1];
+        tris[i].v[0].z = tri_vec[i].p0[2];
+        tris[i].v[1].x = tri_vec[i].p1[0];
+        tris[i].v[1].y = tri_vec[i].p1[1];
+        tris[i].v[1].z = tri_vec[i].p1[2];
+        tris[i].v[2].x = tri_vec[i].p2[0];
+        tris[i].v[2].y = tri_vec[i].p2[1];
+        tris[i].v[2].z = tri_vec[i].p2[2];
+    }
+    return tris;
+}
+
+} // extern "C"
diff --git a/test/load_obj.h b/test/load_obj.h
index 830f78dc..97a24446 100644
--- a/test/load_obj.h
+++ b/test/load_obj.h
@@ -1,6 +1,22 @@
 #ifndef LOAD_OBJ_H
 #define LOAD_OBJ_H
 
+#include <bvh/v2/c_api/bvh.h>
+
+#include <stddef.h>
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+struct tri { struct bvh_vec3f v[3]; };
+struct tri* load_obj(const char*, size_t*);
+
+#ifdef __cplusplus
+}
+#endif
+
+#ifdef __cplusplus
 #include <bvh/v2/tri.h>
 
 #include <vector>
@@ -8,5 +24,6 @@
 
 template <typename T>
 std::vector<bvh::v2::Tri<T, 3>> load_obj(const std::string& file);
+#endif
 
 #endif
diff --git a/test/scenes/cornell_box.obj b/test/scenes/cornell_box.obj
new file mode 100644
index 00000000..6309495d
--- /dev/null
+++ b/test/scenes/cornell_box.obj
@@ -0,0 +1,168 @@
+# The original Cornell Box in OBJ format.
+# Note that the real box is not a perfect cube, so
+# the faces are imperfect in this data set.
+#
+# Created by Guedis Cardenas and Morgan McGuire at Williams College, 2011
+# Released into the Public Domain.
+#
+# http://graphics.cs.williams.edu/data
+# http://www.graphics.cornell.edu/online/box/data.html
+#
+
+mtllib cornell_box.mtl
+
+## Object floor
+v  -1.01  0.00   0.99
+v   1.00  0.00   0.99
+v   1.00  0.00  -1.04
+v  -0.99  0.00  -1.04
+
+g floor
+usemtl floor
+f -4 -3 -2 -1
+
+## Object ceiling 
+v  -1.02  1.99   0.99  
+v  -1.02  1.99  -1.04
+v   1.00  1.99  -1.04
+v   1.00  1.99   0.99
+
+g ceiling
+usemtl ceiling
+f -4 -3 -2 -1
+
+## Object backwall  
+v  -0.99  0.00  -1.04 
+v   1.00  0.00  -1.04
+v   1.00  1.99  -1.04
+v  -1.02  1.99  -1.04
+
+g backWall
+usemtl backWall
+f -4 -3 -2 -1
+
+## Object rightwall 
+v	1.00  0.00  -1.04   
+v	1.00  0.00   0.99
+v	1.00  1.99   0.99
+v	1.00  1.99  -1.04
+
+g rightWall
+usemtl rightWall
+f -4 -3 -2 -1
+
+## Object leftWall   
+v  -1.01  0.00   0.99
+v  -0.99  0.00  -1.04
+v  -1.02  1.99  -1.04
+v  -1.02  1.99   0.99
+
+g leftWall
+usemtl leftWall
+f -4 -3 -2 -1
+
+## Object shortBox 
+usemtl shortBox
+
+# Top Face 
+v	0.53  0.60   0.75 
+v	0.70  0.60   0.17  
+v	0.13  0.60   0.00
+v  -0.05  0.60   0.57
+f -4 -3 -2 -1
+
+# Left Face 
+v  -0.05  0.00   0.57
+v  -0.05  0.60   0.57
+v   0.13  0.60   0.00 
+v   0.13  0.00   0.00
+f -4 -3 -2 -1
+
+# Front Face
+v	0.53  0.00   0.75
+v	0.53  0.60   0.75
+v  -0.05  0.60   0.57
+v  -0.05  0.00   0.57
+f -4 -3 -2 -1
+
+# Right Face
+v	0.70  0.00   0.17
+v	0.70  0.60   0.17
+v	0.53  0.60   0.75
+v	0.53  0.00   0.75
+f -4 -3 -2 -1
+
+# Back Face
+v	0.13  0.00   0.00
+v	0.13  0.60   0.00
+v	0.70  0.60   0.17
+v	0.70  0.00   0.17
+f -4 -3 -2 -1 
+
+# Bottom Face
+v	0.53  0.00   0.75 
+v	0.70  0.00   0.17  
+v	0.13  0.00   0.00
+v  -0.05  0.00   0.57
+f -12 -11 -10 -9
+
+g shortBox
+usemtl shortBox
+
+## Object tallBox 
+usemtl tallBox
+
+# Top Face
+v	-0.53  1.20   0.09
+v	 0.04  1.20  -0.09
+v	-0.14  1.20  -0.67
+v	-0.71  1.20  -0.49
+f -4 -3 -2 -1
+
+# Left Face
+v	-0.53  0.00   0.09
+v	-0.53  1.20   0.09
+v	-0.71  1.20  -0.49
+v	-0.71  0.00  -0.49
+f -4 -3 -2 -1
+
+# Back Face
+v	-0.71  0.00  -0.49
+v	-0.71  1.20  -0.49
+v	-0.14  1.20  -0.67
+v	-0.14  0.00  -0.67
+f -4 -3 -2 -1
+
+# Right Face
+v	-0.14  0.00  -0.67
+v	-0.14  1.20  -0.67
+v	 0.04  1.20  -0.09
+v	 0.04  0.00  -0.09
+f -4 -3 -2 -1
+
+# Front Face
+v	 0.04  0.00  -0.09
+v	 0.04  1.20  -0.09
+v	-0.53  1.20   0.09
+v	-0.53  0.00   0.09
+f -4 -3 -2 -1
+
+# Bottom Face
+v	-0.53  0.00   0.09
+v	 0.04  0.00  -0.09
+v	-0.14  0.00  -0.67
+v	-0.71  0.00  -0.49
+f -8 -7 -6 -5
+
+g tallBox
+usemtl tallBox
+
+## Object light 
+v	-0.24  1.98   0.16
+v	-0.24  1.98  -0.22
+v	 0.23  1.98  -0.22
+v	 0.23  1.98   0.16
+
+g light
+usemtl light
+f -4 -3 -2 -1