Location>code7788 >text

games101 Assignment 6 Detailed SAH

Popularity:346 ℃/2024-08-15 22:15:46

games101 Assignment 6 Detailed SAH

References available:/3ed-2018/Primitives_and_Intersection_Acceleration/Bounding_Volume_Hierarchies

code analysis

The code framework for Assignment 6 has been further enriched and improved from that of Assignment 5
The first step was to create the scene, add the objects and the lights, this time the only object was a 3d model of a rabbit that was read with the objloader, but when converting the model to a meshtriangle data structure, a BVH (Bounding volume hierarchies) acceleration structure was also created for the triangles in the model. here. The BVHs are also created for the objects in the scene, but there is only one object in the scene at the moment.

MeshTriangle bunny(Utils::PathFromAsset("model/bunnyAssignment6/"));

(&bunny);
(std::make_unique<Light>(Vector3f(-20, 70, 20), 1));
(std::make_unique<Light>(Vector3f(20, 70, 20), 1));
();

MeshTriangle(const std::string& filename){
.......
std::vector<Object*> ptrs;
for (auto& tri : triangles)
    ptrs.push_back(&tri);

bvh = new BVHAccel(ptrs);
}

BVHAccel::BVHAccel(std::vector<Object*> p, int maxPrimsInNode,
                   SplitMethod splitMethod)
    : maxPrimsInNode(std::min(255, maxPrimsInNode)), splitMethod(splitMethod),
      primitives(std::move(p))
{
    time_t start, stop;
    time(&start);
    if (())
        return;

    root = recursiveBuild(primitives);

    time(&stop);
    double diff = difftime(stop, start);
    int hrs = (int)diff / 3600;
    int mins = ((int)diff / 60) - (hrs * 60);
    int secs = (int)diff - (hrs * 3600) - (mins * 60);

    printf(
        "\rBVH Generation complete: \nTime Taken: %i hrs, %i mins, %i secs\n\n",
        hrs, mins, secs);
}

The process of building a BVH is the recursive construction of a tree The method of dividing the BVH nodes is to sort all the triangles in the longest dimension to find out the position of the triangle in the middle The left node is the left node The right node is the right node, and the termination condition for the end of the division is that there is only one object.
The PowerPoint that corresponds to the lesson:
img
pbrt illustration:
img

BVHBuildNode* BVHAccel::recursiveBuild(std::vector<Object*> objects)
{
    BVHBuildNode* node = new BVHBuildNode();

    // Compute bounds of all primitives in BVH node
    Bounds3 bounds;
    for (int i = 0; i < (); ++i)
        bounds = Union(bounds, objects[i]->getBounds());
    if (() == 1) {
        // Create leaf _BVHBuildNode_
        node->bounds = objects[0]->getBounds();
        node->object = objects[0];
        node->left = nullptr;
        node->right = nullptr;
        return node;
    }
    else if (() == 2) {
        node->left = recursiveBuild(std::vector{objects[0]});
        node->right = recursiveBuild(std::vector{objects[1]});

        node->bounds = Union(node->left->bounds, node->right->bounds);
        return node;
    }
    else {
        //Use each of the triangles corresponding to theboxmidpoint of a triangle put sth. togetherbox Then find out what thisboxLongest dimension Segmentation on this dimension
        Bounds3 centroidBounds;
        for (int i = 0; i < (); ++i)
            centroidBounds =
                Union(centroidBounds, objects[i]->getBounds().Centroid());
        int dim = ();
        switch (dim) {
        case 0:
            std::sort((), (), [](auto f1, auto f2) {
                return f1->getBounds().Centroid().x <
                       f2->getBounds().Centroid().x;
            });
            break;
        case 1:
            std::sort((), (), [](auto f1, auto f2) {
                return f1->getBounds().Centroid().y <
                       f2->getBounds().Centroid().y;
            });
            break;
        case 2:
            std::sort((), (), [](auto f1, auto f2) {
                return f1->getBounds().Centroid().z <
                       f2->getBounds().Centroid().z;
            });
            break;
        }

        auto beginning = ();
        auto middling = () + (() / 2);
        auto ending = ();

        auto leftshapes = std::vector<Object*>(beginning, middling);
        auto rightshapes = std::vector<Object*>(middling, ending);

        assert(() == (() + ()));

        node->left = recursiveBuild(leftshapes);
        node->right = recursiveBuild(rightshapes);

        node->bounds = Union(node->left->bounds, node->right->bounds);
    }

    return node;
}

