C++ BBox::SurfaceArea方法代码示例

TEST_F(BBoxTest, SurfaceAreaDegenerateBBoxWorks) {
	BBox b (Point(0, 0, 0), Point (0, 0, 0));

	float area = b.SurfaceArea();

	EXPECT_FLOAT_EQ (0.f, area);
}

BVHNode* BVH::RecursiveBuild(uint32_t start, uint32_t end, uint32_t depth)
{
    maxDepth = fmaxf(depth, maxDepth);
    totalNodes++;
    BVHNode* node = new BVHNode;

    //compute bounds of all primitives in BVH node
    BBox bbox;
    for(auto i = start; i < end; ++i)
        bbox = Union(bbox, workList[i].bounds);

    uint32_t nPrims = end - start;
    //if number of primitives are less than threshold, create leaf node
    if(nPrims <= MAX_LEAF_PRIM_NUM)
    {
        uint32_t firstPrimOffset = orderedPrims.size();
        for(auto i = start; i < end; ++i)
        {
            auto pIdx = workList[i].pIdx;
            orderedPrims.push_back(mesh.faces[pIdx]);
        }
        node->InitLeaf(firstPrimOffset, nPrims, bbox);
    }
    else
    {
        //compute bound of primitive centroids, choose split dimension
        BBox centroidBounds;
        for(auto i = start; i < end; ++i)
            centroidBounds = Union(centroidBounds, workList[i].bounds.bcenter);

        // split along max span axis
        int dim = centroidBounds.MaxExtent();

        //partition primitives into two sets and build children
        uint32_t mid = (end + start) / 2;
        // if max span axis is too small, create a leaf node
        if((centroidBounds.bmax[dim] - centroidBounds.bmin[dim]) < 1e-4)
        {
            uint32_t firstPrimOffset = orderedPrims.size();
            for(auto i = start; i < end; ++i)
            {
                auto pIdx = workList[i].pIdx;
                orderedPrims.push_back(mesh.faces[pIdx]);
            }
            node->InitLeaf(firstPrimOffset, nPrims, bbox);

            return node;
        }

        //partition primitives based on SAH
        std::vector<BucketInfo> buckets(nBuckets);
        float extent = centroidBounds.bmax[dim] - centroidBounds.bmin[dim];
        for(auto i = start; i < end; ++i)
        {
            uint32_t b = nBuckets * ((workList[i].bounds.bcenter[dim] - centroidBounds.bmin[dim]) / extent);
            if(b == nBuckets) b -= 1;
            buckets[b].count++;
            buckets[b].bounds = Union(buckets[b].bounds, workList[i].bounds);
        }

        //compute costs for splitting after each bucket
        float cost[nBuckets - 1];
        for(auto i = 0; i < nBuckets - 1; ++i)
        {
            BBox b0, b1;
            int count0 = 0, count1 = 0;

            for(auto j = 0; j <= i; ++j)
            {
                b0 = Union(b0, buckets[j].bounds);
                count0 += buckets[j].count;
            }
            for(auto j = i + 1; j < nBuckets; ++j)
            {
                b1 = Union(b1, buckets[j].bounds);
                count1 += buckets[j].count;
            }

            cost[i] = (count0 * b0.SurfaceArea() + count1 * b1.SurfaceArea()) / bbox.SurfaceArea();
        }

        //find best split
        float minCost = cost[0];
        uint32_t bestSplit = 0;
        for(auto i = 1; i < nBuckets - 1; ++i)
        {
            if(cost[i] < minCost)
            {
                minCost = cost[i];
                bestSplit = i;
            }
        }

        //either create leaf or split at selected SAH bucket
        if(nPrims > MAX_LEAF_PRIM_NUM || minCost < nPrims)
        {
            auto compare = [&](BVHPrimitiveInfo& p) {
                auto b = nBuckets * ((p.bounds.bcenter[dim] - centroidBounds.bmin[dim]) / extent);
                b = (b == nBuckets) ? (b - 1) : b;
                return b <= bestSplit;
//.........这里部分代码省略.........

//.........这里部分代码省略.........
            if (mid != start && mid != end)
                // for lots of prims with large overlapping bounding boxes, this
                // may fail to partition; in that case don't break and fall through
                // to SPLIT_EQUAL_COUNTS
                break;
        }
        case SPLIT_EQUAL_COUNTS: {
            // Partition primitives into equally-sized subsets
            mid = (start + end) / 2;
            std::nth_element(&buildData[start], &buildData[mid],
                             &buildData[end-1]+1, ComparePoints(dim));
            break;
        }
        case SPLIT_SAH: default: {
            // Partition primitives using approximate SAH
            if (nPrimitives <= 4) {
                // Partition primitives into equally-sized subsets
                mid = (start + end) / 2;
                std::nth_element(&buildData[start], &buildData[mid],
                                 &buildData[end-1]+1, ComparePoints(dim));
            }
            else {
                // Allocate _BucketInfo_ for SAH partition buckets
                const int nBuckets = 12;
                struct BucketInfo {
                    BucketInfo() { count = 0; }
                    int count;
                    BBox bounds;
                };
                BucketInfo buckets[nBuckets];

                // Initialize _BucketInfo_ for SAH partition buckets
                for (uint32_t i = start; i < end; ++i) {
                    int b = nBuckets *
                        ((buildData[i].centroid[dim] - centroidBounds.pMin[dim]) /
                         (centroidBounds.pMax[dim] - centroidBounds.pMin[dim]));
                    if (b == nBuckets) b = nBuckets-1;
                    Assert(b >= 0 && b < nBuckets);
                    buckets[b].count++;
                    buckets[b].bounds = Union(buckets[b].bounds, buildData[i].bounds);
                }

                // Compute costs for splitting after each bucket
                float cost[nBuckets-1];
                for (int i = 0; i < nBuckets-1; ++i) {
                    BBox 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] = .125f + (count0*b0.SurfaceArea() + count1*b1.SurfaceArea()) /
                              bbox.SurfaceArea();
                }

                // Find bucket to split at that minimizes SAH metric
                float minCost = cost[0];
                uint32_t minCostSplit = 0;
                for (int i = 1; i < nBuckets-1; ++i) {
                    if (cost[i] < minCost) {
                        minCost = cost[i];
                        minCostSplit = i;
                    }
                }

                // Either create leaf or split primitives at selected SAH bucket
                if (nPrimitives > maxPrimsInNode ||
                    minCost < nPrimitives) {
                    BVHPrimitiveInfo *pmid = std::partition(&buildData[start],
                        &buildData[end-1]+1,
                        CompareToBucket(minCostSplit, nBuckets, dim, centroidBounds));
                    mid = pmid - &buildData[0];
                }
                
                else {
                    // Create leaf _BVHBuildNode_
                    uint32_t firstPrimOffset = orderedPrims.size();
                    for (uint32_t i = start; i < end; ++i) {
                        uint32_t primNum = buildData[i].primitiveNumber;
                        orderedPrims.push_back(primitives[primNum]);
                    }
                    node->InitLeaf(firstPrimOffset, nPrimitives, bbox);
                    return node;
                }
            }
            break;
        }
        }
        node->InitInterior(dim,
                           recursiveBuild(buildArena, buildData, start, mid,
                                          totalNodes, orderedPrims),
                           recursiveBuild(buildArena, buildData, mid, end,
                                          totalNodes, orderedPrims));
    }
    return node;
}

