|
Today I was reading some code for atmospheric lighting. Glancing over the code, I noticed this: (I did rename a couple of things).
+ Show Spoiler +// As the y tex coord goes from 0 to 1, the angle goes from 0 to 180 degrees float fCos = 1 - (nAngle + nAngle) / (float)mHeight; float fAngle = std::acos(fCos); vec3 ray(std::sin(fAngle), std::cos(fAngle), 0); for(int nHeight=0; nHeight < nSize; nHeight++) { // As the x tex coord goes from 0 to 1, the height goes from the bottom of the atmosphere to the top float fHeight = delta + innerRadius + ((outerRadius - innerRadius) * nHeight) / mWidth; // The position of the camera vec3 eyepos(0, fHeight, 0);
// If the ray from vPos heading in the vRay direction intersects the inner // radius (i.e. the planet), then this spot is not visible from the viewpoint float B = 2.0f * dot(eyepos, ray); float Bsq = B * B; float Cpart = dot(eyepos, eyepos); float C = Cpart - innerRadius*innerRadius; float discr = Bsq - 4.0f * C; bool bVisible = (fDet < 0 || (0.5f * (-B - std::sqrt(discr)) <= 0) && (0.5f * (-B + std::sqrt(discr)) <= 0));
What it does is pretty obvious, first it finds some cosine and sine to a particular angle, then applies the quadratic formula to solve the intersection of a ray with the planet. If you know something about trigonometry, this will bother you a lot :
float fAngle = acosf(fCos); vec3 ray(std::sin(fAngle), std::cos(fAngle), 0);
It does seem a little strange to call an inverse cosine only to take it's cosine again the next instruction, and of course we can replace it by fCos. Secondly, we can eliminate the acos call entirely by using the pythagorean identity :
you can easily see how you can solve this to get sin(angle)
float fCos = 1 - (nAngle + nAngle) / (float)mHeight; vec3 ray(std::sqrt(1 - fCos*fCos), fCos, 0);
secondly is the quadratic formula, to solve the intersection of the ray with the sphere. What most people forget to do is to remove all the 2, 4, 0.5 coefficients which will cancel each other anyway. we get:
float B = dot(eyepos, ray); float Bsq = B * B; float Cpart = dot(eyepos, eyepos); float C = Cpart - innerRadius*innerRadius; float discr = Bsq - C; bool bVisible = discr < 0 || -B - std::sqrt(discr) <= 0) && -B + std::sqrt(discr) <= 0; Now you might look at the bVisible boolean expression and wonder if it also can't be simplified any further, and yes indeed it is possible. In particular, the factors between the && signs: We start by letting D = sqrt(discr) and multiplying both expressions by -1 and we get:
D + B >= 0 && -D + B >= 0
suppose the case that the left hand side is false, e.g D + B <= 0, if we subtract 2D from both sides of this equation we get
-D + B <= -2D <= 0
,which means the right hand side will also be false. since -D + B < D + B, when the right hand side IS true, the lhs will also be true.
this allows us to write the expression as :
bVisible = discr < 0 || B >= std::sqrt(discr);
Hope you enjoyed reading it and learned something
|
I think you're smarter than me.
|
28076 Posts
On March 08 2014 16:28 -Kaiser- wrote: I think you're smarter than me.
Makes two of us
|
On March 08 2014 10:35 fmod wrote:float fCos = 1 - (nAngle + nAngle) / (float)mHeight; Make that a static_cast<float>(mHeight) and I am happy
Btw replacing any sine/cosine functions from the standard library is always a good idea. I’ve read that they are pretty buggy if the argument gets too high, but I haven’t checked myself.
|
I'm a computer science major and I have no idea what's going on here *hangs head in shame*
|
I'd be interested in seeing how much of a difference this makes in terms of performance in addition to cleaning up the code. You save yourself several calculations here by making the code simpler.
Did you profile the code at all to compare its performance before and after improving it? It would be cool to see that by applying some more math to the problem, you're not only able to improve the code readability and clarity, but also the performance.
|
the thing that interested me most was how you managed to impelement TeX code here. Then I found out it was a pasted pic Edit: not that the other stuff is uninteresting. math is a useful skill toi have.
|
Since you already have Bsq you could save yourself one more call to std::sqrt() by writing the last expression as
bVisible = discr < 0 || Bsq >= discr;
|
On March 09 2014 03:13 Mr. Wiggles wrote: I'd be interested in seeing how much of a difference this makes in terms of performance in addition to cleaning up the code. You save yourself several calculations here by making the code simpler.
Did you profile the code at all to compare its performance before and after improving it? It would be cool to see that by applying some more math to the problem, you're not only able to improve the code readability and clarity, but also the performance.
We can estimate the theoretical relative performance decrease/increase by looking into the computational complexity of each of the operations and functions. For example, suppose the functions and operations like addition, sin, square root, etc., follow the complexity linked here:
http://en.wikipedia.org/wiki/Computational_complexity_of_mathematical_operations
In the original code, the first block that we looked into called 3 functions/operations: arccos, sin, cos. This has a combined time complexity of 3 * O(log(n)^2 * M(n)). In complexity arithmetic, we don't care about constants so the block has a complexity of O(log(n)^2 * M(n)). Don't worry about what M(n) is; we're just comparing, not computing.
In the modified code, the first block now calls several functions/operations with equivalent complexity: subtraction x2, addition, division, and square root. This has a combined time complexity of O(M(n)).
Not even knowing which multiplication algorithm the program/compiler is using, we can then argue that the modified code is faster by a factor of log(n)^2. This is quite small I would say. However, this is under A LOT of assumptions, chiefly that the compiler is so dumb such that it couldn't figure this optimization out on its own and that the algorithms used for trigonometric functions is equivalent to that used by the program/compiler.
|
On March 09 2014 04:38 Cyx. wrote:Since you already have Bsq you could save yourself one more call to std::sqrt() by writing the last expression as bVisible = discr < 0 || Bsq >= discr; That'll only work for positive B though. The square root is assumed to be positive always. If not the early out behaviour of || will assure it never gets there,
BirdKiller you are right that these are very 'micro' optimizations. I mostly posted this because it interested me and maybe others as well. One thing though is that by avoiding (inverse) trig functions the result will usually be more accurate, apart from being faster.
|
|
|
|