Jump to content
  • Advertisement

Calculating Irradiance Map

Recommended Posts

I am trying to calculate irradiance map based on this article: http://www.codinglabs.net/article_physically_based_rendering.aspx

I understand the idea, you must place an hemisphere over the normal, then you need to sum every incoming radiance. Then write the radiance back at the corresponding pixel in the cube map. However I didn't understand how author iterated over the hemisphere and the direct translation from hlsl code to glsl didn't worked as I expected. It smooth the original cube map a lot. As if I dropped resolution from 1024px to 32px.


I changed the hemisphere iteration as below:


    normal = normalize(normal);
    vec3 up = vec3(0.0, 1.0, 0.0);
    vec3 right = normalize(cross(normal, up));
    int index = 0;
    vec3 irradiance = vec3(0.0 ,0.0 ,0.0);
    for (float longi = 0.0; longi <= 90.0; longi += 3.0)
     	mat4 trl = rotationMatrix(right, radians(longi));
        for (float azi = 0.0; azi <= 360.0; azi += 3.0)
         	mat4 tra = rotationMatrix(normal, radians(azi));
            vec3 sampleVec = (tra * trl * vec4(normal, 1.0)).xyz;
            irradiance += texture(iChannel0, sampleVec).rgb * dot(sampleVec, normal);

    fragColor = vec4((PI * irradiance / float(index)), 1.0);

Generated irradiance map seems to me too bright. Also I don't understand why we are averaging the summed radiance by dividing it to "index" ? Aren't we after the total incoming radiance to a point ?

Here is the link to the shader.



Edited by afraidofdark

Share this post

Link to post
Share on other sites
6 hours ago, afraidofdark said:

I didn't understand how author iterated over the hemisphere


Now I understand how the author traversed over the hemisphere. However there is still a problem. In monte carlo integration, you consider the integration out put as the average of random samples. This is why we divide the sum by "index". However, this part isn't clear to me yet. If I want to approximate x^2 integral over 0 to 5 period, does this method suggest me to take for example 100 random number in that interval and average them, and expect this outcome as my integration result ?

Share this post

Link to post
Share on other sites
37 minutes ago, afraidofdark said:

If I want to approximate x^2 integral over 0 to 5 period, does this method suggest me to take for example 100 random number in that interval and average them

I have lived enough enlightenment today thanks to scratchapixel ! I've seen the light ! (Literally) Since I used a cube map to look up incoming radiance, I tie my hands to use an approximation of the actual integral. Because for incoming radiance, I don't have a function which I can integrate over spherical coordinates !

Share this post

Link to post
Share on other sites

So the full monte carlo formula is essentially (1 / NumSamples) * Sum(f(x) / p(x), NumSamples), where f(x) is the function you're integrating, and p(x) is probability distribution function (PDF) evaluated for x. The PDF is basically the likelihood of a particular sample being chosen. For "uniform" sampling schemes where samples are distributed evenly across the whole domain (for instance, the surface of a sphere) the PDF is the same for all samples, and so you can pull it out of the sum.

Now for integrating irradiance for a point with normal N, you'll need integrate over the surface of the hemisphere that surrounds N. Uniformly sampling a hemisphere oriented around Z = 1 (the canonical "upper" hemisphere) is pretty simple:

// Returns a direction on the hemisphere around z = 1
Float3 SampleDirectionHemisphere(float u1, float u2)
    float z = u1;
    float r = std::sqrt(std::max(0.0f, 1.0f - z * z));
    float phi = 2 * Pi * u2;
    float x = r * std::cos(phi);
    float y = r * std::sin(phi);

    return Float3(x, y, z);

As input you need u1 and u2, which are two random variables in the range [0, 1]. These can be generated from a random number generator like rand(), or can be generated by a sequence that produces more a optimal sample distribution that avoids clustering (this is known as Quasi-Monte Carlo sampling).

To make this work for our case, we need to sample the hemisphere surrounding the surface normal. This means we need a transformation that can go from tangent space -> world space. If you're doing normal mapping then you probably already have such a matrix, and if not you can build one by starting with the normal as your Z basis and then generating a perpendicular vector to use as the X basis (after that you can use a cross product to generate the Y basis). Then you just transform the hemisphere sample direction by this matrix, and you now have a world space direction.Now if you recall from earlier, we need to divide samples by the PDF of each sample. For uniform hemisphere sampling the PDF is constant for all samples, and it ends up being one over the surface area of a unit hemisphere (1 / (2 * Pi)). So our algorithm ends up looking like this:

