提交 4d2e5636 编写于 作者: M Matt Pharr

Update from book source. No functional changes.

上级 2bec4dc4
......@@ -116,11 +116,9 @@ pstd::optional<BSDFSample> DielectricBxDF::Sample_f(
} else {
// Sample rough dielectric BSDF
// Sample half-angle vector for outgoing direction and compute Fresnel factor
Vector3f wh = mfDistrib.Sample_wm(wo, u);
Float R = FrDielectric(Dot(wo, wh), eta);
Vector3f wm = mfDistrib.Sample_wm(wo, u);
Float R = FrDielectric(Dot(wo, wm), eta);
Float T = 1 - R;
// Compute probabilities _pr_ and _pt_ for sampling reflection and transmission
Float pr = R, pt = T;
if (!(sampleFlags & BxDFReflTransFlags::Reflection))
......@@ -131,27 +129,27 @@ pstd::optional<BSDFSample> DielectricBxDF::Sample_f(
return {};
if (uc < pr / (pr + pt)) {
// Sample reflection at non-delta dielectric interface
Vector3f wi = Reflect(wo, wh);
// Sample reflection at rough dielectric interface
Vector3f wi = Reflect(wo, wm);
if (!SameHemisphere(wo, wi))
return {};
// Compute PDF of direction $\wi$ for microfacet reflection
Float pdf = mfDistrib.PDF(wo, wh) / (4 * AbsDot(wo, wh)) * pr / (pr + pt);
Float pdf = mfDistrib.PDF(wo, wm) / (4 * AbsDot(wo, wm)) * pr / (pr + pt);
CHECK(!IsNaN(pdf));
// Evaluate BRDF and return _BSDFSample_ for dielectric microfacet reflection
Float cosTheta_o = CosTheta(wo), cosTheta_i = CosTheta(wi);
if (cosTheta_i == 0 || cosTheta_o == 0)
return {};
SampledSpectrum f(mfDistrib.D(wh) * mfDistrib.G(wo, wi) * R /
SampledSpectrum f(mfDistrib.D(wm) * mfDistrib.G(wo, wi) * R /
(4 * cosTheta_i * cosTheta_o));
return BSDFSample(f, wi, pdf, BxDFFlags::GlossyReflection);
} else {
// Sample transmission at non-delta dielectric interface
// Sample transmission at rough dielectric interface
Float etap;
Vector3f wi;
bool tir = !Refract(wo, (Normal3f)wh, eta, &etap, &wi);
bool tir = !Refract(wo, (Normal3f)wm, eta, &etap, &wi);
CHECK_RARE(1e-5f, tir);
if (SameHemisphere(wo, wi))
return {};
......@@ -159,17 +157,17 @@ pstd::optional<BSDFSample> DielectricBxDF::Sample_f(
return {};
// Evaluate BSDF
Float denom = Sqr(Dot(wo, wh) + etap * Dot(wi, wh));
Float denom = Sqr(Dot(wo, wm) + etap * Dot(wi, wm));
Float factor = (mode == TransportMode::Radiance) ? Sqr(1 / etap) : 1;
SampledSpectrum f(T * factor *
std::abs(mfDistrib.D(wh) * mfDistrib.G(wo, wi) *
Dot(wi, wh) * Dot(wo, wh) /
std::abs(mfDistrib.D(wm) * mfDistrib.G(wo, wi) *
Dot(wi, wm) * Dot(wo, wm) /
(CosTheta(wi) * CosTheta(wo) * denom)));
Float dwh_dwi =
/*Sqr(etap) * */ AbsDot(wi, wh) / denom;
Float pdf = mfDistrib.PDF(wo, wh) * dwh_dwi * pt / (pr + pt);
Float dwm_dwi =
/*Sqr(etap) * */ AbsDot(wi, wm) / denom;
Float pdf = mfDistrib.PDF(wo, wm) * dwm_dwi * pt / (pr + pt);
CHECK(!IsNaN(pdf));
return BSDFSample(f, wi, pdf, BxDFFlags::GlossyTransmission, etap);
......@@ -178,88 +176,86 @@ pstd::optional<BSDFSample> DielectricBxDF::Sample_f(
}
SampledSpectrum DielectricBxDF::f(Vector3f wo, Vector3f wi, TransportMode mode) const {
// Return zero for unsupported dielectric configurations
if (eta == 1 || mfDistrib.EffectivelySmooth())
return SampledSpectrum(0.f);
// Evaluate rough dielectric BSDF
// Compute generalized half vector
Float cosTheta_o = CosTheta(wo), cosTheta_i = CosTheta(wi);
if (cosTheta_i == 0 || cosTheta_o == 0)
return SampledSpectrum(0.f);
// Compute generalized half vector
return {};
bool reflect = cosTheta_i * cosTheta_o > 0;
float etap = 1;
if (!reflect)
etap = cosTheta_o > 0 ? eta : (1 / eta);
Vector3f wh = wi * etap + wo;
CHECK_RARE(1e-5f, LengthSquared(wh) == 0);
if (LengthSquared(wh) == 0)
Vector3f wm = wi * etap + wo;
CHECK_RARE(1e-5f, LengthSquared(wm) == 0);
if (LengthSquared(wm) == 0)
return {};
wh = FaceForward(Normalize(wh), Normal3f(0, 0, 1));
wm = FaceForward(Normalize(wm), Normal3f(0, 0, 1));
if (reflect) {
// Compute reflection at non-delta dielectric interface
Float F = FrDielectric(Dot(wi, wh), eta);
return SampledSpectrum(mfDistrib.D(wh) * mfDistrib.G(wo, wi) * F /
// Compute reflection at rough dielectric interface
Float F = FrDielectric(Dot(wi, wm), eta);
return SampledSpectrum(mfDistrib.D(wm) * mfDistrib.G(wo, wi) * F /
std::abs(4 * cosTheta_i * cosTheta_o));
} else {
// Compute transmission at non-delta dielectric interface
// Return no transmission if _wi_ and _wo_ are on the same side of _wh_
if (Dot(wi, wh) * Dot(wo, wh) > 0)
// Compute transmission at rough dielectric interface
// Return no transmission if _wi_ and _wo_ are on the same side of _wm_
if (Dot(wi, wm) * Dot(wo, wm) > 0)
return {};
// Evaluate BTDF for transmission through microfacet interface
Float F = FrDielectric(Dot(wo, wh), eta);
Float denom = Sqr(Dot(wo, wh) + etap * Dot(wi, wh));
Float F = FrDielectric(Dot(wo, wm), eta);
Float denom = Sqr(Dot(wo, wm) + etap * Dot(wi, wm));
Float factor = (mode == TransportMode::Radiance) ? Sqr(1 / etap) : 1;
return SampledSpectrum(
(1 - F) * factor *
std::abs(mfDistrib.D(wh) * mfDistrib.G(wo, wi) * Dot(wi, wh) * Dot(wo, wh) /
std::abs(mfDistrib.D(wm) * mfDistrib.G(wo, wi) * Dot(wi, wm) * Dot(wo, wm) /
(cosTheta_i * cosTheta_o * denom)));
}
}
Float DielectricBxDF::PDF(Vector3f wo, Vector3f wi, TransportMode mode,
BxDFReflTransFlags sampleFlags) const {
// Return zero PDF for unsupported dielectric configurations
if (eta == 1 || mfDistrib.EffectivelySmooth())
return 0.f;
// Evaluate sampling PDF of rough dielectric BSDF
// Compute generalized half vector
Float cosTheta_o = CosTheta(wo), cosTheta_i = CosTheta(wi);
if (cosTheta_i == 0 || cosTheta_o == 0)
return 0.f;
// Compute generalized half vector
return {};
bool reflect = cosTheta_i * cosTheta_o > 0;
float etap = 1;
if (!reflect)
etap = cosTheta_o > 0 ? eta : (1 / eta);
Vector3f wh = wi * etap + wo;
CHECK_RARE(1e-5f, LengthSquared(wh) == 0);
if (LengthSquared(wh) == 0)
Vector3f wm = wi * etap + wo;
CHECK_RARE(1e-5f, LengthSquared(wm) == 0);
if (LengthSquared(wm) == 0)
return {};
wh = FaceForward(Normalize(wh), Normal3f(0, 0, 1));
wm = FaceForward(Normalize(wm), Normal3f(0, 0, 1));
if (reflect) {
// Return PDF for non-delta dielectric reflection
// Return PDF for rough dielectric reflection
if (!(sampleFlags & BxDFReflTransFlags::Reflection))
return 0;
// Compute Fresnel factor and probabilities for dielectric reflection PDF
Float F = FrDielectric(Dot(wo, wh), eta);
Float F = FrDielectric(Dot(wo, wm), eta);
CHECK_RARE(1e-5f, F == 0);
Float pr = F, pt = 1 - F;
if (!(sampleFlags & BxDFReflTransFlags::Transmission))
pt = 0;
return mfDistrib.PDF(wo, wh) / (4 * AbsDot(wo, wh)) * pr / (pr + pt);
return mfDistrib.PDF(wo, wm) / (4 * AbsDot(wo, wm)) * pr / (pr + pt);
} else {
// Return PDF for non-delta dielectric transmission
// Return PDF for rough dielectric transmission
if (!(sampleFlags & BxDFReflTransFlags::Transmission))
return 0;
if (Dot(wi, wh) * Dot(wo, wh) > 0)
if (Dot(wi, wm) * Dot(wo, wm) > 0)
return 0;
// Compute Fresnel factor and probabilities for dielectric transmission PDF
Float F = FrDielectric(Dot(wo, wh), eta);
Float F = FrDielectric(Dot(wo, wm), eta);
Float pr = F, pt = 1 - F;
if (pt == 0)
return 0;
......@@ -267,8 +263,8 @@ Float DielectricBxDF::PDF(Vector3f wo, Vector3f wi, TransportMode mode,
if (!(sampleFlags & BxDFReflTransFlags::Reflection))
pr = 0;
Float dwh_dwi = AbsDot(wi, wh) / Sqr(Dot(wo, wh) + etap * Dot(wi, wh));
return mfDistrib.PDF(wo, wh) * dwh_dwi * pt / (pr + pt);
Float dwm_dwi = AbsDot(wi, wm) / Sqr(Dot(wo, wm) + etap * Dot(wi, wm));
return mfDistrib.PDF(wo, wm) * dwm_dwi * pt / (pr + pt);
}
}
......
......@@ -217,7 +217,6 @@ class DielectricBxDF {
PBRT_CPU_GPU
SampledSpectrum f(Vector3f wo, Vector3f wi, TransportMode mode) const;
PBRT_CPU_GPU
Float PDF(Vector3f wo, Vector3f wi, TransportMode mode,
BxDFReflTransFlags sampleFlags = BxDFReflTransFlags::All) const;
......@@ -336,26 +335,25 @@ class ConductorBxDF {
return BSDFSample(f, wi, 1, BxDFFlags::SpecularReflection);
}
// Sample rough conductor BRDF
// Sample microfacet orientation $\wh$ and reflected direction $\wi$
// Sample microfacet normal $\wm$ and reflected direction $\wi$
if (wo.z == 0)
return {};
Vector3f wh = mfDistrib.Sample_wm(wo, u);
Vector3f wi = Reflect(wo, wh);
Vector3f wm = mfDistrib.Sample_wm(wo, u);
Vector3f wi = Reflect(wo, wm);
if (!SameHemisphere(wo, wi))
return {};
// Compute PDF of _wi_ for microfacet reflection
Float pdf = mfDistrib.PDF(wo, wh) / (4 * AbsDot(wo, wh));
Float pdf = mfDistrib.PDF(wo, wm) / (4 * AbsDot(wo, wm));
Float cosTheta_o = AbsCosTheta(wo), cosTheta_i = AbsCosTheta(wi);
if (cosTheta_i == 0 || cosTheta_o == 0)
return {};
// Evaluate Fresnel factor _F_ for conductor BRDF
Float frCosTheta_i = AbsDot(wi, wh);
SampledSpectrum F = FrComplex(frCosTheta_i, eta, k);
SampledSpectrum F = FrComplex(AbsDot(wo, wm), eta, k);
SampledSpectrum f =
mfDistrib.D(wh) * mfDistrib.G(wo, wi) * F / (4 * cosTheta_i * cosTheta_o);
mfDistrib.D(wm) * mfDistrib.G(wo, wi) * F / (4 * cosTheta_i * cosTheta_o);
return BSDFSample(f, wi, pdf, BxDFFlags::GlossyReflection);
}
......@@ -366,20 +364,19 @@ class ConductorBxDF {
if (mfDistrib.EffectivelySmooth())
return {};
// Evaluate rough conductor BRDF
// Compute cosines and $\wh$ for conductor BRDF
// Compute cosines and $\wm$ for conductor BRDF
Float cosTheta_o = AbsCosTheta(wo), cosTheta_i = AbsCosTheta(wi);
if (cosTheta_i == 0 || cosTheta_o == 0)
return {};
Vector3f wh = wi + wo;
if (wh.x == 0 && wh.y == 0 && wh.z == 0)
Vector3f wm = wi + wo;
if (wm.x == 0 && wm.y == 0 && wm.z == 0)
return {};
wh = Normalize(wh);
wm = Normalize(wm);
// Evaluate Fresnel factor _F_ for conductor BRDF
Float frCosTheta_i = AbsDot(wi, wh);
SampledSpectrum F = FrComplex(frCosTheta_i, eta, k);
SampledSpectrum F = FrComplex(AbsDot(wo, wm), eta, k);
return mfDistrib.D(wh) * mfDistrib.G(wo, wi) * F / (4 * cosTheta_i * cosTheta_o);
return mfDistrib.D(wm) * mfDistrib.G(wo, wi) * F / (4 * cosTheta_i * cosTheta_o);
}
PBRT_CPU_GPU
......
......@@ -135,10 +135,10 @@ FilterSampler::FilterSampler(Filter filter, Allocator alloc)
f(int(32 * filter.Radius().x), int(32 * filter.Radius().y), alloc),
distrib(alloc) {
// Tabularize unnormalized filter function in _f_
for (int y = 0; y < f.ySize(); ++y)
for (int x = 0; x < f.xSize(); ++x) {
for (int y = 0; y < f.YSize(); ++y)
for (int x = 0; x < f.XSize(); ++x) {
Point2f p =
domain.Lerp(Point2f((x + 0.5f) / f.xSize(), (y + 0.5f) / f.ySize()));
domain.Lerp(Point2f((x + 0.5f) / f.XSize(), (y + 0.5f) / f.YSize()));
f(x, y) = filter.Evaluate(p);
}
......
......@@ -1164,7 +1164,7 @@ PortalImageInfiniteLight::PortalImageInfiniteLight(
});
// Initialize sampling distribution for portal image infinite light
auto duv_dw = [&](const Point2f &p) {
auto duv_dw = [&](Point2f p) {
Float duv_dw;
(void)RenderFromImage(p, &duv_dw);
return duv_dw;
......
......@@ -687,7 +687,7 @@ class PortalImageInfiniteLight : public LightBase {
}
PBRT_CPU_GPU
pstd::optional<Bounds2f> ImageBounds(const Point3f &p) const {
pstd::optional<Bounds2f> ImageBounds(Point3f p) const {
pstd::optional<Point2f> p0 = ImageFromRender(Normalize(portal[0] - p));
pstd::optional<Point2f> p1 = ImageFromRender(Normalize(portal[2] - p));
if (!p0 || !p1)
......
......@@ -83,7 +83,7 @@ class PowerLightSampler {
Float PDF(Light light) const {
if (!aliasTable.size())
return 0;
return aliasTable.PDF(lightToIndex[light]);
return aliasTable.PMF(lightToIndex[light]);
}
PBRT_CPU_GPU
......
......@@ -48,8 +48,9 @@ RGBSigmoidPolynomial RGBToSpectrumTable::operator()(const RGB &rgb) const {
i = j;
// Compute floating-point offsets into polynomial coefficient table
float z = rgb[i], x = rgb[(i + 1) % 3] * (res - 1) / z,
y = rgb[(i + 2) % 3] * (res - 1) / z;
float z = rgb[i];
float x = rgb[(i + 1) % 3] * (res - 1) / z;
float y = rgb[(i + 2) % 3] * (res - 1) / z;
// Compute integer indices and offsets for coefficient interpolation
constexpr int nCoeffs = 3;
......
......@@ -213,7 +213,7 @@ class Array2D {
Array2D(int nx, int ny, T def, allocator_type allocator = {})
: Array2D({{0, 0}, {nx, ny}}, def, allocator) {}
Array2D(const Array2D &a, allocator_type allocator = {})
: Array2D(a.begin(), a.end(), a.xSize(), a.ySize(), allocator) {}
: Array2D(a.begin(), a.end(), a.XSize(), a.YSize(), allocator) {}
~Array2D() {
int n = extent.Area();
......@@ -282,9 +282,9 @@ class Array2D {
PBRT_CPU_GPU
int size() const { return extent.Area(); }
PBRT_CPU_GPU
int xSize() const { return extent.pMax.x - extent.pMin.x; }
int XSize() const { return extent.pMax.x - extent.pMin.x; }
PBRT_CPU_GPU
int ySize() const { return extent.pMax.y - extent.pMin.y; }
int YSize() const { return extent.pMax.y - extent.pMin.y; }
PBRT_CPU_GPU
iterator begin() { return values; }
......@@ -804,8 +804,8 @@ class SampledGrid {
}
size_t BytesAllocated() const { return values.size() * sizeof(T); }
int xSize() const { return nx; }
int ySize() const { return ny; }
int XSize() const { return nx; }
int YSize() const { return ny; }
int zSize() const { return nz; }
const_iterator begin() const { return values.begin(); }
......
......@@ -17,8 +17,8 @@ TEST(Array2D, Basics) {
const int nx = 5, ny = 9;
Array2D<Float> a(nx, ny);
EXPECT_EQ(nx, a.xSize());
EXPECT_EQ(ny, a.ySize());
EXPECT_EQ(nx, a.XSize());
EXPECT_EQ(ny, a.YSize());
EXPECT_EQ(nx * ny, a.size());
for (int y = 0; y < ny; ++y)
......@@ -34,8 +34,8 @@ TEST(Array2D, Bounds) {
Bounds2i b(Point2i(-4, 3), Point2i(10, 7));
Array2D<Point2f> a(b);
EXPECT_EQ(b.pMax.x - b.pMin.x, a.xSize());
EXPECT_EQ(b.pMax.y - b.pMin.y, a.ySize());
EXPECT_EQ(b.pMax.x - b.pMin.x, a.XSize());
EXPECT_EQ(b.pMax.y - b.pMin.y, a.YSize());
for (Point2i p : b)
a[p] = Point2f(p.y, p.x);
......
......@@ -369,7 +369,7 @@ Point2f RejectionSampleDisk(RNG &rng) {
do {
p.x = 1 - 2 * rng.Uniform<Float>();
p.y = 1 - 2 * rng.Uniform<Float>();
} while (p.x * p.x + p.y * p.y > 1);
} while (Sqr(p.x) + Sqr(p.y) > 1);
return p;
}
......@@ -557,7 +557,7 @@ AliasTable::AliasTable(pstd::span<const Float> weights, Allocator alloc)
Float sum = std::accumulate(weights.begin(), weights.end(), 0.);
CHECK_GT(sum, 0);
for (size_t i = 0; i < weights.size(); ++i)
bins[i].pdf = weights[i] / sum;
bins[i].p = weights[i] / sum;
// Create alias table work lists
struct Outcome {
......@@ -567,7 +567,7 @@ AliasTable::AliasTable(pstd::span<const Float> weights, Allocator alloc)
std::vector<Outcome> under, over;
for (size_t i = 0; i < bins.size(); ++i) {
// Add outcome _i_ to an alias table work list
Float pHat = bins[i].pdf * bins.size();
Float pHat = bins[i].p * bins.size();
if (pHat < 1)
under.push_back(Outcome{pHat, i});
else
......@@ -608,16 +608,16 @@ AliasTable::AliasTable(pstd::span<const Float> weights, Allocator alloc)
}
}
int AliasTable::Sample(Float u, Float *pdfOut, Float *uRemapped) const {
int AliasTable::Sample(Float u, Float *pmf, Float *uRemapped) const {
// Compute alias table _offset_ and remapped random sample _up_
int offset = std::min<int>(u * bins.size(), bins.size() - 1);
Float up = std::min<Float>(u * bins.size() - offset, OneMinusEpsilon);
if (up < bins[offset].q) {
// Return sample for alias table at _offset_
DCHECK_GT(bins[offset].pdf, 0);
if (pdfOut)
*pdfOut = bins[offset].pdf;
DCHECK_GT(bins[offset].p, 0);
if (pmf)
*pmf = bins[offset].p;
if (uRemapped)
*uRemapped = std::min<Float>(up / bins[offset].q, OneMinusEpsilon);
return offset;
......@@ -626,9 +626,9 @@ int AliasTable::Sample(Float u, Float *pdfOut, Float *uRemapped) const {
// Return sample for alias table at _alias[offset]_
int alias = bins[offset].alias;
DCHECK_GE(alias, 0);
DCHECK_GT(bins[alias].pdf, 0);
if (pdfOut)
*pdfOut = bins[alias].pdf;
DCHECK_GT(bins[alias].p, 0);
if (pmf)
*pmf = bins[alias].p;
if (uRemapped)
*uRemapped = std::min<Float>((up - bins[offset].q) / (1 - bins[offset].q),
OneMinusEpsilon);
......@@ -639,7 +639,7 @@ int AliasTable::Sample(Float u, Float *pdfOut, Float *uRemapped) const {
std::string AliasTable::ToString() const {
std::string s = "[ AliasTable bins: [ ";
for (const auto &b : bins)
s += StringPrintf("[ Bin q: %f pdf: %f alias: %d ] ", b.q, b.pdf, b.alias);
s += StringPrintf("[ Bin q: %f p: %f alias: %d ] ", b.q, b.p, b.alias);
return s + "] ]";
}
......
......@@ -303,21 +303,21 @@ PBRT_CPU_GPU inline Float InvertSmoothStepSample(Float x, Float a, Float b) {
return (P(x) - P(a)) / (P(b) - P(a));
}
PBRT_CPU_GPU inline Point2f SampleUniformDiskPolar(const Point2f &u) {
PBRT_CPU_GPU inline Point2f SampleUniformDiskPolar(Point2f u) {
Float r = std::sqrt(u[0]);
Float theta = 2 * Pi * u[1];
return {r * std::cos(theta), r * std::sin(theta)};
}
PBRT_CPU_GPU
inline Point2f InvertUniformDiskPolarSample(const Point2f &p) {
inline Point2f InvertUniformDiskPolarSample(Point2f p) {
Float phi = std::atan2(p.y, p.x);
if (phi < 0)
phi += 2 * Pi;
return Point2f(Sqr(p.x) + Sqr(p.y), phi / (2 * Pi));
}
PBRT_CPU_GPU inline Point2f SampleUniformDiskConcentric(const Point2f &u) {
PBRT_CPU_GPU inline Point2f SampleUniformDiskConcentric(Point2f u) {
// Map _u_ to $[-1,1]^2$ and handle degeneracy at the origin
Point2f uOffset = 2 * u - Vector2f(1, 1);
if (uOffset.x == 0 && uOffset.y == 0)
......@@ -336,7 +336,7 @@ PBRT_CPU_GPU inline Point2f SampleUniformDiskConcentric(const Point2f &u) {
}
PBRT_CPU_GPU
inline Point2f InvertUniformDiskConcentricSample(const Point2f &p) {
inline Point2f InvertUniformDiskConcentricSample(Point2f p) {
Float theta = std::atan2(p.y, p.x); // -pi -> pi
Float r = std::sqrt(Sqr(p.x) + Sqr(p.y));
......@@ -365,7 +365,7 @@ inline Point2f InvertUniformDiskConcentricSample(const Point2f &p) {
return {(uo.x + 1) / 2, (uo.y + 1) / 2};
}
PBRT_CPU_GPU inline Vector3f SampleUniformHemisphere(const Point2f &u) {
PBRT_CPU_GPU inline Vector3f SampleUniformHemisphere(Point2f u) {
Float z = u[0];
Float r = SafeSqrt(1 - Sqr(z));
Float phi = 2 * Pi * u[1];
......@@ -376,16 +376,16 @@ PBRT_CPU_GPU inline Float UniformHemispherePDF() {
return Inv2Pi;
}
PBRT_CPU_GPU inline Point2f InvertUniformHemisphereSample(const Vector3f &v) {
Float phi = std::atan2(v.y, v.x);
PBRT_CPU_GPU inline Point2f InvertUniformHemisphereSample(Vector3f w) {
Float phi = std::atan2(w.y, w.x);
if (phi < 0)
phi += 2 * Pi;
return Point2f(v.z, phi / (2 * Pi));
return Point2f(w.z, phi / (2 * Pi));
}
PBRT_CPU_GPU inline Vector3f SampleUniformSphere(const Point2f &u) {
PBRT_CPU_GPU inline Vector3f SampleUniformSphere(Point2f u) {
Float z = 1 - 2 * u[0];
Float r = SafeSqrt(1 - z * z);
Float r = SafeSqrt(1 - Sqr(z));
Float phi = 2 * Pi * u[1];
return {r * std::cos(phi), r * std::sin(phi), z};
}
......@@ -394,14 +394,14 @@ PBRT_CPU_GPU inline Float UniformSpherePDF() {
return Inv4Pi;
}
PBRT_CPU_GPU inline Point2f InvertUniformSphereSample(const Vector3f &v) {
Float phi = std::atan2(v.y, v.x);
PBRT_CPU_GPU inline Point2f InvertUniformSphereSample(Vector3f w) {
Float phi = std::atan2(w.y, w.x);
if (phi < 0)
phi += 2 * Pi;
return Point2f((1 - v.z) / 2, phi / (2 * Pi));
return Point2f((1 - w.z) / 2, phi / (2 * Pi));
}
PBRT_CPU_GPU inline Vector3f SampleCosineHemisphere(const Point2f &u) {
PBRT_CPU_GPU inline Vector3f SampleCosineHemisphere(Point2f u) {
Point2f d = SampleUniformDiskConcentric(u);
Float z = SafeSqrt(1 - Sqr(d.x) - Sqr(d.y));
return Vector3f(d.x, d.y, z);
......@@ -411,25 +411,24 @@ PBRT_CPU_GPU inline Float CosineHemispherePDF(Float cosTheta) {
return cosTheta * InvPi;
}
PBRT_CPU_GPU inline Point2f InvertCosineHemisphereSample(const Vector3f &v) {
return InvertUniformDiskConcentricSample({v.x, v.y});
PBRT_CPU_GPU inline Point2f InvertCosineHemisphereSample(Vector3f w) {
return InvertUniformDiskConcentricSample({w.x, w.y});
}
PBRT_CPU_GPU inline Float UniformConePDF(Float cosThetaMax) {
return 1 / (2 * Pi * (1 - cosThetaMax));
}
PBRT_CPU_GPU inline Vector3f SampleUniformCone(const Point2f &u, Float cosThetaMax) {
PBRT_CPU_GPU inline Vector3f SampleUniformCone(Point2f u, Float cosThetaMax) {
Float cosTheta = (1 - u[0]) + u[0] * cosThetaMax;
Float sinTheta = SafeSqrt(1 - Sqr(cosTheta));
Float phi = u[1] * 2 * Pi;
return SphericalDirection(sinTheta, cosTheta, phi);
}
PBRT_CPU_GPU inline Point2f InvertUniformConeSample(const Vector3f &v,
Float cosThetaMax) {
Float cosTheta = v.z;
Float phi = SphericalPhi(v);
PBRT_CPU_GPU inline Point2f InvertUniformConeSample(Vector3f w, Float cosThetaMax) {
Float cosTheta = w.z;
Float phi = SphericalPhi(w);
return {(cosTheta - 1) / (cosThetaMax - 1), phi / (2 * Pi)};
}
......@@ -574,9 +573,8 @@ class WeightedReservoirSampler {
PBRT_CPU_GPU
void Merge(const WeightedReservoirSampler &wrs) {
DCHECK_LE(weightSum + wrs.WeightSum(), 1e80);
if (wrs.HasSample()) {
if (wrs.HasSample())
Add(wrs.reservoir, wrs.weightSum);
}
}
std::string ToString() const {
......@@ -645,7 +643,7 @@ class PiecewiseConstant1D {
PBRT_CPU_GPU
Float Integral() const { return funcInt; }
PBRT_CPU_GPU
size_t size() const { return func.size(); }
size_t Size() const { return func.size(); }
PBRT_CPU_GPU
Float Sample(Float u, Float *pdf = nullptr, int *offset = nullptr) const {
......@@ -665,7 +663,7 @@ class PiecewiseConstant1D {
*pdf = (funcInt > 0) ? func[o] / funcInt : 0;
// Return $x$ corresponding to sample
return Lerp((o + du) / size(), min, max);
return Lerp((o + du) / Size(), min, max);
}
PBRT_CPU_GPU
......@@ -699,10 +697,10 @@ class PiecewiseConstant2D {
: PiecewiseConstant2D(data, nx, ny, Bounds2f(Point2f(0, 0), Point2f(1, 1)),
alloc) {}
explicit PiecewiseConstant2D(const Array2D<Float> &data, Allocator alloc = {})
: PiecewiseConstant2D(pstd::span<const Float>(data), data.xSize(), data.ySize(),
: PiecewiseConstant2D(pstd::span<const Float>(data), data.XSize(), data.YSize(),
alloc) {}
PiecewiseConstant2D(const Array2D<Float> &data, Bounds2f domain, Allocator alloc = {})
: PiecewiseConstant2D(pstd::span<const Float>(data), data.xSize(), data.ySize(),
: PiecewiseConstant2D(pstd::span<const Float>(data), data.XSize(), data.YSize(),
domain, alloc) {}
PBRT_CPU_GPU
......@@ -717,7 +715,7 @@ class PiecewiseConstant2D {
PBRT_CPU_GPU
Point2i Resolution() const {
return {int(pConditionalV[0].size()), int(pMarginal.size())};
return {int(pConditionalV[0].Size()), int(pMarginal.Size())};
}
std::string ToString() const {
......@@ -751,8 +749,7 @@ class PiecewiseConstant2D {
Float Integral() const { return pMarginal.Integral(); }
PBRT_CPU_GPU
Point2f Sample(const Point2f &u, Float *pdf = nullptr,
Point2i *offset = nullptr) const {
Point2f Sample(Point2f u, Float *pdf = nullptr, Point2i *offset = nullptr) const {
Float pdfs[2];
Point2i uv;
Float d1 = pMarginal.Sample(u[1], &pdfs[1], &uv[1]);
......@@ -765,11 +762,11 @@ class PiecewiseConstant2D {
}
PBRT_CPU_GPU
Float PDF(const Point2f &pr) const {
Float PDF(Point2f pr) const {
Point2f p = Point2f(domain.Offset(pr));
int iu =
Clamp(int(p[0] * pConditionalV[0].size()), 0, pConditionalV[0].size() - 1);
int iv = Clamp(int(p[1] * pMarginal.size()), 0, pMarginal.size() - 1);
Clamp(int(p[0] * pConditionalV[0].Size()), 0, pConditionalV[0].Size() - 1);
int iv = Clamp(int(p[1] * pMarginal.Size()), 0, pMarginal.Size() - 1);
return pConditionalV[iv].func[iu] / pMarginal.Integral();
}
......@@ -804,18 +801,18 @@ class AliasTable {
AliasTable(pstd::span<const Float> weights, Allocator alloc = {});
PBRT_CPU_GPU
int Sample(Float u, Float *pdf = nullptr, Float *uRemapped = nullptr) const;
int Sample(Float u, Float *pmf = nullptr, Float *uRemapped = nullptr) const;
std::string ToString() const;
PBRT_CPU_GPU
size_t size() const { return bins.size(); }
PBRT_CPU_GPU
Float PDF(int index) const { return bins[index].pdf; }
Float PMF(int index) const { return bins[index].p; }
private:
// AliasTable Private Members
struct Bin {
Float q, pdf;
Float q, p;
int alias;
};
pstd::vector<Bin> bins;
......@@ -827,17 +824,17 @@ class SummedAreaTable {
// SummedAreaTable Public Methods
SummedAreaTable(Allocator alloc) : sum(alloc) {}
SummedAreaTable(const Array2D<Float> &values, Allocator alloc = {})
: sum(values.xSize(), values.ySize(), alloc) {
: sum(values.XSize(), values.YSize(), alloc) {
sum(0, 0) = values(0, 0);
// Compute sums along first scanline and column
for (int x = 1; x < sum.xSize(); ++x)
// Compute sums along first row and column
for (int x = 1; x < sum.XSize(); ++x)
sum(x, 0) = values(x, 0) + sum(x - 1, 0);
for (int y = 1; y < sum.ySize(); ++y)
for (int y = 1; y < sum.YSize(); ++y)
sum(0, y) = values(0, y) + sum(0, y - 1);
// Compute sums for the remainder of the entries
for (int y = 1; y < sum.ySize(); ++y)
for (int x = 1; x < sum.xSize(); ++x)
for (int y = 1; y < sum.YSize(); ++y)
for (int x = 1; x < sum.XSize(); ++x)
sum(x, y) =
(values(x, y) + sum(x - 1, y) + sum(x, y - 1) - sum(x - 1, y - 1));
}
......@@ -848,7 +845,7 @@ class SummedAreaTable {
(double)Lookup(extent.pMin.x, extent.pMax.y)) +
((double)Lookup(extent.pMin.x, extent.pMin.y) -
(double)Lookup(extent.pMax.x, extent.pMin.y)));
return std::max<Float>(s / (sum.xSize() * sum.ySize()), 0);
return std::max<Float>(s / (sum.XSize() * sum.YSize()), 0);
}
std::string ToString() const;
......@@ -858,8 +855,8 @@ class SummedAreaTable {
PBRT_CPU_GPU
Float Lookup(Float x, Float y) const {
// Rescale $(x,y)$ to table resolution and compute integer coordinates
x *= sum.xSize();
y *= sum.ySize();
x *= sum.XSize();
y *= sum.YSize();
int x0 = (int)x, y0 = (int)y;
// Bilinearly interpolate between surrounding table values
......@@ -877,8 +874,8 @@ class SummedAreaTable {
return 0;
// Reindex $(x,y)$ and return actual stored value
x = std::min(x - 1, sum.xSize() - 1);
y = std::min(y - 1, sum.ySize() - 1);
x = std::min(x - 1, sum.XSize() - 1);
y = std::min(y - 1, sum.YSize() - 1);
return sum(x, y);
}
......@@ -912,11 +909,11 @@ class WindowedPiecewiseConstant2D {
// Sample marginal windowed function in $x$
Point2f p;
p.x = SampleBisection(Px, u[0], b.pMin.x, b.pMax.x, func.xSize());
p.x = SampleBisection(Px, u[0], b.pMin.x, b.pMax.x, func.XSize());
// Sample conditional windowed function in $y$
// Compute 2D bounds _bCond_ for conditional sampling
int nx = func.xSize();
int nx = func.XSize();
Bounds2f bCond(Point2f(pstd::floor(p.x * nx) / nx, b.pMin.y),
Point2f(pstd::ceil(p.x * nx) / nx, b.pMax.y));
if (bCond.pMin.x == bCond.pMax.x)
......@@ -933,7 +930,7 @@ class WindowedPiecewiseConstant2D {
by.pMax.y = y;
return sat.Integral(by) / condIntegral;
};
p.y = SampleBisection(Py, u[1], b.pMin.y, b.pMax.y, func.ySize());
p.y = SampleBisection(Py, u[1], b.pMin.y, b.pMax.y, func.YSize());
// Compute PDF and return point sampled from windowed function
*pdf = Eval(p) / bInt;
......@@ -971,8 +968,8 @@ class WindowedPiecewiseConstant2D {
PBRT_CPU_GPU
Float Eval(Point2f p) const {
Point2i pi(std::min<int>(p[0] * func.xSize(), func.xSize() - 1),
std::min<int>(p[1] * func.ySize(), func.ySize() - 1));
Point2i pi(std::min<int>(p[0] * func.XSize(), func.XSize() - 1),
std::min<int>(p[1] * func.YSize(), func.YSize() - 1));
return func[pi];
}
......
......@@ -215,12 +215,12 @@ TEST(Sobol, IntervalToIndexRandoms) {
TEST(PiecewiseConstant1D, Continuous) {
PiecewiseConstant1D dist({1.f, 1.f, 2.f, 4.f, 8.f});
EXPECT_EQ(5, dist.size());
EXPECT_EQ(5, dist.Size());
Float pdf;
int offset;
EXPECT_EQ(0., dist.Sample(0., &pdf, &offset));
EXPECT_FLOAT_EQ(dist.size() * 1. / 16., pdf);
EXPECT_FLOAT_EQ(dist.Size() * 1. / 16., pdf);
EXPECT_EQ(0, offset);
// Right at the boundary between the 4 and the 8 segments.
......@@ -228,7 +228,7 @@ TEST(PiecewiseConstant1D, Continuous) {
// Middle of the 8 segment
EXPECT_FLOAT_EQ(.9, dist.Sample(0.75, &pdf, &offset));
EXPECT_FLOAT_EQ(dist.size() * 8. / 16., pdf);
EXPECT_FLOAT_EQ(dist.Size() * 8. / 16., pdf);
EXPECT_EQ(4, offset);
EXPECT_FLOAT_EQ(0., dist.Sample(0., &pdf));
......@@ -1212,8 +1212,8 @@ TEST(AliasTable, Varying) {
TEST(SummedArea, Constant) {
Array2D<Float> v(4, 4);
for (int y = 0; y < v.ySize(); ++y)
for (int x = 0; x < v.xSize(); ++x)
for (int y = 0; y < v.YSize(); ++y)
for (int x = 0; x < v.XSize(); ++x)
v(x, y) = 1;
SummedAreaTable sat(v);
......@@ -1228,25 +1228,25 @@ TEST(SummedArea, Constant) {
TEST(SummedArea, Rect) {
Array2D<Float> v(8, 4);
for (int y = 0; y < v.ySize(); ++y)
for (int x = 0; x < v.xSize(); ++x)
for (int y = 0; y < v.YSize(); ++y)
for (int x = 0; x < v.XSize(); ++x)
v(x, y) = x + y;
SummedAreaTable sat(v);
// All boxes that line up with boundaries exactly
for (int y0 = 0; y0 < v.ySize(); ++y0)
for (int x0 = 0; x0 < v.xSize(); ++x0)
for (int y1 = y0; y1 < v.ySize(); ++y1)
for (int x1 = x0; x1 < v.xSize(); ++x1) {
for (int y0 = 0; y0 < v.YSize(); ++y0)
for (int x0 = 0; x0 < v.XSize(); ++x0)
for (int y1 = y0; y1 < v.YSize(); ++y1)
for (int x1 = x0; x1 < v.XSize(); ++x1) {
Float mySum = 0;
for (int y = y0; y < y1; ++y)
for (int x = x0; x < x1; ++x)
mySum += v(x, y);
Bounds2f b(Point2f(Float(x0) / v.xSize(), Float(y0) / v.ySize()),
Point2f(Float(x1) / v.xSize(), Float(y1) / v.ySize()));
EXPECT_EQ(mySum / (v.xSize() * v.ySize()), sat.Integral(b));
Bounds2f b(Point2f(Float(x0) / v.XSize(), Float(y0) / v.YSize()),
Point2f(Float(x1) / v.XSize(), Float(y1) / v.YSize()));
EXPECT_EQ(mySum / (v.XSize() * v.YSize()), sat.Integral(b));
}
}
......@@ -1256,24 +1256,24 @@ TEST(SummedArea, Randoms) {
for (const auto d : dims) {
Array2D<Float> v(d[0], d[1]);
for (int y = 0; y < v.ySize(); ++y)
for (int x = 0; x < v.xSize(); ++x)
for (int y = 0; y < v.YSize(); ++y)
for (int x = 0; x < v.XSize(); ++x)
v(x, y) = rng.Uniform<int>(32);
SummedAreaTable sat(v);
for (int i = 0; i < 100; ++i) {
Bounds2i bi({rng.Uniform<int>(v.xSize()), rng.Uniform<int>(v.ySize())},
{rng.Uniform<int>(v.xSize()), rng.Uniform<int>(v.ySize())});
Bounds2f bf(Point2f(Float(bi.pMin.x) / Float(v.xSize()),
Float(bi.pMin.y) / Float(v.ySize())),
Point2f(Float(bi.pMax.x) / Float(v.xSize()),
Float(bi.pMax.y) / Float(v.ySize())));
Bounds2i bi({rng.Uniform<int>(v.XSize()), rng.Uniform<int>(v.YSize())},
{rng.Uniform<int>(v.XSize()), rng.Uniform<int>(v.YSize())});
Bounds2f bf(Point2f(Float(bi.pMin.x) / Float(v.XSize()),
Float(bi.pMin.y) / Float(v.YSize())),
Point2f(Float(bi.pMax.x) / Float(v.XSize()),
Float(bi.pMax.y) / Float(v.YSize())));
double ref = 0;
for (Point2i p : bi)
ref += v[p];
ref /= v.xSize() * v.ySize();
ref /= v.XSize() * v.YSize();
double s = sat.Integral(bf);
if (ref != s)
......@@ -1289,8 +1289,8 @@ TEST(SummedArea, NonCellAligned) {
for (const auto d : dims) {
Array2D<Float> v(d[0], d[1]);
for (int y = 0; y < v.ySize(); ++y)
for (int x = 0; x < v.xSize(); ++x)
for (int y = 0; y < v.YSize(); ++y)
for (int x = 0; x < v.XSize(); ++x)
v(x, y) = rng.Uniform<int>(32);
SummedAreaTable sat(v);
......@@ -1302,7 +1302,7 @@ TEST(SummedArea, NonCellAligned) {
int nSamples = 100000;
for (Point2f u : Hammersley2D(nSamples)) {
Point2f p = b.Lerp(u);
Point2i pi(p.x * v.xSize(), p.y * v.ySize());
Point2i pi(p.x * v.XSize(), p.y * v.YSize());
sampledSum += v[pi];
}
Float sampledResult = sampledSum * b.Area() / nSamples;
......
Markdown is supported
0% .
You are about to add 0 people to the discussion. Proceed with caution.
先完成此消息的编辑!
想要评论请 注册