After that, it's ray tracing, castRay, intersection, and then the coloring calculation, similar to Job 5.
But this intersection is going to be accelerated using the BVH acceleration structure that we've built, and that's part of what we're going to accomplish.

Analysis and resolution

Intersect with BVH

The main purpose here is to traverse the constructed BVH accelerated structure and intersect the boundingbox of a node. If there is no intersection, return. If there is, check whether the node is a leaf node or not. If it is not, then continue to traverse the node's children. If it is a leaf node, then test the intersection of the objects in it to find the nearest intersection point and return:


Intersection BVHAccel::getIntersection(BVHBuildNode* node, const Ray& ray) const
{
    // TODO Traverse the BVH to find intersection
    Intersection hit1;
    Intersection hit2;
    std::array<int, 3> dirIsNeg{ 0,0,0 };
    for (int i = 0; i < 3; i++)
    {
        if ([i] < 0) {
            dirIsNeg[i] = 1;
        }
    }
    if (!node->(ray, ray.direction_inv, std::move(dirIsNeg)))
    {
        return Intersection{};
    }
    if (node->left == nullptr && node->right == nullptr) {
        return node->object->getIntersection(ray);
    }
    hit1 = getIntersection(node->left, ray);
    hit2 = getIntersection(node->right, ray);
    
    return  >  ? hit2 : hit1;
}

boundingbox intersection

The boundingbox is a rectangle and has faces parallel to the xyz-axis, so all we need to do is find the tenter and texit of each set of planes, and find the maximum value of enter and the minimum value of exit, which prevents us from incorrectly locating the intersection point on the extended face of the rectangle, and also takes into account that the dir line is negative, so that it will intersect with planes of smaller coordinates first, and then we'll find the intersection point at the end of the line. This will prevent you from incorrectly locating the intersection point on the extension plane of the rectangle:

inline bool Bounds3::IntersectP(const Ray& ray, const Vector3f& invDir,
                                const std::array<int, 3>& dirIsNeg) const
{
    // invDir: ray direction(x,y,z), invDir=(1.0/x,1.0/y,1.0/z), use this because Multiply is faster that Division
    // dirIsNeg: ray direction(x,y,z), dirIsNeg=[int(x>0),int(y>0),int(z>0)], use this to simplify your logic
    // TODO test if ray bound intersects
    float txmin = (pMin - ).x * ;
    float txmax = (pMax - ).x * ;
    float tymin = (pMin - ).y * ;
    float tymax = (pMax - ).y * ;
    float tzmin = (pMin - ).z * ;
    float tzmax = (pMax - ).z * ;
    //If the direction is negative in one dimension grouping together withpMaxFaces in corresponding dimensions intersect first
    if (dirIsNeg[0]) {
        std::swap(txmin, txmax);
    }
    if (dirIsNeg[1]) {
        std::swap(tymin, tymax);
    }
    if (dirIsNeg[2]) {
        std::swap(tzmin, tzmax);
    }
    float tenter = std::max(std::max(txmin, tymin), tzmin);
    float texit = std::min(std::min(txmax, tymax), tzmax);
    return tenter < texit && texit >= 0;
}

Intersection of triangles

Similar to the implementation in Assignment 5
Here you need to return the properties of the intersection, such as material, normal, and so on:

inline Intersection Triangle::getIntersection(Ray ray)
{
    Intersection inter;
    //If the normal is on the opposite side of the ray This means that the light is coming in from the back side
    if (dotProduct(, normal) > 0)
        return inter;
    double u, v, t_tmp = 0;
    // D E2 = S1
    Vector3f pvec = crossProduct(, e2);
    // S1 . E1
    double det = dotProduct(e1, pvec);
    if (fabs(det) < EPSILON)
        return inter;

    double det_inv = 1. / det;
    //S
    Vector3f tvec = - v0;
    u = dotProduct(tvec, pvec) * det_inv;
    if (u < 0 || u > 1)
        return inter;
    //S2
    Vector3f qvec = crossProduct(tvec, e1);
    v = dotProduct(, qvec) * det_inv;
    if (v < 0 || u + v > 1)
        return inter;
    t_tmp = dotProduct(e2, qvec) * det_inv;

    // TODO find ray triangle intersection
    if (t_tmp < 0) {
        return inter;
    }
     = true;
     = ray(t_tmp);
     = normal;
     = t_tmp;
     = this;
     = this->m;


    return inter;
}

SAH

