Coding a Lafortune specular shader, some questions though...

Hi everybody,
I started to code a Lafortune shader this week, since I’ve got some time.(I’m not sure I’m in the good section). And I’m encountering some various problems :slight_smile:

Before explaining them, I’ll try to give an overview of what is going on with this shading technique (otherwise, it’s hard to understand the problems I’m fighting against :smiley: )

Everything is explained in this paper : http://www.graphics.cornell.edu/pubs/1997/LFTG97.pdf
However, if you don’t want to read it, here’s the concept : the Phong specular shader computes the cosine between the view vector (v) and the vector of the incident light (l), mirrored around the surface normal vector (n). The cosine of this angle (a on the picture below) is equal to the dot product between r and v, provided those vectors are normalized, of course.

http://pix.nofrag.com/c/2/2/c6672ed949c32a1df21f34e17440a.jpg

Then, one elevates this cosine to a chosen power, depending on the wanted roughness of the material to be rendered : the higher this power (called n in the following lines) is, the harsher the highlight will look. This is the hardness parameter of Blender (hardness is an integer between 1 and 512 though, but it’s basically the same). One can also multiply this result by a number between 0 and 1, which could be the albedo of the material (in Blender, the “Spec” factor is between 0 and 200%… :D)
It gives this : f(r,v)=albedocos(r,v)^n = albedo(xrxv + yryv + zr*zv)^n

The mirroring around the normal n can be done using a householder matrix (mentioned in the paper). I do not understand what it is… We even didn’t learn anything about matrices in school, so, I’m not used to handle with them. I only know the very basics :slight_smile:

What’s new, is that, by replacing the householder matrix by a general matrix, some calcul lines later, you’ve got the dot product l.v multiplied by a diagonal matrix 33, whose coefficients can be chosen. It means you can give a different weight to the terms of the dot product, by varying a Cx, a Cy, and a Cz coefficient. The l and v vectors must be expressed in local coordinates, with z axis aligned to the surface normal n, and x/y axis aligned to the principal directions of anisotropy (I won’t be able to code anisotropy I think… It’s a lot of geometry, in space what’s more… :slight_smile: )
To sum up, here is the formula : f(r,v)=albedo
(Cx * xlxv + Cy * ylyv + Cz * zl*zv)^n
(as I said, I won’t code anisotropy, so Cx = Cy : the vectors are equally deformed along x and y axis).
This shading method allows to get a more realistic rendering of the highlights, especially at grazing angles, for which the ordinary specular shaders aren’t convincing at all :

http://pix.nofrag.com/9/7/1/20b52c5822af6e228a9bae3e3eb7c.jpg

