Skip to content

Commit

Permalink
a better importance sampling for material micrograin
Browse files Browse the repository at this point in the history
  • Loading branch information
Cchen-77 committed Jul 27, 2024
1 parent 2b491d3 commit 3c79793
Show file tree
Hide file tree
Showing 7 changed files with 140 additions and 138 deletions.
14 changes: 5 additions & 9 deletions scenes/micrograin-teapot/scene.json
Original file line number Diff line number Diff line change
Expand Up @@ -6,17 +6,13 @@
"type":"porousLayerMicrograin",
"bulkMaterial":
{
"name":"conductor",
"type":"conductor",
"albedo":"textures/rustyMetal.jpg",
"eta":[2,2,2],
"k":[0,0,0],
"roughness":0.3,
"distribution":"ggx"
"name":"plastic",
"albedo":0.9,
"type":"plastic"
},
"micrograinType":"conductor",
"tau0":0.3,
"beta":0.8,
"tau0":0.27,
"beta":0.1,
"R0":[0.588,0.294,0.001],
"k":[0,0,0]
},
Expand Down
67 changes: 32 additions & 35 deletions src/FunctionLayer/Material/BxDF/ConductorBxDF.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2,80 +2,77 @@
#include "Fresnel.h"
#include "CoreLayer/Geometry/Frame.h"


