ImportanceSampling.h

// Ported from https://github.com/jMonkeyEngine/jmonkeyengine/blob/PBRisComing/jme3-core/src/main/java/com/jme3/environment/generation/PrefilteredEnvMapFaceGenerator.java
// License: MIT
// To use this with GLSL set #define GLSL
// To use this with OpenCL set #define OpenCL
#define WRAP 0
#define STRETCH 1
#define FixSeamsMethod int

#ifdef GLSL
  // Convert source from OpenCL to GLSL
    #define toFloat(X) float(X)
    #define toFloat2(X,Y) float2(X,Y)
    #define toFloat3(X,Y,Z) float3(X,Y,Z)
    #define toFloat4(X,Y,Z,W) float4(X,Y,Z,W)
    #define toInt(X) int(X)
    #define toInt2(X,Y) int2(X,Y)
    #define toInt3(X,Y,Z) int3(X,Y,Z)
    #define toInt4(X,Y,Z,W) int4(X,Y,Z,W)
    #define fabs(X) abs(X)
    #define float2 vec2
    #define float3 vec3
    #define float4 vec4
    #define int2 ivec2
    #define int3 ivec3
    #define int4 ivec4
    #define __inline 
    #define __global
    #define __read_only
    #define __write_only
    #define __kernel
    #define f(x) x
    #define  __p(x) x
    #define  __in in
    #define  __out out
    #define p(x) x
    #define r(x) x
#endif

#ifdef OpenCL
    #define toFloat(X) (float)(X)
    #define toFloat2(X,Y) (float2)(X,Y)
    #define toFloat3(X,Y,Z) (float3)(X,Y,Z)
    #define toFloat4(X,Y,Z,W) (float4)(X,Y,Z,W)
    #define toInt(X) (int)(X)
    #define toInt2(X,Y) (int2)(X,Y)
    #define toInt3(X,Y,Z) (int3)(X,Y,Z)
    #define toInt4(X,Y,Z,W) (int4)(X,Y,Z,W)
    #define f(x) x##f
    #define  __p(x) *x
    #define  __in 
    #define  __out 
    #define  p(x) (*x)
    #define r(x) &x
#endif

#ifndef PI
#define PI f(3.14159265358979323846264)
#endif

__inline float3 importanceSampleGGX(__in float4 xi, __in float a2, __in float3 normal) {
      float cosTheta = sqrt((f(1.0) - xi.x) / (f(1.0) + (a2 - f(1.0)) * xi.x));
      float sinTheta = sqrt(f(1.0) - cosTheta * cosTheta);

      float sinThetaCosPhi = sinTheta * xi.z;  // xi.z is cos(phi)
      float sinThetaSinPhi = sinTheta * xi.w;  // xi.w is sin(phi)

      float3 upVector = toFloat3(f(1.0), f(0.0), f(0.0));

      if (fabs(normal.z) < f(0.999)) {
        upVector = toFloat3(f(0.0), f(1.0), f(0.0));
      }

      float3 tangentX = normalize(cross(upVector,normal));
      float3 tangentY = cross(normal, tangentX);

      // Tangent to world space
      tangentX *= (sinThetaCosPhi);
      tangentY *= (sinThetaSinPhi);
      float3 vect5 = normal * cosTheta;

      // Tangent to world space
      float3 store = (tangentX) + (tangentY) + (vect5);

      return store;
}

// From http://holger.dammertz.org/stuff/notes_HammersleyOnHemisphere.html
 __inline float radicalInverse_VdC(__in uint bits) {
     bits = (bits << 16u) | (bits >> 16u);
     bits = ((bits & 0x55555555u) << 1u) | ((bits & 0xAAAAAAAAu) >> 1u);
     bits = ((bits & 0x33333333u) << 2u) | ((bits & 0xCCCCCCCCu) >> 2u);
     bits = ((bits & 0x0F0F0F0Fu) << 4u) | ((bits & 0xF0F0F0F0u) >> 4u);
     bits = ((bits & 0x00FF00FFu) << 8u) | ((bits & 0xFF00FF00u) >> 8u);
     return toFloat(bits) * 2.3283064365386963e-10; // / 0x100000000
 }