SAH, or Surface Area Heuristic, is a method of partitioning during the construction of a BVH. We have previously mentioned a simple partitioning method. The previous method was used in the case shown in the following figure.
img
If you choose a division like the blue line in the diagram, you'll end up with a partial overlap of the bounding box, which reduces the efficiency of light intersection.
SAH proposes a function to quantitatively estimate the cost of dividing a region Take the example of dividing a boundingbox into two ABs:
img
where t-trav is the cost of traversing the tree nodes, here fixed to a constant 0.125 pbrt is explained by the fact that we set the cost of light intersection to 1, and the traversal efficiency is roughly 1/8th that of intersection, it is the relative relationship between them that is important, not the numerical value.
The t-isect represents the cost of intersection, which is set to the number of tuples in each region.
pa pb represents the probability that a ray of light will pass through a child node using the ratio of the boundingbox area of A to the total boundingbox area.

The specific process is summarized below:
1. First we find the longest dimension as in the previous division algorithm.
2. Here we define a certain number of buckets to be divided on the overall boundingbox and calculate the index of the bucket where each object is located.
3. Each bucket represents a partitioning method. we iterate over the cost of each method and find the bucket index corresponding to the smallest cost:
img
4. Split the objects into two piles with the help of minimum indexing

Code: reference /mmp/pbrt-v3/blob/master/src/accelerators/

BVHBuildNode* BVHAccel::SAHrecursiveBuild(std::vector<Object*> objects)
{
    BVHBuildNode* node = new BVHBuildNode();

    // Compute bounds of all primitives in BVH node
    Bounds3 bounds;
    for (int i = 0; i < (); ++i)
        bounds = Union(bounds, objects[i]->getBounds());
   //If the number of objects is too small Just use the original division directly unnecessary
    if (() <= 2) {
        node = recursiveBuild(objects);
        return node;
    }
    else {
        Bounds3 centroidBounds;
        for (int i = 0; i < (); ++i)
            centroidBounds =
            Union(centroidBounds, objects[i]->getBounds().Centroid());
        int dim = ();
        //Calculation Index
        constexpr int nBuckets = 8;
        BucketInfo buckets[nBuckets];
        for (int i = 0; i < (); ++i) {
            int b = nBuckets *
                (
                    objects[i]->getBounds().Centroid())[dim];
            if (b == nBuckets) b = nBuckets - 1;
            buckets[b].count++;
            buckets[b].bounds =
                Union(buckets[b].bounds, objects[i]->getBounds());
        }
        //Calculate the cost of each division
        float cost[nBuckets - 1];
        for (int i = 0; i < nBuckets - 1; ++i) {
            Bounds3 b0, b1;
            int count0 = 0, count1 = 0;
            for (int j = 0; j <= i; ++j) {
                b0 = Union(b0, buckets[j].bounds);
                count0 += buckets[j].count;
            }
            for (int j = i + 1; j < nBuckets; ++j) {
                b1 = Union(b1, buckets[j].bounds);
                count1 += buckets[j].count;
            }
            cost[i] = 0.125f +
                (count0 * () +
                    count1 * ()) /
                ();
        }
        //Finding the Minimum Cost
        float minCost = cost[0];
        int minCostSplitBucket = 0;
        for (int i = 1; i < nBuckets - 1; ++i) {
            if (cost[i] < minCost) {
                minCost = cost[i];
                minCostSplitBucket = i;
            }
        }
        float leafCost = ();
        int mid = 0;
        //If too many objects are included Or the cost of narrowing the division Then divide
        if (() > maxPrimsInNode || minCost < leafCost) {
            for (int i = 0; i < (); ++i) {
                int b = nBuckets *
                    (
                        objects[i]->getBounds().Centroid())[dim];
                if (b == nBuckets) b = nBuckets - 1;

                // The elements that satisfy the condition are placed in the first half of the array.
                if (b <= minCostSplitBucket) {
                    std::swap(objects[i], objects[mid]);
                    mid++;
                }
            }
        }
       
        auto beginning = ();
        auto middling = () + mid;
        auto ending = ();

        auto leftshapes = std::vector<Object*>(beginning, middling);
        auto rightshapes = std::vector<Object*>(middling, ending);

        assert(() == (() + ()));

        node->left = SAHrecursiveBuild(leftshapes);
        node->right = SAHrecursiveBuild(rightshapes);

        node->bounds = Union(node->left->bounds, node->right->bounds);
    }

    return node;
}

Results Showcase

BVH:
img
BVH-SAH:
img
img