void BVHAccel::FindBestSplit(std::vector<BVHAccelTreeNode *> &list, unsigned int begin, unsigned int end, float *splitValue, unsigned int *bestAxis) {
	if (end - begin == 2) {
		// Trivial case with two elements
		*splitValue = (list[begin]->bbox.pMax[0] + list[begin]->bbox.pMin[0] +
				list[end - 1]->bbox.pMax[0] + list[end - 1]->bbox.pMin[0]) / 2;
		*bestAxis = 0;
	} else {
		// Calculate BBs mean center (times 2)
		Point mean2(0, 0, 0), var(0, 0, 0);
		for (unsigned int i = begin; i < end; i++)
			mean2 += list[i]->bbox.pMax + list[i]->bbox.pMin;
		mean2 /= static_cast<float>(end - begin);

		// Calculate variance
		for (unsigned int i = begin; i < end; i++) {
			Vector v = list[i]->bbox.pMax + list[i]->bbox.pMin - mean2;
			v.x *= v.x;
			v.y *= v.y;
			v.z *= v.z;
			var += v;
		}
		// Select axis with more variance
		if (var.x > var.y && var.x > var.z)
			*bestAxis = 0;
		else if (var.y > var.z)
			*bestAxis = 1;
		else
			*bestAxis = 2;

		if (costSamples > 1) {
			BBox nodeBounds;
			for (unsigned int i = begin; i < end; i++)
				nodeBounds = Union(nodeBounds, list[i]->bbox);

			Vector d = nodeBounds.pMax - nodeBounds.pMin;
			const float invTotalSA = 1.f / nodeBounds.SurfaceArea();

			// Sample cost for split at some points
			float increment = 2 * d[*bestAxis] / (costSamples + 1);
			float bestCost = INFINITY;
			for (float splitVal = 2 * nodeBounds.pMin[*bestAxis] + increment; splitVal < 2 * nodeBounds.pMax[*bestAxis]; splitVal += increment) {
				int nBelow = 0, nAbove = 0;
				BBox bbBelow, bbAbove;
				for (unsigned int j = begin; j < end; j++) {
					if ((list[j]->bbox.pMax[*bestAxis] + list[j]->bbox.pMin[*bestAxis]) < splitVal) {
						nBelow++;
						bbBelow = Union(bbBelow, list[j]->bbox);
					} else {
						nAbove++;
						bbAbove = Union(bbAbove, list[j]->bbox);
					}
				}
				const float pBelow = bbBelow.SurfaceArea() * invTotalSA;
				const float pAbove = bbAbove.SurfaceArea() * invTotalSA;
				float eb = (nAbove == 0 || nBelow == 0) ? emptyBonus : 0.f;
				float cost = traversalCost + isectCost * (1.f - eb) * (pBelow * nBelow + pAbove * nAbove);
				// Update best split if this is lowest cost so far
				if (cost < bestCost) {
					bestCost = cost;
					*splitValue = splitVal;
				}
			}
		} else {
			// Split in half around the mean center
			*splitValue = mean2[*bestAxis];
		}
	}
}

void KDTree::BuildTree(int nodeNum, const BBox &nodeBounds,
        const vector<BBox> &allPrimBounds, uint32_t *primNums,
        int nPrimitives, int depth, BoundEdge *edges[3],
        uint32_t *prims0, uint32_t *prims1, int badRefines) {
    Assert(nodeNum == nextFreeNode);
	totalNodes++;
    // Get next free node from _nodes_ array
    if (nextFreeNode == nAllocedNodes) 
	{
        int nAlloc = max(2 * nAllocedNodes, 512);
        KdNode *n = AllocAligned<KdNode>(nAlloc);
        if (nAllocedNodes > 0) 
		{
            memcpy(n, nodes, nAllocedNodes * sizeof(KdNode));
            FreeAligned(nodes);
        }
        nodes = n;
        nAllocedNodes = nAlloc;
    }
    ++nextFreeNode;

    // Initialize leaf node if termination criteria met

	if (nPrimitives <= maxPrims || depth == 0) {
       
        nodes[nodeNum].initLeaf(primNums, nPrimitives, arena);
        return;
    }

    // Initialize interior node and continue recursion

    // Choose split axis position for interior node
    int bestAxis = -1, bestOffset = -1;
    float bestCost = INFINITY;
    float oldCost = isectCost * float(nPrimitives);
    float totalSA = nodeBounds.SurfaceArea();
    float invTotalSA = 1.f / totalSA;
    Vector d = nodeBounds.pMax - nodeBounds.pMin;

    // Choose which axis to split along
    uint32_t axis = nodeBounds.MaximumExtent();
    int retries = 0;
    bool retrySplit = true;
	while (retrySplit) {
		retrySplit=false;
		// Initialize edges for _axis_
		for (int i = 0; i < nPrimitives; ++i) {
			int pn = primNums[i];
		    const BBox &bbox = allPrimBounds[pn];
		 edges[axis][2*i] =   BoundEdge(bbox.pMin[axis], pn, true);
		 edges[axis][2*i+1] = BoundEdge(bbox.pMax[axis], pn, false);
		}
		sort(&edges[axis][0], &edges[axis][2*nPrimitives]);

		// Compute cost of all splits for _axis_ to find best
		int nBelow = 0, nAbove = nPrimitives;
		for (int i = 0; i < 2*nPrimitives; ++i) {
		 if (edges[axis][i].type == BoundEdge::END) --nAbove;
		 float edget = edges[axis][i].t;
		 if (edget > nodeBounds.pMin[axis] &&
		     edget < nodeBounds.pMax[axis]) {
				// Compute cost for split at _i_th edge
				uint32_t otherAxis0 = (axis + 1) % 3, otherAxis1 = (axis + 2) % 3;
				float belowSA = 2 * (d[otherAxis0] * d[otherAxis1] +
                                 (edget - nodeBounds.pMin[axis]) *
                                 (d[otherAxis0] + d[otherAxis1]));
				float aboveSA = 2 * (d[otherAxis0] * d[otherAxis1] +
                                 (nodeBounds.pMax[axis] - edget) *
                                 (d[otherAxis0] + d[otherAxis1]));
				float pBelow = belowSA * invTotalSA;
				float pAbove = aboveSA * invTotalSA;
				float eb = (nAbove == 0 || nBelow == 0) ? emptyBonus : 0.f;
				float cost = traversalCost +
                         isectCost * (1.f - eb) * (pBelow * nBelow + pAbove * nAbove);

				// Update best split if this is lowest cost so far
				if (cost < bestCost)  {
					bestCost = cost;
					bestAxis = axis;
					bestOffset = i;
				}
         }
         if (edges[axis][i].type == BoundEdge::START) ++nBelow;
		}
		Assert(nBelow == nPrimitives && nAbove == 0);

		// Create leaf if no good splits were found
		 if (bestAxis == -1 && retries < 2) {
			++retries;
			axis = (axis+1) % 3;
			retrySplit=true;
		}
	}
    if (bestCost > oldCost) ++badRefines;
    if ((bestCost > 4.f * oldCost && nPrimitives < 16) ||
        bestAxis == -1 || badRefines == 3) {
        nodes[nodeNum].initLeaf(primNums, nPrimitives, arena);
        return;
    }

//.........这里部分代码省略.........