(top left, Phong specular shader, top right, raytraced soft reflections, bottom left, Lafortune specular shader… The weakness of the Phong shader is obvious here. Even with a Ward isotropic shader, or a Cook Torrance shader, the shape of the highlights is badly described… :frowning:

Thas was the theory, here’s my code now. I commented everything to be commented, for my code to be clear. I just show the shading function : the others lines of code I added do not need to be shown (adding buttons in the interface, that kind of thing)

/* Lafortune isotropic spec */
static float Lafortune_Spec( float *n, float *l, float *v, int hard , int tangent, float Cxy, float Cz )
{
	float rslt ;
	float x[3] ;
	float y[3] ;
	
	//find x axis perpendicular to surface normal : x[0]*n[0]+x[1]*n[1]+x[2]*n[2]=0 //
	if (n[0] == 0.0)
	{
		x[0] = 0.0 ;
		x[1] = n[2] ;
		x[2] = -n[1] ;
	}
	else if (n[1] == 0.0)
	{
		x[0] = n[2] ;
		x[1] = 0.0 ;
		x[2] = -n[0] ;
	}
	else
	{
		x[0] = n[1] ;
		x[1] = -n[0] ;
		x[2] = 0.0 ;
	}
	Normalize (x);
	
	//find y axis perpendicular to x axis and surface normal, with cross product between n and x//
	y[0] = n[1] * x[2] - n[2] * x[1] ;
	y[1] = n[2] * x[0] - n[0] * x[2] ;
	y[2] = n[0] * x[1] - n[1] * x[0] ;
	Normalize (y) ;
		
	//get view vector in new local coordinates system//
	v[0] = v[0]*x[0] + v[1]*x[1] + v[2]*x[2] ; // = v.x = cos (v, x) since x and v are normalized
	v[1] = v[0]*y[0] + v[1]*y[1] + v[2]*y[2] ; // = v.y = cos (v, y)
	v[2] = v[0]*n[0] + v[1]*n[1] + v[2]*n[2] ; // = v.n = cos (v, n)
	Normalize (v) ;
	
	//get light vector in new local coordinates system//
	l[0] = x[0]*l[0] + x[1]*l[1] + x[2]*l[2] ; // = l.x = cos (l, x) since x and l are normalized.
	l[1] = y[0]*l[0] + y[1]*l[1] + y[2]*l[2] ; // = l.y = cos (l, y)
	l[2] = n[0]*l[0] + n[1]*l[1] + n[2]*l[2] ; // = l.n = cos (l, n)
	Normalize (l) ;
	
	rslt = Cxy*l[0]*v[0] + Cxy*l[1]*v[1] + Cz*l[2]*v[2] ;
	// weighted dot product between light vector and view vector //
	if (tangent) rslt = sasqrt( 1.0f - rslt*rslt ) ; // in case tangent shading is turned on
	if( rslt > 0.0f ) rslt= spec(rslt, hard); // take the hardness parameter into account
	else rslt = 0.0f; // clamp negative values to 0
	return rslt;
	
}

To sum up : I first find the local x axis which has to be perpendicular to the normal n (the code works quite simply : imagine n is : (3, 4, 5), it will give x (-4, 3, 0), so that n.x = 0). Then, thanks to the cross product between x and n, I find the y axis perpendicular to both x and n. I can then project l and v on x and y founded before, and on n. Doing so, I’ve got l and v in local coordinates, just like in the paper. Finally, I compute the weighted dot product between l and v, in which Cxy weights the x and y terms of the dot product, and Cz weights the z term of the dot product. My coordinates transform is quite strange, but by intuition, it should work. What do you guys think of it ? It could and would surely be useful to use matrices for this, in order for the code not to be burdened…

Now, here are my problems (finally :slight_smile: ). I would say I’ve got 3 problems, and they are related I think…

1- I wonder how vectors are… handled in blender. For instance : I setup a scene with a perfectly vertical Sun, toward the ground, and a perfectly horizontal plane, directed toward the sky, so, its normal is vertical, toward the sky. The Sun lamp throws rays from an infinite distance, so, they are all parallel ; the light vector is thus constant.
I decide to display xn, yn, zn, xl, yl, and zl in the console, and here’s what I get :
> n(-0.010817, 0.895343, 0.445245)
> l(0.010817, -0.895343, -0.445245)
It’s reassuring and strange at the same time… Reassuring because the vectors are actually opposite, but strange because I can’t figure out how vertical vectors can have such weird xyz components…

2- I wonder if my coordinates change is operational. Logically, yes, but the point is that at the moment, what I have IS NOT a lafortune shader… :smiley: It must be related to problem #1
Example : we are in a local coordinates system, with n = z axis. Which means, by choosing Cxy = -1 and Cz = 1, we should get a Phong shader : setting -1 for Cxy will mean mirroring the l vector around the normal, and this is what Phong does before computating the cosine between the 2 vectors. And, with those precise parameters, I get some nice speculars, but not a Phong shader, without a doubt.

And yet, some renders that I get give nice things, which look like what the Lafortune shader is capable of :

http://pix.nofrag.com/a/d/2/13378efd36a8c3481920adefceadb.jpg

So, I was wondering whether my coordinates change is valid or not.

3- The third problem is completely different from the 2 others : the specular highlights do not add correctly. With images :

http://pix.nofrag.com/2/7/b/abc22028312c05557da00a50ecaab.jpg

As we see, when shadows are activated, it’s better, but still hardly bugged :smiley: When an area is shadowed, blender seems to consider that speculars cannot be seen anyway, so it allows to show the highlights resulting from others lamps (picture on the right). But when the shadows are turned off, even an intensity of 0 for one specular highlight is enough to overwhelm the others. I’m not sure that it’s clear, that’s why I showed those images^^. (By the way, here we can see that the shading isn’t properly done… On the left image, it looks like tangent shading is turned on…)
I first thought that it didn’t come from the shading function itself, but rather from a wrong declaration somewhere, that kind of thing. But actually, it comes from the inside of the shading function I show here : I copied and pasted the code for the Phong shader into my own function, and everything worked perfectly… So, I have absolutely no idea of what is happening wrong.

That’s all folks, if you came at the end of this post, congrats :slight_smile: I hope I’ll find a response to my problems, 'cause the few promising results I got give me the will to go on… :smiley:

and for fun, a useless video with my shading function : http://www.vimeo.com/1441982
oh, btw, I posted exactly the same thing on a french forum, and I was asked why I’m not using PyNodes instead, which fit my needs exactly… The problem is that I don’t know anything about Python :slight_smile:

There is a good tutorial showing how to write basic diffuse shader (Lambert) using PyNodes at http://gallery.hiddenworlds.org/index.php?action=page&aid=17&page=2 . You might find it useufl. Another resource to check out is the cookbook, http://wiki.blender.org/index.php/Resources/PyNode_Cookbook .

Thanks BeBraw :slight_smile: I already took a look at those pages though… I even wrote the Lambert shader. I understood what I was coding, but I couldn’t invent the functions myself :smiley: Because I don’t know anything about Python, neither do I know how python scripts are integrated to blender.

Hey you convinced me, instead of having headaches on my C code, I’ll try it with a PyNode tomorrow :slight_smile: There is nothing really hard to code indeed (additions, multiplications, conditions with if/elseif/else… it cannot be so different in Python).

Pixelvore, for you as a programmer to get a quick overview of Python, I can recommend the online book A Byte of Python and for a handy reference on one page for the most important functions the Python Quick Reference.

Thanks, that was exactly what I was looking for :slight_smile:

For the first part of your code you could have used VecOrthoBasisf. Just do this to get x and y from n:

VecOrthoBasisf(n, x, y);

Regarding 1), all rendering calculations are done in camera space. That means, objects are transformed by the camera matrix, so the positions and vectors will depend on the camera position and rotation.

