 |
|
|
|
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Can someone check what's wrong with the following code? It's supposed to
implement the formula (4) from section 2,1 of this paper:
http://graphics.stanford.edu/papers/bssrdf/
The #if'ed-out section implements the approximation given by the first formula
of section 2.4 of the same paper, and seems to work perfectly.
I've been reviewing the code over and over, comparing it to the formula given in
the paper, another implementation, and various tutorials and papers further
exploring on the subject, but still can't make out why it doesn't seem to work.
Note that the paper doesn't specify a formula for "alpha_prime", which is
supposed to be the "reduced albedo". In most sources I found this defined as
"sigma_prime_s / sigma_prime_t". Some source, however, specifies it as
"sigma_prime_s / sigma_a", which makes the code yield *similar* results as the
approximation does with the "standard" definition of alpha_prime. However, I
can't imagine that the BSSRDF paper uses two concurring definitions, and anyway
the results are still a deal off the mark, and also seem to get closer or
further off depending on material parameters. (See [1] in my code)
Also note that formula (4) seems asymmetric regarding "zr" and "zv"; in fact,
various other sources add "zr" as an additional factor to the first term in the
angular brackets, making the formula symmetric. Strangely enough, this factor
does not seem to have any noticeable effect. (See [3] in my code)
One source also had a factor of 2 instead of 3 inside the root in "sigma_t_r".
Again, I did not see any noticeable impact. (See [2] in my code)
So, here's my C++ implementation - anyone able to spot the difference?
-----
double sigma_prime_t = sigma_prime_s + sigma_a;
double F_dr = (-1.440 / Sqr(mu) + 0.710 / mu + 0.668 + 0.0636 * mu);
double Aconst = ((1 + F_dr) / (1 - F_dr));
double dist;
double Rd;
double phi_in, phi_out;
VDist(dist, in.IPoint, out.IPoint);
dist *= sceneData->mmPerUnit;
#if 1
// full BSSRDF model
// TODO FIXME - why doesn't this give roughly similar results as the BRDF
approximation below??
double Dconst = 1 / (3 * sigma_prime_t);
double alpha_prime = sigma_prime_s / sigma_a; // [1]
double sigma_t_r = sqrt(3 * sigma_a * sigma_prime_t); // [2]
double zr = 1.0 / sigma_prime_t;
double dr = sqrt(Sqr(zr) + Sqr(dist));
double zv = zr + 4 * Aconst * Dconst;
double dv = sqrt (Sqr(zv) + Sqr(dist));
// calculate Rd
double m = alpha_prime / (4 * M_PI);
double t1 = zr * (sigma_t_r + 1.0/dr); // [3]
t1 *= exp(-sigma_t_r * dr);
t1 /= (sigma_prime_t * Sqr(dr));
double t2 = zv * (sigma_t_r + 1.0/dv);
t2 *= exp(-sigma_t_r * dv);
t2 /= (sigma_prime_t * Sqr(dv));
Rd = m * (t1 + t2);
#else
// uniform illumination BRDF approximation
double alpha_prime = sigma_prime_s / sigma_prime_t;
// calculate Rd
double root = sqrt(3.0 * (1.0 - alpha_prime));
Rd = (alpha_prime / 2.0) * (1 + exp((-4.0/3.0) * Aconst * root)) *
exp(-root);
#endif
VDot(phi_in, *vIn, in.INormal);
VDot(phi_out, *vOut, out.INormal);
phi_in = acos(phi_in);
phi_out = acos(phi_out);
*sd = 1 / M_PI;
*sd *= ComputeFt(phi_in, mu);
*sd *= Rd;
*sd *= ComputeFt(phi_out, mu);
-----
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
"clipka" <nomail@nomail> wrote:
> http://graphics.stanford.edu/papers/bssrdf/
I just worked through a good share of this the other day. You'd sure be in luck
if I had a better grasp on it...
> Note that the paper doesn't specify a formula for "alpha_prime", which is
> supposed to be the "reduced albedo". In most sources I found this defined as
> "sigma_prime_s / sigma_prime_t". Some source, however, specifies it as
> "sigma_prime_s / sigma_a", which makes the code yield *similar* results as the
> approximation does with the "standard" definition of alpha_prime.
You can find an implementation (PBRT) here. According to their code,
Spectrum ap = (sigma_sp / sigma_tp);
And they get good results...
> So, here's my C++ implementation - anyone able to spot the difference?
I'll try my best.
> double sigma_prime_t = sigma_prime_s + sigma_a;
> double F_dr = (-1.440 / Sqr(mu) + 0.710 / mu + 0.668 + 0.0636 * mu);
eta? Fair enough. Where did you get the sigma_prime_s and sigma_a from? Does
this require evaluation of the mean cosine, or do you assume isotropy? Should
these be vector parameters for the rgb values?
> double Aconst = ((1 + F_dr) / (1 - F_dr));
These are all promoted correctly to double precision, but I try to avoid getting
lazy here since I've done some pretty stupid things in the past...
> double dist;
> double Rd;
> double phi_in, phi_out;
>
> VDist(dist, in.IPoint, out.IPoint);
> dist *= sceneData->mmPerUnit;
>
> #if 1
> // full BSSRDF model
>
> // TODO FIXME - why doesn't this give roughly similar results as the BRDF
> approximation below??
>
> double Dconst = 1 / (3 * sigma_prime_t);
> double alpha_prime = sigma_prime_s / sigma_a; // [1]
See note above regarding choices.
> double sigma_t_r = sqrt(3 * sigma_a * sigma_prime_t); // [2]
>
> double zr = 1.0 / sigma_prime_t;
> double dr = sqrt(Sqr(zr) + Sqr(dist));
"dr = ||x-xr|| is the distance from x to the real source, and dv = ||x-xv|| is
the distance from x to the virtual source." I need to take a closer look at
this section. You're sure of this formula?
> double zv = zr + 4 * Aconst * Dconst;
> double dv = sqrt (Sqr(zv) + Sqr(dist));
Same as above.
> // calculate Rd
> double m = alpha_prime / (4 * M_PI);
> double t1 = zr * (sigma_t_r + 1.0/dr); // [3]
> t1 *= exp(-sigma_t_r * dr);
> t1 /= (sigma_prime_t * Sqr(dr));
For a second I thought there was an error in your simplification, but except for
the zr, this looks good. This does seem much more symmetric, so whether or not
it works, this form would seem more logical.
> double t2 = zv * (sigma_t_r + 1.0/dv);
> t2 *= exp(-sigma_t_r * dv);
> t2 /= (sigma_prime_t * Sqr(dv));
>
> Rd = m * (t1 + t2);
The code that you've provided looks good for the most part, but I'm assuming
that you have correct input and output. This doesn't look easy to debug since
there aren't many sanity checks for intermediate values.
- Ricky
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
"triple_r" <nomail@nomail> wrote:
> You can find an implementation (PBRT) here. According to their code,
Um... where's "here"? ;)
> > double sigma_prime_t = sigma_prime_s + sigma_a;
> > double F_dr = (-1.440 / Sqr(mu) + 0.710 / mu + 0.668 + 0.0636 * mu);
>
> eta? Fair enough. Where did you get the sigma_prime_s and sigma_a from? Does
> this require evaluation of the mean cosine, or do you assume isotropy? Should
> these be vector parameters for the rgb values?
Both sigma_prime_s and sigma_a are "input parameters" of the whole smash; R, G
and B values are processed separately, so they're just plain floating-point
values at this stage.
> > double sigma_t_r = sqrt(3 * sigma_a * sigma_prime_t); // [2]
> >
> > double zr = 1.0 / sigma_prime_t;
> > double dr = sqrt(Sqr(zr) + Sqr(dist));
>
> "dr = ||x-xr|| is the distance from x to the real source, and dv = ||x-xv|| is
> the distance from x to the virtual source." I need to take a closer look at
> this section. You're sure of this formula?
In that particular section of the paper, I assume that "x" should actually read
"xo" (xi doesn't make much sense, as the distances would then be zv and zr
respectively; nor would xv or xr be fitting, and there are no other points of
interest in the model). Enter Pythagoras.
> > // calculate Rd
> > double m = alpha_prime / (4 * M_PI);
> > double t1 = zr * (sigma_t_r + 1.0/dr); // [3]
> > t1 *= exp(-sigma_t_r * dr);
> > t1 /= (sigma_prime_t * Sqr(dr));
>
> For a second I thought there was an error in your simplification, but except for
> the zr, this looks good. This does seem much more symmetric, so whether or not
> it works, this form would seem more logical.
It also matches formulae (9) and (10) in this paper:
http://graphics.ucsd.edu/~henrik/papers/fast_bssrdf/
> The code that you've provided looks good for the most part, but I'm assuming
> that you have correct input and output. This doesn't look easy to debug since
> there aren't many sanity checks for intermediate values.
Right now I'm assuming that I have *something* wrong with the choice and
weighting of sample points.
Theoretically(!), if I'd monte-carlo-integrate over all surface points xi of the
object, "ad infinitum" if it is a plane, and with uniform density, the weighting
of each sample should be 1/numberOfSamples, right? (Just for argument's sake; I
know it's not possible on infinite objects)
So if I'd monte-carlo-integrate instead with a sample density proportional to
some function f(), each sample should be weighted by f()/SUM(f()), or am I
mistaken?
Assuming that I'm dealing with a "locally flat" geometry, I pick sample points
on the surface in a random direction at a distance of
r = fabs(log(rnd()))/sigma_tr
This should give me a sample density proportional to
(sigma_tr * exp(-sigma_tr * r) / sigma_tr) / (2*pi * r)
right? (The first term being due to the "radial density" and the other due to
the "angular density")
The further processing is obviously fine, as with the uniform illumination BRDF
approximation (plus single-scattering term) it gives results almost perfectly
matching the diffuse reflectance values presented in the paper.
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Looks like what I need is an additional factor of something like
1.0 / (sigma_a * sigma_prime_t)
.... but *why*? (And, more importantly, what's the *exact* factor I need? It
still seems slightly too dark; but the order of magnitude is a surprisingly
good match)
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
|
|