float3 irradiance = 0.0f;
for(uint i = 0; i < NumSamples; ++i)
    float u1 = RandomFloat();
    float u2 = RandomFloat();
    float3 sampleDirInTangentSpace = SampleDirectionHemisphere(u1, u2);
    float3 sampleDirInWorldSpace = mul(sampleDirInTangentSpace, tangentToWorld);
    float3 radiance = SampleEnvironment(sampleDirInWorldSpace);
    irradiance += radiance * saturate(sampleDirInWorldSpace, normalInWorldSpace);

float hemispherePDF = 1.0f / (2 * Pi);
irradiance /= (NumSamples * hemispherePDF);

Just be aware that if you use this result to light a surface using a Lambertian BRDF, make sure that you include the 1 / Pi term that's part of that BRDF. A lot of people like to bake the 1/ Pi into their light sources which is fine, you just have to be careful to make sure that you're doing it somewhere otherwise your diffuse lighting can be too bright.

Once you have this working, you can improve the quality by being smarter about choosing your sample directions. The first way to do that is to use QMC techniques like I mentioned earlier, for instance using Halton sequence or using stratified sampling. You can also just pick values that are evenly spaced out over the [0, 1] domain, which is basically stratified sampling without any jitter (this is essentially what you're doing in the code you posted above). The other way to is to importance sample your function by choosing samples that match the shape of the function being sampled. For irradiance, the common way to do this is to choose sample directions that are proportional to the cosine of the angle between that direction and the surface normal.

Share this post

Link to post
Share on other sites

Create an account or sign in to comment

You need to be a member in order to leave a comment

Create an account

Sign up for a new account in our community. It's easy!

Register a new account

Sign in

Already have an account? Sign in here.

Sign In Now

  • Advertisement
  • Advertisement
  • Popular Tags

  • Similar Content

    • By Colm
      Hey all, I've been trying to work out how LittleBigPlanet handles its objects for a while now.
      For those unaware, LittleBigPlanet has a building component where you can build 2D-ish (there are 2 - 16 2D layers that you can build on) objects. There are a number of shaped brushes to do this with, from basic squares and circles to teardrops and eye shapes. There's a decent video showing this off, actually.
      Anyways, I've been trying to work out how this works for a while now. My current thought is that it might be along the lines of storing a list of object corners and then drawing an object within those bounds - this makes the most sense to me because the engine has a corner editor for making more advanced shapes, and because some of the restrictions in the engine are based around corners. Of course, that could also be completely wrong and it's something else entirely.
      What are your thoughts?
    • By ethancodes
      I'm wondering if anyone has any examples or tips/ideas on how to handle pickups in a game. The game is an arkanoid style game. I'm going to have at least 5 different pick up types, and they are going to be in a queue, where only one can be active at a time. Once one pick up is expended, the next one should automatically start up. Some of the pick ups have an immediate effect, such as ball speed. Others will activate when the pickup is hit, but doesn't actually do anything until the ball hits the paddle. Those type of pick ups have a limited number of shots instead of a time limit. What I'm trying to figure out is what kind of structure I should have for a system like this? I'm just curious how these things are handled in other games, especially if anyone has any examples that would be great. Thank you!
    • By Mihumihu
      I'm trying to solve a problem where I can get all colors from image. 
      I see only one: walk through a loop at raster data and collect all bytes, I'm sure there is a better way to get colors from image.  I'm thinking about some sort  of collecting colors in result texture...
      Is it ordinary situation, could you help me, I didn find anithing on the internet...
    • By Thomas Harper
      Hi Gamedev team,
      I am trying to re-implement the fluid simulation from the following paper . http://movement.stanford.edu/courses/cs448-01-spring/papers/foster.pdf. 
      However, I am stuck with the level set method as it generates very strange resultt. The water should move much faster than the level set actually propagates, so I think there is a bug in my levelSetPropagation. However, I tried now for several months, but I just can't get it right. In the lower image you can see the problem. The scene is a breaking dam. While the inserted particles move quite fast , my level set seems to stick to the wall and does not propagate correctly. As the particles are moving correctly, I assume the error lies in my level set propagation.  The level set sticks to the wall very long which looks really strange and not correct. I also inserted my level set propagation code below. If anybody could give me a hint on what I am doing wrong please let me know. If you require more information just let me know and I will reply as soon as possible. thank you for your support!

    • By scullsold
      Hi I read some tutorials on STL as many people on this forum say it's faster than the most selfwritten container classes and somehow i can't even declare a list in VC .net 2003...the compiler says: Compiling... VertexManager.cpp c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(33) : error C2143: syntax error : missing ';' before '<' c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(33) : error C2501: 'CVertexManager::list' : missing storage-class or type specifiers c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(33) : error C2238: unexpected token(s) preceding ';' c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(34) : error C2143: syntax error : missing ';' before '<' c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(34) : error C2501: 'CVertexManager::list' : missing storage-class or type specifiers c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(34) : error C2238: unexpected token(s) preceding ';' c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(35) : error C2143: syntax error : missing ';' before '<' c:\Dokumente und Einstellungen\gregor\Eigene Dateien\nGin\VertexManager.h(35) : error C2501: 'CVertexManager::list' : missing storage-class or type specifiers my code: class CVertexManager { private: list<VertListEntry> VertList; list<VertGroup> VertGroup; list<int> TextureChange; CVertexManager(); public: ~CVertexManager(); void addEntry(VertListEntry Entry); static CVertexManager& Instance() { static CVertexManager TheOneAndOnly; return CVertexManager; } }; btw what does the list.insert function want as the first parameter? it says something with iterator...what is that? can i just pass an int as the number where i want to have the new entry? regards, m4gnus  
    • By Nikita Sidorenko
      I'm making render just for fun (c++, opengl)
      Want to add decals support. Here what I found
      A couple of slides from doom
      http://advances.realtimerendering.com/s2016/Siggraph2016_idTech6.pdf Decals but deferred 
      http://martindevans.me/game-development/2015/02/27/Drawing-Stuff-… space-Decals/
      No implementation details here
      As I see there should be a list of decals for each tile same as for light sources. But what to do next?
      Let assume that all decals are packed into a spritesheet. Decal will substitute diffuse and normal.
      - What data should be stored for each decal on the GPU? 
      - Articles above describe decals as OBB. Why OBB if decals seem to be flat?
      - How to actually render a decal during object render pass (since it's forward)? Is it projected somehow? Don't understand this part completely.
      Are there any papers for this topic?
    • By Ming-Lun "Allen" Chou
      Here is the original blog post.
      Edit: Sorry, I can't get embedded LaTeX to display properly.
      The pinned tutorial post says I have to do it in plain HTML without embedded images?
      I actually tried embedding pre-rendered equations and they seemed fine when editing, 
      but once I submit the post it just turned into a huge mess.
      So...until I can find a proper way to fix this, please refer to the original blog post for formatted formulas.
      I've replaced the original LaTex mess in this post with something at least more readable.
      Any advice on fixing this is appreciated.
      This post is part of my Game Math Series.
      Source files are on GitHub.
      Shortcut to sterp implementation.
      Shortcut to code used to generate animations in this post.
      An Alternative to Slerp
      Slerp, spherical linear interpolation, is an operation that interpolates from one orientation to another, using a rotational axis paired with the smallest angle possible.
      Quick note: Jonathan Blow explains here how you should avoid using slerp, if normalized quaternion linear interpolation (nlerp) suffices. Long store short, nlerp is faster but does not maintain constant angular velocity, while slerp is slower but maintains constant angular velocity; use nlerp if you’re interpolating across small angles or you don’t care about constant angular velocity; use slerp if you’re interpolating across large angles and you care about constant angular velocity. But for the sake of using a more commonly known and used building block, the remaining post will only mention slerp. Replacing all following occurrences of slerp with nlerp would not change the validity of this post.
      In general, slerp is considered superior over interpolating individual components of Euler angles, as the latter method usually yields orientational sways.
      But, sometimes slerp might not be ideal. Look at the image below showing two different orientations of a rod. On the left is one orientation, and on the right is the resulting orientation of rotating around the axis shown as a cyan arrow, where the pivot is at one end of the rod.

      If we slerp between the two orientations, this is what we get:

      Mathematically, slerp takes the “shortest rotational path”. The quaternion representing the rod’s orientation travels along the shortest arc on a 4D hyper sphere. But, given the rod’s elongated appearance, the rod’s moving end seems to be deviating from the shortest arc on a 3D sphere.
      My intended effect here is for the rod’s moving end to travel along the shortest arc in 3D, like this:

      The difference is more obvious if we compare them side-by-side:

      This is where swing-twist decomposition comes in.
      Swing-Twist Decomposition
      Swing-Twist decomposition is an operation that splits a rotation into two concatenated rotations, swing and twist. Given a twist axis, we would like to separate out the portion of a rotation that contributes to the twist around this axis, and what’s left behind is the remaining swing portion.
      There are multiple ways to derive the formulas, but this particular one by Michaele Norel seems to be the most elegant and efficient, and it’s the only one I’ve come across that does not involve any use of trigonometry functions. I will first show the formulas now and then paraphrase his proof later:
      Given a rotation represented by a quaternion R = [W_R, vec{V_R}] and a twist axis vec{V_T}, combine the scalar part from R the projection of vec{V_R} onto vec{V_T} to form a new quaternion: T = [W_R, proj_{vec{V_T}}(vec{V_R})]. We want to decompose R into a swing component and a twist component. Let the S denote the swing component, so we can write R = ST. The swing component is then calculated by multiplying R with the inverse (conjugate) of T: S= R T^{-1} Beware that S and T are not yet normalized at this point. It's a good idea to normalize them before use, as unit quaternions are just cuter. Below is my code implementation of swing-twist decomposition. Note that it also takes care of the singularity that occurs when the rotation to be decomposed represents a 180-degree rotation. public static void DecomposeSwingTwist ( Quaternion q, Vector3 twistAxis, out Quaternion swing, out Quaternion twist ) { Vector3 r = new Vector3(q.x, q.y, q.z); // singularity: rotation by 180 degree if (r.sqrMagnitude < MathUtil.Epsilon) { Vector3 rotatedTwistAxis = q * twistAxis; Vector3 swingAxis = Vector3.Cross(twistAxis, rotatedTwistAxis); if (swingAxis.sqrMagnitude > MathUtil.Epsilon) { float swingAngle = Vector3.Angle(twistAxis, rotatedTwistAxis); swing = Quaternion.AngleAxis(swingAngle, swingAxis); } else { // more singularity: // rotation axis parallel to twist axis swing = Quaternion.identity; // no swing } // always twist 180 degree on singularity twist = Quaternion.AngleAxis(180.0f, twistAxis); return; } // meat of swing-twist decomposition Vector3 p = Vector3.Project(r, twistAxis); twist = new Quaternion(p.x, p.y, p.z, q.w); twist = Normalize(twist); swing = q * Quaternion.Inverse(twist); } Now that we have the means to decompose a rotation into swing and twist components, we need a way to use them to interpolate the rod’s orientation, replacing slerp.
      Swing-Twist Interpolation
      Replacing slerp with the swing and twist components is actually pretty straightforward. Let the Q_0 and Q_1 denote the quaternions representing the rod's two orientations we are interpolating between. Given the interpolation parameter t, we use it to find "fractions" of swing and twist components and combine them together. Such fractiona can be obtained by performing slerp from the identity quaternion, Q_I, to the individual components. So we replace: Slerp(Q_0, Q_1, t) with: Slerp(Q_I, S, t) Slerp(Q_I, T, t) From the rod example, we choose the twist axis to align with the rod's longest side. Let's look at the effect of the individual components Slerp(Q_I, S, t) and Slerp(Q_I, T, t) as t varies over time below, swing on left and twist on right:
      And as we concatenate these two components together, we get a swing-twist interpolation that rotates the rod such that its moving end travels in the shortest arc in 3D. Again, here is a side-by-side comparison of slerp (left) and swing-twist interpolation (right):

      I decided to name my swing-twist interpolation function sterp. I think it’s cool because it sounds like it belongs to the function family of lerp and slerp. Here’s to hoping that this name catches on.
      And here’s my code implementation:
      public static Quaternion Sterp ( Quaternion a, Quaternion b, Vector3 twistAxis, float t ) { Quaternion deltaRotation = b * Quaternion.Inverse(a); Quaternion swingFull; Quaternion twistFull; QuaternionUtil.DecomposeSwingTwist ( deltaRotation, twistAxis, out swingFull, out twistFull ); Quaternion swing = Quaternion.Slerp(Quaternion.identity, swingFull, t); Quaternion twist = Quaternion.Slerp(Quaternion.identity, twistFull, t); return twist * swing; } Proof
      Lastly, let’s look at the proof for the swing-twist decomposition formulas. All that needs to be proven is that the swing component S does not contribute to any rotation around the twist axis, i.e. the rotational axis of S is orthogonal to the twist axis. Let vec{V_{R_para}} denote the parallel component of vec{V_R} to vec{V_T}, which can be obtained by projecting vec{V_R} onto vec{V_T}: vec{V_{R_para}} = proj_{vec{V_T}}(vec{V_R}) Let vec{V_{R_perp}} denote the orthogonal component of vec{V_R} to vec{V_T}: vec{V_{R_perp}} = vec{V_R} - vec{V_{R_para}} So the scalar-vector form of T becomes: T = [W_R, proj_{vec{V_T}}(vec{V_R})] = [W_R, vec{V_{R_para}}] Using the quaternion multiplication formula, here is the scalar-vector form of the swing quaternion: S = R T^{-1} = [W_R, vec{V_R}] [W_R, -vec{V_{R_para}}] = [W_R^2 - vec{V_R} ‧ (-vec{V_{R_para}}), vec{V_R} X (-vec{V_{R_para}}) + W_R vec{V_R} + W_R (-vec{V_{R_para}})] = [W_R^2 - vec{V_R} ‧ (-vec{V_{R_para}}), vec{V_R} X (-vec{V_{R_para}}) + W_R (vec{V_R} -vec{V_{R_para}})] = [W_R^2 - vec{V_R} ‧ (-vec{V_{R_para}}), vec{V_R} X (-vec{V_{R_para}}) + W_R vec{V_{R_perp}}] Take notice of the vector part of the result: vec{V_R} X (-vec{V_{R_para}}) + W_R vec{V_{R_perp}} This is a vector parallel to the rotational axis of S. Both vec{V_R} X(-vec{V_{R_para}}) and vec{V_{R_perp}} are orthogonal to the twist axis vec{V_T}, so we have shown that the rotational axis of S is orthogonal to the twist axis. Hence, we have proven that the formulas for S and T are valid for swing-twist decomposition. Conclusion
      That’s all.
      Given a twist axis, I have shown how to decompose a rotation into a swing component and a twist component.
      Such decomposition can be used for swing-twist interpolation, an alternative to slerp that interpolates between two orientations, which can be useful if you’d like some point on a rotating object to travel along the shortest arc.
      I like to call such interpolation sterp.
      Sterp is merely an alternative to slerp, not a replacement. Also, slerp is definitely more efficient than sterp. Most of the time slerp should work just fine, but if you find unwanted orientational sway on an object’s moving end, you might want to give sterp a try.
    • By Hanseul Shin
      Thanx to @Randy Gaul, I succesfully implemented cube/cube collision detection and response.
      1- substract the center of each AABB = 3d vector a.
      2- if |x| of a is the biggest, this represents a face on each AABB.
      3- if x is pointing at the same(or exact opposte) direction of the normal(of a face), two AABB are colliding on those faces.
      But these steps only work if two colliders are cubes, because the size of each half-lengths are different in a right square prism.
      I'd like to check which faces are collided with two right square prism, please help!
      Thank you!
    • By StepperDox
      I've been digging around online and can't seem to find any formulas for 3D mesh simplification. I'm not sure where to start but I generally want to know how I could make a function that takes in an array of vertices, indices, and a float/double for the decimation rate. And could I preserve the general shape of the object too?
      Thanks for the help!
      P.S. I was hoping to do something with Quadric Error / Quadric Edge Collapse if that's possible.
    • By Armando Neto
      INTRODUCTION - It's the following, for a university chair I'm making an application that will be an accessory to assist in face-to-face RPG tables, it has rooms in which the master and players can manage tokens, I just want to do one more thing, a map where all the players in the room could see and move their "miniatures" in it the master could put their monsters, but without animation, it would be just to replace the paper map that is normally used in face-to-face sessions.
      PROBLEM - To do the map explained in the introduction, I looked for some tools but I did not find any good, so I would ask to be told APIs to do this, if someone has already played table RPG knows more or less how I want to do, I will leave an image of the roll 20 that is more or less similar, and some more images of how I want it to stay on the cell phone for those who never played understand , I just need to know a good tool to do this.
      I want the master himself to put the image that will be used when starting the map (photo or maybe tiled map).
      It can be solutions for web to and call the website into application, literally anything that does what I want will help.
      If you can help me, I'll be grateful, sorry for the awful English.

  • Advertisement
  • Popular Now

  • Forum Statistics

    • Total Topics
    • Total Posts

Important Information

By using GameDev.net, you agree to our community Guidelines, Terms of Use, and Privacy Policy.

Participate in the game development conversation and more when you create an account on GameDev.net!

Sign me up!