Regarding 2), I’m not sure your reasoning is correct, but I might be wrong. You say setting Cxy to -1 will mirror the view vector around the normal but I’m not sure that is true? For that, the x vector would have to point along the direction of that mirroring, but it’s actually chosen rather arbitrarily in the first part of your code. For that to work it should be computed more like this (I may have gotten the Cross product order wrong, and if v == n I’m not sure what happens):

tmp = Cross(v, n);
x = Normalize(Cross(tmp, n))

Separate Cx/Cy weights are only for anisotropic results, so for an isotropic shader the x/y you have should work.

VecOrthoBasisf(n, x, y);

Great, it will simplifie the code somehow. I took a look at the function, and basically, it does the same as mine, but at least, I’m sure that there won’t be bug coming from this part of the process. Thanks :slight_smile:

Regarding 1), all rendering calculations are done in camera space. That means, objects are transformed by the camera matrix, so the positions and vectors will depend on the camera position and rotation.

Okay. Is there a way to obtain the “raw” vectors for the surface normal, light, and so on ? [edit] even though it won’t change anything to the angle between L and V…

About what you say then, I think I understood what you mean. I took a 2nd look at the paper, and it is said at one moment : "the original cosine lobe model is obtained by choosing -Cx = -Cy = Cz = n-root(Cs), with n being the Phong exponent (hardness) (I’m not sure the word “n-root” exists in english.) Cs is a normalization factor = (n+2)/2*pi. This coefficient must be useful for physically plausibles renderers, where it is necessary to ensure energy conservation. But here, I just didn’t take into account this coefficient, so for me, Cxy = -1 and Cz = 1 should give the Phong model I think.

Another question… What does the “f” stand for in all calculs ? For example : if(rslt > 0.0f) … and so on. I’m very new to the source code of Blender… And to programming in general :slight_smile:

Thanks for your explanations anyway, they are very useful :slight_smile:

f means just that you are dealing with a float number.

f means just that you are dealing with a float number.

Lol, shamefull. I think now you realize how new I am to the programming world :slight_smile:

Hi again everybody,
I think it’s okay now. Someone found what was wrong :smiley:

when doing this…

//get view vector in new local coordinates system//
	v[0] = v[0]*x[0] + v[1]*x[1] + v[2]*x[2] ; // = v.x = cos (v, x) since x and v are normalized
	v[1] = v[0]*y[0] + v[1]*y[1] + v[2]*y[2] ; // = v.y = cos (v, y)
	v[2] = v[0]*n[0] + v[1]*n[1] + v[2]*n[2] ; // = v.n = cos (v, n)

… I did a noobish mistake. I forgot that, to compute v[1] for example, I need v[0] BEFORE v[0] is modified… I hope you know what I mean. So, here is the new code :

//get view vector in new local coordinates system//
	tmp_v[0] = v[0]*x[0] + v[1]*x[1] + v[2]*x[2] ; // = v.x = cos (v, x) since x and v are normalized
	tmp_v[1] = v[0]*y[0] + v[1]*y[1] + v[2]*y[2] ; // = v.y = cos (v, y)
	tmp_v[2] = v[0]*n[0] + v[1]*n[1] + v[2]*n[2] ; // = v.n = cos (v, n)
	Normalize (tmp_v) ;
	v[0] = tmp_v[0] ;
	v[1] = tmp_v[1] ;
	v[2] = tmp_v[2] ;