__inline float4 getHammersleyPoint(__in uint ui, __in int nbrSample) {
    float4 store;
    float phi;
    store.x = (toFloat(ui) / toFloat(nbrSample));
    store.y = radicalInverse_VdC(ui);
    phi = f(2.0) * PI * store.y;
    store.z = (cos(phi));
    store.w = (sin(phi));
    return store;
}

 __inline  int cubemapVectorToUVFace(__in float3 texelVect, __out float2 __p(store)) {
        float u =f(0.0);
        float v = f(0.0);
        float  bias = f(0.0);
        int face;
        float absX = fabs(texelVect.x);
        float absY = fabs(texelVect.y);
        float absZ = fabs(texelVect.z);
        float m = max(max(absX, absY), absZ);
        if (m == absX) {
            face = texelVect.x > f(0.0) ? 0 : 1;
        } else if (m == absY) {
            face = texelVect.y > f(0.0) ? 2 : 3;
        } else {
            face = texelVect.z > f(0.0) ? 4 : 5;
        }

        //compute vector depending on the face
        // Code from Nvtt : http://code.google.com/p/nvidia-texture-tools/source/browse/trunk/src/nvtt/CubeSurface.cpp	
        switch (face) {
            case 0:
                bias = f(1.0) / texelVect.x;
                u = -texelVect.z;
                v = -texelVect.y;
                break;
            case 1:
                bias = f(-1.0) / texelVect.x;
                u = texelVect.z;
                v = -texelVect.y;
                break;
            case 2:
                bias = f(1.0) / texelVect.y;
                u = texelVect.x;
                v = texelVect.z;
                break;
            case 3:
                bias = f(-1.0) / texelVect.y;
                u = texelVect.x;
                v = -texelVect.z;
                break;
            case 4:
                bias = f(1.0) / texelVect.z;
                u = texelVect.x;
                v = -texelVect.y;
                break;
            case 5:
                bias = f(-1.0) / texelVect.z;
                u = -texelVect.x;
                v = -texelVect.y;
                break;
        }
        u *= bias;
        v *= bias;

        p(store).xy=toFloat2(u,v);
        return face;
    }

__inline  float3 cubemapUvFaceToVector(__in float2 uv, __in int face) {
	  float u=uv.x;
    float v=uv.y;
		float3 store;
		//compute vector depending on the face
		// Code from Nvtt : http://code.google.com/p/nvidia-texture-tools/source/browse/trunk/src/nvtt/CubeSurface.cpp	
		switch(face){
			case 0:
				store.xyz=toFloat3(f(1.0),-v,-u);
				break;
			case 1:
				store.xyz=toFloat3(-f(1.0),-v,u);
				break;
			case 2:
				store.xyz=toFloat3(u,f(1.0),v);
				break;
			case 3:
				store.xyz=toFloat3(u,-f(1.0),-v);
				break;
			case 4:
			    store.xyz=toFloat3(u,-v,f(1.0));
				break;
			case 5:
				store.xyz=toFloat3(-u,-v,-f(1.0));
				break;            
		}

		store.xyz=normalize(store.xyz);
        return store;
	}

PrefilteredEnvMapGenerator.frag

// License: MIT
// Requires MTR.
#define GLSL
#include "ImportanceSampling.h"

uniform sampler2D m_Face1;
uniform sampler2D m_Face2;
uniform sampler2D m_Face3;
uniform sampler2D m_Face4;
uniform sampler2D m_Face5;
uniform sampler2D m_Face6;

uniform float m_Roughness;
uniform int m_Samples;

out vec4 outFace1;
out vec4 outFace2;
out vec4 outFace3;
out vec4 outFace4;
out vec4 outFace5;
out vec4 outFace6;

in vec2 texCoord;

void kernel(in int faceId,out vec4 outFace){
    float roughness=m_Roughness;
    int samples=m_Samples;

    float a2 = roughness * roughness;
    a2 *= a2;
    a2 *= 10.;

    vec3 N=cubemapUvFaceToVector(texCoord, faceId);

    float total_weight = 0.0;
    vec4 prefiltered_color = vec4(0,0,0,1);
    for (uint i = 0; i < samples; i++) {
        vec4 Xi = getHammersleyPoint(i, samples);
        vec3 H = importanceSampleGGX(Xi, a2, N);

        H = normalize(H);


        float NoH = dot(N, H);
        vec3 L = (H * (NoH * 2.)) - N;
        float NoL = clamp(dot(N, L), 0., 1.);

        if (NoL > 0.) {
            vec2 luv ; // Selected Uv inside the envemap
            int lface = cubemapVectorToUVFace(L, luv); // Selected Face inside the envmap

            vec3 c ;
            if (lface == 0) c = texture(m_Face1,  luv).xyz;
            else if (lface == 1) c = texture(m_Face2,   luv).xyz;
            else if (lface == 2) c = texture(m_Face3,  luv).xyz;
            else if (lface == 3) c = texture(m_Face4,   luv).xyz;
            else if (lface == 4) c = texture(m_Face5,   luv).xyz;
            else if (lface == 5) c = texture(m_Face6,   luv).xyz;
            else {
                outFace=vec4(1,1,0,1);
                return;
            }

            prefiltered_color.x = (prefiltered_color.x + c.r * NoL);
            prefiltered_color.y = (prefiltered_color.y + c.g * NoL);
            prefiltered_color.z = (prefiltered_color.z + c.b * NoL);
            total_weight += NoL;
        }    

    }
    prefiltered_color.rgb /= total_weight;
    outFace=prefiltered_color;
}

void main(){   
    kernel(0,outFace1);
    kernel(1,outFace2);
    kernel(2,outFace3);
    kernel(3,outFace4);
    kernel(4,outFace5);
    kernel(5,outFace6);
}

PrefilteredEnvMapGenerator.vert

in vec4 inPosition;
in vec2 inTexCoord;

out vec2 texCoord;

void main() {
    vec2 pos = inPosition.xy * 2.0 - 1.0;      
    gl_Position = vec4(pos, 0.0, 1.0);
    texCoord = inTexCoord;
}