ConductorBxDF::ConductorBxDF(Vec3d _eta,Vec3d _k, Spectrum _albedo)
:eta(_eta),k(_k),albedo(_albedo){

ConductorBxDF::ConductorBxDF(Vec3d _eta, Vec3d _k, Spectrum _albedo)
: eta(_eta), k(_k), albedo(_albedo) {
}


Spectrum ConductorBxDF::f(const Vec3d & out, const Vec3d & in) const {
Spectrum ConductorBxDF::f(const Vec3d &out, const Vec3d &in) const {
return {};
}

double ConductorBxDF::pdf(const Vec3d & out, const Vec3d & in) const {
double ConductorBxDF::pdf(const Vec3d &out, const Vec3d &in) const {
return 0;
}

BxDFSampleResult ConductorBxDF::sample(const Vec3d & out, const Point2d & sample) const {
BxDFSampleResult ConductorBxDF::sample(const Vec3d &out, const Point2d &sample) const {
BxDFSampleResult result;
if(out.z<0)
if (out.z < 0)
return result;
result.directionIn = Frame::reflect(out);
result.pdf = 1;
auto F = Fresnel::conductorReflectance(eta, k, out.z);
result.s = albedo * F / (abs(result.directionIn.z));
result.s = albedo * F / (abs(result.directionIn.z));
result.bxdfSampleType = BXDFType(BXDF_REFLECTION | BXDF_SPECULAR);
return result;
}

RoughConductorBxDF::RoughConductorBxDF(Vec3d _eta, Vec3d _k, Spectrum _albedo, double _uRoughness, double _vRoughness,
std::shared_ptr < MicrofacetDistribution > _distrib) :
eta(_eta),k(_k),albedo(_albedo),distrib(std::move(_distrib))
{
alphaXY = Vec2d(distrib->roughnessToAlpha(_uRoughness),distrib->roughnessToAlpha(_vRoughness));
}
std::shared_ptr<MicrofacetDistribution> _distrib) : eta(_eta), k(_k), albedo(_albedo), distrib(std::move(_distrib)) {
alphaXY = Vec2d(distrib->roughnessToAlpha(_uRoughness), distrib->roughnessToAlpha(_vRoughness));
}

Spectrum RoughConductorBxDF::f(const Vec3d & out, const Vec3d & in) const {
if(out.z<0 || in.z<0)
Spectrum RoughConductorBxDF::f(const Vec3d &out, const Vec3d &in) const {
if (out.z < 0 || in.z < 0)
return {0.};
double cosThetaO = AbsCosTheta(out), cosThetaI = AbsCosTheta(in);
Vec3d wh = in + out;
// Handle degenerate cases for microfacet reflection
if (cosThetaI == 0 || cosThetaO == 0) return {0.};
if (wh.x == 0 && wh.y == 0 && wh.z == 0) return {0.};
wh = normalize(wh);
double cosI = dot(out,wh.z>0?wh:-wh);
auto F = Fresnel::conductorReflectance(eta,k,cosI);
return albedo * (distrib->D(wh,alphaXY) * F * distrib->G(out, in,alphaXY)/ (4 * cosThetaI * cosThetaO));
if (dot(wh, out) < 0) {
return 0.;
}
double cosI = dot(out, wh.z > 0 ? wh : -wh);
auto F = Fresnel::conductorReflectance(eta, k, cosI);
return albedo * (distrib->D(wh, alphaXY) * F * distrib->G(out, in, alphaXY) / (4 * cosThetaI * cosThetaO));
}

double RoughConductorBxDF::pdf(const Vec3d & out, const Vec3d & in) const {
if ( CosTheta(out) < 0 || CosTheta(in) < 0) return 0;
double RoughConductorBxDF::pdf(const Vec3d &out, const Vec3d &in) const {
if (CosTheta(out) < 0 || CosTheta(in) < 0) return 0;
Vec3d wh = normalize(out + in);
return distrib->Pdf(out, wh,alphaXY) / (4 * dot(out, wh));
return distrib->Pdf(out, wh, alphaXY) / (4 * dot(out, wh));
}

BxDFSampleResult RoughConductorBxDF::sample(const Vec3d & out, const Point2d & sample) const {
BxDFSampleResult RoughConductorBxDF::sample(const Vec3d &out, const Point2d &sample) const {
BxDFSampleResult result{};
if(CosTheta(out)<0){
if (CosTheta(out) < 0) {
return result;
}
Vec3d wh = distrib->Sample_wh(out,sample,alphaXY);
Vec3d in = Frame::reflect(out,wh);
if (dot(wh,out) < 0 || CosTheta(in)<0)
{
result.s=Spectrum(0);
Vec3d wh = distrib->Sample_wh(out, sample, alphaXY);
Vec3d in = Frame::reflect(out, wh);
if (dot(wh, out) < 0 || CosTheta(in) < 0) {
result.s = Spectrum(0);
result.directionIn = in;
return result;
}
double mPdf =distrib->Pdf(out,wh,alphaXY);
auto woDoth = dot(out,wh);
double pdf = mPdf*0.25f/woDoth;
double mPdf = distrib->Pdf(out, wh, alphaXY);
auto woDoth = dot(out, wh);
double pdf = mPdf * 0.25f / woDoth;
Spectrum F = Fresnel::conductorReflectance(eta, k, woDoth);
result.pdf = pdf;
result.s = ( F * distrib->G(out,in,alphaXY)*distrib->D(wh,alphaXY))/abs(4 * out.z * in.z);
result.s = (F * distrib->G(out, in, alphaXY) * distrib->D(wh, alphaXY)) / abs(4 * out.z * in.z);
result.directionIn = in;
result.bxdfSampleType = BXDFType(BXDF_REFLECTION | BXDF_GLOSSY);
return result;
}

135 changes: 65 additions & 70 deletions src/FunctionLayer/Material/BxDF/GlintBxDF.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2,40 +2,38 @@
#include "Fresnel.h"
#include "CoreLayer/Geometry/Frame.h"

GlintBxDF::GlintBxDF(Vec3d _eta, Vec3d _k, Spectrum _albedo, double _uRoughness, double _vRoughness,Intersection _it,
std::shared_ptr < MicrofacetDistribution > _distrib) :
eta(_eta),k(_k),albedo(_albedo),it(_it),distrib(std::move(_distrib))
{
alphaXY = Vec2d(distrib->roughnessToAlpha(_uRoughness),distrib->roughnessToAlpha(_vRoughness));
}
GlintBxDF::GlintBxDF(Vec3d _eta, Vec3d _k, Spectrum _albedo, double _uRoughness, double _vRoughness, Intersection _it,
std::shared_ptr<MicrofacetDistribution> _distrib) : eta(_eta), k(_k), albedo(_albedo), it(_it), distrib(std::move(_distrib)) {
alphaXY = Vec2d(distrib->roughnessToAlpha(_uRoughness), distrib->roughnessToAlpha(_vRoughness));
}

BxDFSampleResult GlintBxDF::sample(const Vec3d & out, const Point2d & sample) const {
BxDFSampleResult GlintBxDF::sample(const Vec3d &out, const Point2d &sample) const {
BxDFSampleResult result{};
if(CosTheta(out)<0){
if (CosTheta(out) < 0) {
return result;
}
Vec3d wh = distrib->Sample_wh(out,sample,alphaXY);
Vec3d in = Frame::reflect(out,wh);
if (dot(wh,out) < 0 || CosTheta(in)<0)
{
result.s=Spectrum(0);
Vec3d wh = distrib->Sample_wh(out, sample, alphaXY);
Vec3d in = Frame::reflect(out, wh);
if (dot(wh, out) < 0 || CosTheta(in) < 0) {
result.s = Spectrum(0);
result.directionIn = in;
return result;
}
double mPdf =distrib->Pdf(out,wh,alphaXY);
auto woDoth = dot(out,wh);
double pdf = mPdf*0.25f/woDoth;
double mPdf = distrib->Pdf(out, wh, alphaXY);
auto woDoth = dot(out, wh);
double pdf = mPdf * 0.25f / woDoth;
Spectrum F = Fresnel::conductorReflectance(eta, k, woDoth);
result.pdf = pdf;
result.s = ( F * distrib->G(out,in,alphaXY)*distrib->D(wh,alphaXY))/abs(4 * out.z * in.z);
result.s = (F * distrib->G(out, in, alphaXY) * distrib->D(wh, alphaXY)) / abs(4 * out.z * in.z);
result.directionIn = in;
result.bxdfSampleType = BXDFType(BXDF_REFLECTION | BXDF_GLOSSY);
return result;
}

double GlintBxDF::pdf(const Vec3d & out, const Vec3d & in) const {
if ( CosTheta(out) < 0 || CosTheta(in) < 0) return 0;
double GlintBxDF::pdf(const Vec3d &out, const Vec3d &in) const {
if (CosTheta(out) < 0 || CosTheta(in) < 0) return 0;
Vec3d wh = normalize(out + in);
return distrib->Pdf(out, wh,alphaXY) / (4 * dot(out, wh));
return distrib->Pdf(out, wh, alphaXY) / (4 * dot(out, wh));
}

// count particles number fall in tex area
Expand All @@ -46,120 +44,117 @@ double GlintBxDF::count_spatial(BoundingBox2d &queryBox) const {
SpatialNode node;
BoundingBox2d clipBox;
double p, w;
//bfs
while (!queue.empty()){
// bfs
while (!queue.empty()) {
node = queue.back();
queue.pop_back();
if(!queryBox.overlap(node.first) || node.second <= 0){
if (!queryBox.overlap(node.first) || node.second <= 0) {
continue;
}
else if(queryBox.contains(node.first)){
} else if (queryBox.contains(node.first)) {
countNum += node.second;
}
else{
} else {
p = countNum / SPATIAL_SAMPLE_NUM;
// approximation
if(sqrt(node.second * p * (1 - p)) < EPSILON * countNum){
if (sqrt(node.second * p * (1 - p)) < EPSILON * countNum) {
clipBox = node.first;
clipBox = BoundingBoxUnionIntersect(queryBox, clipBox);
countNum += (uint32_t)(node.second * clipBox.getSurfaceArea() / node.first.getSurfaceArea());
}
else{
if(queryDepth >= MAX_QUERY_DEPTH){
} else {
if (queryDepth >= MAX_QUERY_DEPTH) {
break;
}
// split box
w = 1 - queryPs[queryDepth][0] - queryPs[queryDepth][1] - queryPs[queryDepth][2];
Point2d center = node.first.getCenter(), pMin = node.first.pMin, pMax = node.first.pMax;
queue.push_back(std::make_pair(
BoundingBox2(pMin, center),
BoundingBox2(pMin, center),
(uint32_t)(node.second * queryPs[queryDepth][0])));
queue.push_back(std::make_pair(
BoundingBox2(Point2d(pMin[0], center[1]), Point2d(center[0],pMax[1])),
BoundingBox2(Point2d(pMin[0], center[1]), Point2d(center[0], pMax[1])),
(uint32_t)(node.second * queryPs[queryDepth][1])));
queue.push_back(std::make_pair(
BoundingBox2(Point2d(center[0], pMin[1]), Point2d(pMax[0],center[1])),
BoundingBox2(Point2d(center[0], pMin[1]), Point2d(pMax[0], center[1])),
(uint32_t)(node.second * queryPs[queryDepth][2])));
queue.push_back(std::make_pair(
BoundingBox2(center, pMax),
BoundingBox2(center, pMax),
(uint32_t)(node.second * w)));
queryDepth++;
}
}
}
assert(countNum < SPATIAL_SAMPLE_NUM);
return (double)countNum / SPATIAL_SAMPLE_NUM;


}

// count particles number fall in cone
double GlintBxDF::count_direction(const Vec3d &wi, const Vec3d & wo) const{
double GlintBxDF::count_direction(const Vec3d &wi, const Vec3d &wo) const {
uint32_t countNum = 0, queryDepth = 1;
std::vector<DirectionNode> queue;
ConicQuery query(wi, wo);
queue.push_back(std::make_pair(DirTriangle(Point2d(0.f, 0.f), Point2d(M_PI / 2, 0.f), Point2d(0.f, M_PI /2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * queryPs[0].x)));
queue.push_back(std::make_pair(DirTriangle(Point2d(M_PI / 2, 0.f), Point2d(M_PI, 0.f), Point2d(0.f, M_PI /2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * queryPs[0].y)));
queue.push_back(std::make_pair(DirTriangle(Point2d(M_PI , 0.f), Point2d(M_PI * 3 / 2, 0.f), Point2d(0.f, M_PI /2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * queryPs[0].z)));
queue.push_back(std::make_pair(DirTriangle(Point2d(M_PI * 3 / 2, 0.f), Point2d(M_PI * 2, 0.f), Point2d(0.f, M_PI /2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * (1 - queryPs[0].x - queryPs[0].y - queryPs[0].z))));
queue.push_back(std::make_pair(DirTriangle(Point2d(0.f, 0.f), Point2d(M_PI / 2, 0.f), Point2d(0.f, M_PI / 2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * queryPs[0].x)));
queue.push_back(std::make_pair(DirTriangle(Point2d(M_PI / 2, 0.f), Point2d(M_PI, 0.f), Point2d(0.f, M_PI / 2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * queryPs[0].y)));
queue.push_back(std::make_pair(DirTriangle(Point2d(M_PI, 0.f), Point2d(M_PI * 3 / 2, 0.f), Point2d(0.f, M_PI / 2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * queryPs[0].z)));
queue.push_back(std::make_pair(DirTriangle(Point2d(M_PI * 3 / 2, 0.f), Point2d(M_PI * 2, 0.f), Point2d(0.f, M_PI / 2)),
(uint32_t)(DIRECTION_SAMPLE_NUM * (1 - queryPs[0].x - queryPs[0].y - queryPs[0].z))));
DirectionNode node;
double p;
//bfs
while (!queue.empty()){
// bfs
while (!queue.empty()) {
node = queue.back();
queue.pop_back();
// not overlap
if((!(query.inConic(node.first[0]) || query.inConic(node.first[1]) ||
query.inConic(node.first[2]))) || node.second <= 0){
if ((!(query.inConic(node.first[0]) || query.inConic(node.first[1]) ||
query.inConic(node.first[2]))) ||
node.second <= 0) {
continue;
}
// contains
else if(query.conicContain(node.first)){
else if (query.conicContain(node.first)) {
countNum += node.second;
}
else{
} else {
p = countNum / SPATIAL_SAMPLE_NUM;
// approximation
// todo: calculate the ratio of area
if(sqrt(node.second * p * (1 - p)) < EPSILON * countNum){
if (sqrt(node.second * p * (1 - p)) < EPSILON * countNum) {
countNum += node.second;
}
else{
if(queryDepth >= MAX_QUERY_DEPTH){
} else {
if (queryDepth >= MAX_QUERY_DEPTH) {
break;
}
// split triangle
queue.push_back(std::make_pair(DirTriangle(node.first.edgeMidPoint(0), node.first.edgeMidPoint(1), node.first.edgeMidPoint(2)),
(uint32_t)(node.second * queryPs[queryDepth].x)));
(uint32_t)(node.second * queryPs[queryDepth].x)));
queue.push_back(std::make_pair(DirTriangle(node.first.edgeMidPoint(0), node.first.edgeMidPoint(2), node.first[0]),
(uint32_t)(node.second * queryPs[queryDepth].y)));
(uint32_t)(node.second * queryPs[queryDepth].y)));
queue.push_back(std::make_pair(DirTriangle(node.first.edgeMidPoint(0), node.first.edgeMidPoint(1), node.first[1]),
(uint32_t)(node.second * queryPs[queryDepth].z)));
(uint32_t)(node.second * queryPs[queryDepth].z)));
queue.push_back(std::make_pair(DirTriangle(node.first.edgeMidPoint(1), node.first.edgeMidPoint(2), node.first[2]),
(uint32_t)(node.second * (1- queryPs[0].x - queryPs[0].y - queryPs[0].z))));
queryDepth++;
(uint32_t)(node.second * (1 - queryPs[0].x - queryPs[0].y - queryPs[0].z))));
queryDepth++;
}
}
}
return (double)countNum / DIRECTION_SAMPLE_NUM;
}

Spectrum GlintBxDF::f(const Vec3d & out, const Vec3d & in) const {
if(out.z<0 || in.z<0)
Spectrum GlintBxDF::f(const Vec3d &out, const Vec3d &in) const {
if (out.z < 0 || in.z < 0)
return {0.};
double cosThetaO = AbsCosTheta(out), cosThetaI = AbsCosTheta(in);
Vec3d wh = in + out;
// Handle degenerate cases for microfacet reflection
if (cosThetaI == 0 || cosThetaO == 0) return {0.};
if (wh.x == 0 && wh.y == 0 && wh.z == 0) return {0.};
wh = normalize(wh);
double cosI = dot(out,wh.z>0?wh:-wh);
auto F = Fresnel::conductorReflectance(eta,k,cosI);
auto G = distrib->G(out, in,alphaXY);
if (dot(wh, out) < 0) {
return 0.;
}
double cosI = dot(out, wh.z > 0 ? wh : -wh);
auto F = Fresnel::conductorReflectance(eta, k, cosI);
auto G = distrib->G(out, in, alphaXY);

double du = (it.dudx + it.dudy) / 2, dv = (it.dvdx + it.dvdy) / 2;
BoundingBox2d texBox;
Expand All @@ -170,8 +165,8 @@ Spectrum GlintBxDF::f(const Vec3d & out, const Vec3d & in) const {

double D = count_direction(in, out);
D *= count_spatial(texBox);
if(D < 1e-8)
if (D < 1e-8)
return {0.};
return albedo * (dot(in,wh) * F * D * G) /
(texBox.getSurfaceArea() * M_PI * (1 - cos(GAMMA)) * cosThetaI * cosThetaO);
return albedo * (dot(in, wh) * F * D * G) /
(texBox.getSurfaceArea() * M_PI * (1 - cos(GAMMA)) * cosThetaI * cosThetaO);
}
Loading

0 comments on commit 3c79793

Please sign in to comment.