Now, it works. When I set Cxy to -1, and set Cz to 1, I get exactly the same as a Phong shader. It’s great someone saw the mistake :slight_smile:

I’ll post more soon, if I solve problem #3.

btw, it shouldn’t be too difficult to allow anisotropy, the only problem is that Blender’s ‘tangent’ option seems to screw things up (replacing the normal with the tangent vector etc).

For testing purposes, I helped mfoxdogg hack around the ‘tangent’ stuff in this patch . Basically the lines like

if ((ma->mode & MA_TANGENT_V)&&(ma->spec_shader != MA_SPEC_WARDISO)) {

disable the side effects of switching tangents on, if that shader is selected.

Then rather than using Vecorthobasisf(), you can just use the tangent vectors as one of your basis vectors, with a simple cross product between the tangent and the normal to give you the 3rd vector. Having a look at that patch may give you more info. Also, INPR(x,y) is a macro for the dot product, a bit easier to read :slight_smile:

Looking forward to seeing where you can go with this!

Hi broken, thanks for the reply. I also thought I have to use the tangent vectors to align x and y axis to the directions of anisotropy. I’ll try if it isn’t too hard (but it will make more buttons in the interface :slight_smile: )

Apart from this, I solved the problem of multiple lamps. It was because in my function, I modified l and v vectors directly whereas I shouldn’t touch them, since Blender probably needs them for others calculations. Every time the function was called (every pixel with a Lafortune shader), the light and view vectors were modified…
So, now, it works, I only modify tmp_l and a tmp_v vectors.

Here’s the code :

/* Lafortune isotropic spec */
static float Lafortune_Spec( float *n, float *l, float *v, int hard , float Cxy, float Cz , int tangent )
{
	float rslt ;
	float x[3] ;
	float y[3] ;
	float tmp_l[3] ;
	float tmp_v[3] ;
	
	VecOrthoBasisf(n, x, y); // find local vector basis, with surface normal = z axis
		
	/*get view vector in new local coordinates system (projecting it on x,y,n)
	(we cannot afford to modify the view vector itself, since Blender needs it for calculs.
	That's why tmp_l is used - this is also true for the light vector in following lines.) */
	tmp_v[0] = v[0]*x[0] + v[1]*x[1] + v[2]*x[2] ; // = v.x = cos (v, x) since x and v are normalized
	tmp_v[1] = v[0]*y[0] + v[1]*y[1] + v[2]*y[2] ; // = v.y = cos (v, y)
	tmp_v[2] = v[0]*n[0] + v[1]*n[1] + v[2]*n[2] ; // = v.n = cos (v, n)
	Normalize (tmp_v) ;

	//get light vector in new local coordinates system (projecting it on x,y,n)//
	tmp_l[0] = x[0]*l[0] + x[1]*l[1] + x[2]*l[2] ; // = l.x = cos (l, x) since x and l are normalized.
	tmp_l[1] = y[0]*l[0] + y[1]*l[1] + y[2]*l[2] ; // = l.y = cos (l, y)
	tmp_l[2] = n[0]*l[0] + n[1]*l[1] + n[2]*l[2] ; // = l.n = cos (l, n)
	Normalize (tmp_l) ;
		
	rslt = Cxy*tmp_l[0]*tmp_v[0] + Cxy*tmp_l[1]*tmp_v[1] + Cz*tmp_l[2]*tmp_v[2] ;
	// weighted dot product between light vector and view vector //
	if (tangent) rslt = sasqrt( 1.0f - rslt*rslt ) ; // in case tangent shading is turned on
	if( rslt > 0.0f )  rslt= spec(rslt, hard); // take the hardness parameter into account
	else rslt = 0.0f; // clamp negative values to 0
	return rslt;
	
}

I also did a quick vid with a sphere and a spinning light : http://www.vimeo.com/1444096

See how nicely stretched the highlight becomes at grazing angles… :slight_smile: It’s far more realistic than Phong and the others I find.

So, here are the source files for those who would like to play with : http://www.savefile.com/files/1701267

Check the directory.txt file to know where to put the files.

By the way, what is the easiest way for this code to become official ? Should I submit a patch ?

With a patch, it’s better : http://www.savefile.com/files/1701289

With a build, it’s even better :smiley: (for vista)
http://www.savefile.com/files/1701425