Skip to content

Instantly share code, notes, and snippets.

@jamornsriwasansak
Last active June 20, 2021 21:02
Show Gist options
  • Save jamornsriwasansak/1875d25ca66e74f25a7e753b660557a1 to your computer and use it in GitHub Desktop.
Save jamornsriwasansak/1875d25ca66e74f25a7e753b660557a1 to your computer and use it in GitHub Desktop.
#include "lightcuts.hpp"
template <class PointLightType, class PointLightNodeType>
class PointLightTree
{
public:
Spectrum samplePosition(Vec3 * wi, bool * isVisible, const Intersection & isect, const Scene & scene, const Sampler & sampler) const { assert(false && "unimplemented"); return Spectrum(0.0); };
Spectrum sampleDirection(Vec3 * wi, Vec3 * position, const Sampler & sampler) const { assert(false && "unimplemented"); return Spectrum(0.0); }
Spectrum eval(const Intersection & isect, const Vec3 & wo, const Scene & scene) const
{
stat.cutSize = 0;
return this->evalBranchExhaust(mRoot, isect, wo, scene);
}
Spectrum evalBranchExhaust(const shared_ptr<PointLightNodeType> & currentNode, const Intersection & isect, const Vec3 & wo, const Scene & scene) const
{
assert(currentNode != nullptr);
stat.cutSize++;
if (currentNode->isLeaf())
{
stat.geometryTermEval++;
return currentNode->eval(isect, wo, scene);
}
else
{
return this->evalBranchExhaust(currentNode->mChild[0], isect, wo, scene) + this->evalBranchExhaust(currentNode->mChild[1], isect, wo, scene);
}
}
Spectrum evalBranch(const shared_ptr<PointLightNodeType> & currentNode, const Intersection & isect, const Vec3 & wo, const Scene & scene) const
{
assert(currentNode != nullptr);
struct HeapNode
{
const PointLightNodeType * node;
Float error;
Float gvTerm;
Spectrum contrib;
HeapNode() {}
};
auto CompareLessHeapNode = [] (const HeapNode & a, const HeapNode & b) { return a.error < b.error; };
std::priority_queue<HeapNode, std::vector<HeapNode>, decltype(CompareLessHeapNode)> heap(CompareLessHeapNode);
uint32_t cutSize = 0;
HeapNode rootHeapNode;
rootHeapNode.node = currentNode->get();
rootHeapNode.gvTerm = currentNode->geometryVisibility(isect, scene);
rootHeapNode.contrib = currentNode->eval(rootHeapNode.gvTerm, isect, wo);
rootHeapNode.error = currentNode->maxError(isect, wo);
while (heap.size() != 0 && cutSize < 1000)
{
const HeapNode top = heap.top();
}
}
void build(const std::vector<shared_ptr<PointLightType>> & pointLights, const Sampler & sampler)
{
const size_t numLights = pointLights.size();
if (numLights == 0) { return; }
// compute constant c that control relative scaling
Float c2 = 0.0f;
{
Aabb bound;
for (const shared_ptr<PointLightType> & light : pointLights)
{
bound = Aabb::Union(light->mPosition, bound);
}
c2 = Aabb::DiagonalLength2(bound);
}
// populate active nodes
std::vector<shared_ptr<PointLightNodeType>> activeNodes(numLights);
for (size_t i = 0; i < numLights; i++)
{
const shared_ptr<PointLightType> & rp = pointLights[i];
activeNodes[i] = make_shared<PointLightNodeType>(*rp);
}
/// TODO:: wait for @kobayashi to unify oriented light and omni light kdtree and
while (activeNodes.size() != 1)
{
const size_t numActiveNodes = activeNodes.size();
std::pair<size_t, size_t> minPair;
Float minDistance = std::numeric_limits<Float>::infinity();
for (size_t ia = 0;ia < numActiveNodes - 1;ia++)
{
const shared_ptr<PointLightNodeType> & rna = activeNodes[ia];
for (size_t ib = ia + 1;ib < numActiveNodes;ib++)
{
const shared_ptr<PointLightNodeType> & rnb = activeNodes[ib];
Float distance = rna->cost(*rnb, c2);
if (distance < minDistance)
{
minDistance = distance;
minPair.first = ia;
minPair.second = ib;
}
}
}
shared_ptr<PointLightNodeType> & minRnA = activeNodes[minPair.first];
shared_ptr<PointLightNodeType> & minRnB = activeNodes[minPair.second];
shared_ptr<PointLightNodeType> newRn = make_shared<PointLightNodeType>(minRnA, minRnB, sampler);
std::swap(activeNodes[minPair.second], activeNodes.back());
activeNodes.pop_back();
std::swap(activeNodes[minPair.first], activeNodes.back());
activeNodes.pop_back();
activeNodes.push_back(newRn);
}
this->mRoot = activeNodes[0];
}
void buildFast(const std::vector<shared_ptr<PointLightType>> & pointLights)
{
Rng rng(81521);
const size_t numPointLights = pointLights.size();
if (numPointLights == 0) { return; }
// compute constant c that control relative scaling
Float diagonal2 = (Float)0.0;
{
Aabb bound;
for (const shared_ptr<PointLightType> & lightPtr : pointLights)
{
bound = Aabb::Union(lightPtr->mPosition, bound);
}
diagonal2 = Aabb::DiagonalLength2(bound);
}
// populate active nodes
std::vector<const PointLightNodeType *> activeNodes(numPointLights);
{
for (size_t i = 0; i < numPointLights; i++)
{
const PointLightType & p = *pointLights[i];
activeNodes[i] = new PointLightNodeType(p);
}
}
// declare kdtree stuffs
using KDT = LcKdTree<PointLightNodeType>;
struct QueueElem
{
struct NodePair
{
shared_ptr<typename KDT::LcKdTreeNode> kdtnode;
const PointLightNodeType *lightNode;
bool obsoleted() const
{
if (kdtnode->deleted()) { return true; }
if (kdtnode->lightNode != lightNode) { return true; }
return false;
}
};
QueueElem(NodePair first, NodePair second, float distance):
first(first), second(second), distance(distance)
{
}
QueueElem(const QueueElem&) = default;
NodePair first, second;
float distance;
bool operator >(const QueueElem &rhs) const { return distance > rhs.distance; }
};
using Queue = std::priority_queue<QueueElem, std::vector<QueueElem>, std::greater<QueueElem>>;
Queue pairs;
KDT kdt(activeNodes, diagonal2);
auto pushNearestToQueue = [&] (shared_ptr<KDT::LcKdTreeNode> p)
{
auto nearest = kdt.findNearest(p);
if (nearest.first)
{
pairs.emplace(
QueueElem::NodePair{p, p->lightNode},
QueueElem::NodePair{nearest.first, nearest.first->lightNode},
nearest.second);
}
};
for (auto p = kdt.begin(); p; p = kdt.next(p))
{
pushNearestToQueue(p);
}
while (!pairs.empty())
{
auto pair = pairs.top(); pairs.pop();
if (pair.first.obsoleted() && pair.second.obsoleted())
{
continue;
}
else if (pair.first.obsoleted())
{
pushNearestToQueue(pair.second.kdtnode);
}
else if (pair.second.obsoleted())
{
pushNearestToQueue(pair.first.kdtnode);
}
else
{
const PointLightNodeType * minPlnA = pair.first.kdtnode->lightNode;
const PointLightNodeType * minPlnB = pair.second.kdtnode->lightNode;
const PointLightNodeType * newPln = PointLightNodeType::MergeNode(minPlnA, minPlnB, rng);
auto remain_node = pair.first.kdtnode;
auto delete_node = pair.second.kdtnode;
if (remain_node->depth > delete_node->depth)
{
std::swap(remain_node, delete_node);
}
remain_node->lightNode = newPln;
delete_node->lightNode = nullptr;
kdt.cleanup(delete_node);
if (!pairs.empty())
{
if (!kdt.findNearest(remain_node).first)
{
kdt.findNearest(remain_node);
}
}
pushNearestToQueue(remain_node);
}
}
auto begin = kdt.begin();
if (!begin)
{
throw std::runtime_error("unreachable! empty tree");
}
assert(begin);
assert(begin->lightNode);
this->mRoot = begin->lightNode;
}
struct
{
mutable std::atomic<int> cutSize = 0;
mutable std::atomic<int> geometryTermEval = 0;
} stat;
protected:
shared_ptr<PointLightNodeType> mRoot;
};
#pragma once
#include <atomic>
#include <queue>
#include <functional>
#include "common\reflectcuts.h"
#include "common\light.h"
#include "common\rng.h"
#include "lightaccels\util\lckdtree.h"
template <typename DerivedType>
class PointLightNode
{
public:
// must be symmetric a.cost(b) = b.cost(a)
virtual Float cost(const DerivedType & b, const Float c2) const = 0;
// for evaluating lightcuts
virtual Float maxError(const Intersection & isect, const Vec3 & wo) const = 0;
virtual Float geometryVisibility(const Intersection & isect, const Scene & scene) const = 0;
// return mRepresent->eval
virtual Spectrum eval(const Intersection & isect, const Vec3 & wo, const Scene & scene) const = 0;
// same as eval but don't need to recompute G term
virtual Spectrum eval(const Float computedG, const Intersection & isect, const Vec3 & wo) const = 0;
// is leaf
virtual bool isLeaf() const = 0;
};
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment