I am currently writing a physically-based rendering engine.  I need to
deviate 3-Dimensional vectors in a random "cone" whose principle axis is
equal to the vector itself.  The size of the "cone" is pregiven in
steradians.

I understand that povray uses a lookup table to get directions for radiosity
samples.   I do not have this option.

My current method uses orthogonal basis vectors in the following way:

[1] First produce a random unit vector that deviates randomly from the
positive y axis <0,1,0>  by some given steradians, r.  This part is easy
and I have already solved it.

[2] Next I produce unit orthogonal basis vectors around a given surface
normal N, such that the "up" vector is equal to N, and the other two are
both orthogonal to N.

[3] I then project the random vector from [1] into the "axis system" in [2].

One would think this algorithm works perfectly.  IT DOES NOT. It produces
banding artifacts right along the parts of the surfaces where the largest
component of the surface Normal changes from one to the next.  Let me give
an example:  Imagine the surface normals on a cylinder morphing smoothly
from <0.6,0,0.7> to <0.7,0,0.6>   Somewhere along the surface, the x
component overtakes the z component of the vector. This causes an
instantaneous switch in the orthogonal basis vectors from part [2].
"Right" changes to "-Right" and so forth.   Random processes along this
boundary go into a little "minishift" that shows up in the rendered image
as an anomaly along the surface.

Does anyone know a different method for randomly deviating a vector from its
original direction into a new one such that the random vector remains within
a given cone of r steradians????

Ideally, your answer will not involve the use of sin() and cos(), for the
obvious reasons.

Thank you for your help.

Post a reply to this message

I am the original poster of this thread. This may be a hint towards the
problem:

http://mathworld.wolfram.com/SpherePointPicking.html

Please let me know if I am on the right track.

Post a reply to this message

Paris wrote:
> I am currently writing a physically-based rendering engine.  I need to
> deviate 3-Dimensional vectors in a random "cone" whose principle axis is
> equal to the vector itself.  The size of the "cone" is pregiven in
> steradians.
> 
> I understand that povray uses a lookup table to get directions for radiosity
> samples.   I do not have this option.

The radiosity samples in POV-Ray are actually distributed on a hemisphere.
they were once generated with a good algorithm randomly, but there is
nothing inherently forcing you to have a fixed table other than speed.

My guess is that ideally you want to use a hemisphere, too, but determined a
cone is simpler? Or is that a cone as in cone-tracing? Either way, although
I am not sure I understand your problem correctly, the best advice would be
to check out: <http://www.boost.org/libs/random/index.html> if you need a
good set of random number generators. It also comes with a mapping for
spheres that is really useful.

	Thorsten

Post a reply to this message

Paris scripsit::
> I am currently writing a physically-based rendering engine.  I need to
> deviate 3-Dimensional vectors in a random "cone" whose principle axis is
> equal to the vector itself.  The size of the "cone" is pregiven in
> steradians.
> 
> I understand that povray uses a lookup table to get directions for radiosity
> samples.   I do not have this option.
> 
> My current method uses orthogonal basis vectors in the following way:
> 
> [1] First produce a random unit vector that deviates randomly from the
> positive y axis <0,1,0>  by some given steradians, r.  This part is easy
> and I have already solved it.
> 

May I ask how is your distribution on the projection disc ? (at
least, at the specification level ?)
Is it uniform in density (deviation to <0,0.5> has same probability
than deviation to <0,1>, in cartesian coordinate) ?
Or is uniform in direction distribution with uniform deviation on
the radius (using polar coordinate, independantly) ?
Note that both appraoches are not identical in result.

Or yet another method ?

Have a look at

http://mathworld.wolfram.com/SpherePointPicking.html

and

http://mathworld.wolfram.com/DiskPointPicking.html


> [2] Next I produce unit orthogonal basis vectors around a given surface
> normal N, such that the "up" vector is equal to N, and the other two are
> both orthogonal to N.
> 

There come your problem!
You have a freedom which must be removed by your code... hence your
banding!
Remember that a sphere cannot have combed hairs!
Your method will always have at least one singularity.

> [3] I then project the random vector from [1] into the "axis system" in [2].

I would have rotated [y] to match [N] but that the same thing.

> Does anyone know a different method for randomly deviating a vector from its
> original direction into a new one such that the random vector remains within
> a given cone of r steradians????

Depend on your density of distribution along the cone axis...

Brute force is always an option: generate a random unit vector using
three random (cartesian coordinates, trap the bogus <0,0,0>, adjust
length after processing), check if it is in the cone, if not do it
again... Do not rotate/play with it, just use as is.

Oh, yes, use three seeds, not a single one for everythings !
(unless you check the random generator for non-correlation in 3D)

To check the cone thing, a simple scalar product of the two vector
should be enough, you just have to precompute the solid angle once
to the right value (as it involves computing portion of a sphere, it
might involve some heavy computation: solving
A = 2 pi h , the surface of a portion of unit sphere of width h,
then computing the radius of the projected disc...)


Other distribution might be different.

> 
> Ideally, your answer will not involve the use of sin() and cos(), for the
> obvious reasons.


-- 
Eifersucht ist die Leidenschaft, die mit Eifer sucht, was Leiden
schafft.

Eco: -8.75 Soc: -6.72
http://www.politicalcompass.org/

Post a reply to this message

Thorsten Froehlich <tho### [at] trfde> wrote:
> The radiosity samples in POV-Ray are actually distributed on a hemisphere.
> they were once generated with a good algorithm randomly, but there is
> nothing inherently forcing you to have a fixed table other than speed.
>
> My guess is that ideally you want to use a hemisphere, too, but determined a
> cone is simpler? Or is that a cone as in cone-tracing? Either way, although


The "cone" only appears if you consider a "bundle" of all the vectors going
towards a circular section of a hemisphere.  Imagine what you would do in
order to simulate "fuzzy reflection".  Given the reflected vector, you
would need to deviate it slightly and take several random samples.  This is
precisely what I need.

Many commercial rendering packages already do this.  The difference is that
when I render using 600+ samples I begin to get banding artifacts on the
surfaces.  This happens because the basis vectors "switch" between positive
and negative on places where the largest component switches from one axis to
the next.  I described this in my post already.


> I am not sure I understand your problem correctly, the best advice would be
> to check out: <http://www.boost.org/libs/random/index.html> if you need a
> good set of random number generators. It also comes with a mapping for
> spheres that is really useful.

I can already map to a sphere just fine.  I even posted a link to the
mathematics for uniformly-distributed spheres.  That's not my problem.  My
problem is: how do I then project this random vector (taken from +Y) into a
 generic surface normal N, such that the new random vector is only so many
steradians deviated from N, maximally?   My method seems to work on paper,
but produces banding artifacts for some esoteric reason having to do with
switching of Right into "negative Right" because of projection issues.

Ideally I would not have to ROTATE the vector into position. This is way too
slow. It takes at least 4 calls to trigonometric functions such as sin()
asin() and tan().  It seems like there should be some simple,
straightforward way to just project the Random-vector-deviated-from-Y+ into
the "axis space" of N.   My method seems to do this just fine on paper, but
it produces banding artifacts in the rendered image. I dont know how else I
can explain this in simpler language for you.

Post a reply to this message

Le Forgeron <jgr### [at] freefr> wrote:
> Paris scripsit::
> > [1] First produce a random unit vector that deviates randomly from the
> > positive y axis <0,1,0>  by some given steradians, r.  This part is easy
> > and I have already solved it.
> >
>
> May I ask how is your distribution on the projection disc ? (at
> least, at the specification level ?)
> Is it uniform in density (deviation to <0,0.5> has same probability
> than deviation to <0,1>, in cartesian coordinate) ?
> Or is uniform in direction distribution with uniform deviation on
> the radius (using polar coordinate, independantly) ?
> Note that both appraoches are not identical in result.
>
> Or yet another method ?

I am not using a disc actually.  That actually does not matter at this time,
because I am getting banding artifacts on surfaces that bounce a reflected
vector in a totally random direction independent of the angle of incidence.
 So for now, just consider that I am using a unit hemisphere.  I have used
several methods, including ones I cooked up on my own, as well as those
mentioned here:
http://mathworld.wolfram.com/SpherePointPicking.html   The banding artifacts
remain the same regardless of which hemisphere method I use.   I have
already concluded that random-distribution is NOT THE PROBLEM.  The problem
arrises when I project the randomly-deviated vector around <0,1,0>  into the
surface normal , N on the surface of the object.  N can be ANY VECTOR
including <0,1,0> but also including <0,-1,0> as well as any other vector
in space.  My "random cone" function produces a uniformly distributed set
of points that are, at maximum, S steradians deviated from Y+ (ie <0,1,0>).
  The bundle of such random unit vectors are in a "solid cone angle" that is
centered about Y+ and S steradians wide.   I then project this random vector
into an "axis system" that is built up around N.  (U, V, W orthogonal basis
vectors).  The final result is N, deviated randomly, but at most maximally
deviated by S steradians.  Repeated over and over, this gives rise to a
bundle of random vectors who are all contained in a "cone" that has N as
its principle axis.


> > [2] Next I produce unit orthogonal basis vectors around a given surface
> > normal N, such that the "up" vector is equal to N, and the other two are
> > both orthogonal to N.
> >
>
> There come your problem!
> You have a freedom which must be removed by your code... hence your
> banding!
> Remember that a sphere cannot have combed hairs!
> Your method will always have at least one singularity.

That's nice. But I was requesting that someone please help me find the
algorithm out there that does this.  My method works on paper, but not in
the real world application.  I mean, it "works" but it produces banding
artifacts.  I know nearly every commercial rendering package can produce
"fuzzy reflections" so this algorithm must be widely-known.  I had googled
for hours before giving up in exasperation and coming here as a last
resort.


> > [3] I then project the random vector from [1] into the "axis system" in [2].
>
> I would have rotated [y] to match [N] but that the same thing.

That requires something like 6 trigonometric function calls.  It is
ridiculously slow.  First you have to take 2 arc-sines in order to obtain
phi and theta. Then you have to perform 2 more trig functions each for each
rotation.  That is just not feasible, don't you agree?

Post a reply to this message

This is me responding to my own post again haha!
Anyways, I beleive I have solved the problem.  I will have rendered images
for you tommorrow that will demonstrate what the "banding" looks like.  It
seems I needed to use a homogenous axis system everywhere, except for two
singularities at the poles of the unit sphere, which are handled as rare
cases in the code.   Your suggestion about singularities really pointed me
in the right direction.

Post a reply to this message

Paris scripsit::
> My "random cone" function produces a uniformly distributed set
> of points that are, at maximum, S steradians deviated from Y+ (ie <0,1,0>).

Rewrite it to take N directly instead of Y!
Or publish it...

>   The bundle of such random unit vectors are in a "solid cone angle" that is
> centered about Y+ and S steradians wide.   I then project this random vector
> into an "axis system" that is built up around N.  (U, V, W orthogonal basis
> vectors).  The final result is N, deviated randomly, but at most maximally
> deviated by S steradians.  Repeated over and over, this gives rise to a
> bundle of random vectors who are all contained in a "cone" that has N as
> its principle axis.
> 
>>> [3] I then project the random vector from [1] into the "axis system" in [2].
>> I would have rotated [y] to match [N] but that the same thing.
> 
> That requires something like 6 trigonometric function calls.  It is
> ridiculously slow.  First you have to take 2 arc-sines in order to obtain
> phi and theta. Then you have to perform 2 more trig functions each for each
> rotation.  That is just not feasible, don't you agree?

Nop! Only one rotation please!
Stop thinging in your fixed referential.

You need the vector product of Y and N, let's call it A. Of course,
you have to trap Y aligned with N... (two trivial cases)
You then need to rotate V, your deviated Y, around A so that Y
become N. That's what vaxis_rotate is made for!


UNTESTED CODE :
/* it DOES produce banding */
/* You might want to allow for a larger zone than strict equality
=0 and =2... replaced with <allowance and >2-allowance */

#macro Change(V,Y,N)
#local delta=vlength(Y-N);
#if (delta=2)
-V
#else
#if (delta=0)
V
#else
#local A = vcross(Y,N); /* check the left/right hand, I did not */
#local Ang= degrees(acos(vdot(Y,N)));
vaxis_rotate(V,A,Ang)
#end
#end
#end

-- 
Eifersucht ist die Leidenschaft, die mit Eifer sucht, was Leiden
schafft.

Eco: -8.75 Soc: -6.72
http://www.politicalcompass.org/

Post a reply to this message

The top image contains the banding. The lower one does not.
http://home.comcast.net/~exspacea/banding_artifacts.png

The solution follows.

It uses precisely one, (I repeat for emphasis ONE) trig function call.  It
is tan() and it is contained in randomcone()  (Not included here).  I have
no idea how many trig functions Axis_Rotate() uses.  No doubt it is
probably something like 4 or 8 trig function calls.

C_Vector Vector_Machine::jittervector( C_Vector *N, double steradians )
{

 //Given some vector, N, deviate it in a random angular direction
 //within a solid cone that is 'steradians' wide.  Return the new
 //vector.   Despite the name of the variable I chose here,
 //steradians should be in DEGREES, not radians.  (..let the groaning
 //commence..)


 C_Vector bU, bV, bW, projIy;
 C_Vector Iy;

 // Get a random vector that deviates around the positive Y axis.
 Iy = randomcone(steradians);

 // Create an axis system where N is the "height" or "Up" direction
 getorthobasisvectors(N, &bU, &bV, &bW );

 // Project our random vector "as if" N were the Y axis.
 projIy = changeaxissystem( &Iy, &bU, &bV, &bW );

 return projIy;
}


void Vector_Machine::getorthobasisvectors( C_Vector *N,C_Vector *U,
C_Vector *V, C_Vector *W )
{
 /*
  Given some vector in space, N, derive an axis system of
  unit orthogonal basis vectors such that V=N.

   U is considered the "right" direction.
   V is considered the "up" direction.
   W is considered the "into" direction.
 */

 C_Vector rightV, intoV, uN;
 C_Vector eastpole, westpole;
 double singularity;

 uN = (*N);
// We will assume that the right/U direction is always X+ ie <1,0,0>
//  This assumption breaks down at two singularities, namely the east and
west
//  poles of the unit hemisphere.
//  For the east singularity, assume right/U is Z+ (<0,0,1>) and into/W is
Y+
// (<0,1,0>)
//  For the west singularity, assume right/U is Y+ and into/W is Z+

 eastpole.x=1.0;  eastpole.y=0; eastpole.z=0;
 westpole.x=-1.0;  westpole.y=0; westpole.z=0;

 singularity = uN * eastpole;
 if( singularity > 0.0 ) {
  singularity = 1.0 - singularity;
  if( singularity < EPSILON ) {
   rightV.x = 0; rightV.y = 0; rightV.z = 1; // Assumed right.
   goto sendreturn;
  }
 }

 singularity = uN * westpole;
 if( singularity > 0.0 ) {
  singularity = 1.0 - singularity;
  if( singularity < EPSILON ) {
   rightV.x = 0; rightV.y =1; rightV.z = 0; // Assumed right.
   goto sendreturn;
  }
 }

 //  ...otherwise, everyone else gets an X+.
 rightV.x = 1; rightV.y =0; rightV.z = 0;

sendreturn:
 intoV = rightV % uN;  // Give me an orthogonal vector. I dont carewhere.
 intoV.Unitize();
 rightV = uN % intoV;  // The 'true' right.
 rightV.Unitize();

 (*U) = rightV;
 (*V) = uN;
 (*W) = intoV;
}


C_Vector Vector_Machine::changeaxissystem(C_Vector *E, C_Vector *U, C_Vector
*V, C_Vector *W )
{
/*
Given some unit axis system U-V-W, project the point E into this
system.  Return this new point.
*/

 C_Vector P;

 P.x = (U->x)*(E->x)  + (V->x)*(E->y) + (W->x)*(E->z);
 P.y = (U->y)*(E->x)  + (V->y)*(E->y) + (W->y)*(E->z);
 P.z = (U->z)*(E->x)  + (V->z)*(E->y) + (W->z)*(E->z);

 return P;
}

Post a reply to this message

Paris wrote:
> The top image contains the banding. The lower one does not.
> http://home.comcast.net/~exspacea/banding_artifacts.png
> 
> The solution follows.
> 
> It uses precisely one, (I repeat for emphasis ONE) trig function call.  It
> is tan() and it is contained in randomcone()  (Not included here).  I have
> no idea how many trig functions Axis_Rotate() uses.  No doubt it is
> probably something like 4 or 8 trig function calls.
...

vaxis_rotate() uses 2 trig. function calls: 1 to sin() and 1 to cos()

Some relevant discussions can be found here:

http://news.povray.org/povray.binaries.programming/thread/%3C3BFFF01C.38CE12DB%40hotmail.com%3E
http://news.povray.org/povray.advanced-users/thread/%3C3C06D0B8.9DE747A0%40hotmail.com%3E

-- 
Tor Olav
http://subcube.com

Post a reply to this message