球谐光照总结与实现
实习需要学习了一下球谐光照,总结如下:
1. 首先是下面的rendering equation
第二部分是一个半球面上的积分。
2. 球谐变换
类比傅里叶变换[采用定义在圆上的三角函数],球谐变换采用定义在球面上的一组球谐基函数。[Spherical Harmonic Lighting: The Gritty Details]有这些函数的介绍。
现在只需要知道
- 一个定义在球面上的函数可以近似分解成一组正交基函数的加权和:F(s) ≈ ∑ fi Yi(s).
- 两个函数f(s),g(s)乘积的积分可以用两个基向量的点积近似表示:∫f(s)g(s)ds ≈∑ fi gi.
- 坐标fi的求解,类别向量的点积,函数的点积为积分:fi=∫F(s)Yi(s)d
这样我们就把一个球面上的积分转换成了计算向量的点积,问题是如何求解坐标fi
3. 蒙特卡罗积分
解fi需要用到神奇的蒙特卡罗积分
-------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
下面可以写代码了,看一下最简单的“dot product lighting”,也就是漫反射模型,对应的渲染方程里面的brdf部分f(p,wo,wi)是一个常数,可以拿出积分号外,现在积分号里面剩下了∫Li(p,wi)*max(n.dot(wi),0)dwi。根据2.2,我们只需要分别投影分解Li(p,wi)和max(n.dot(wi),0)dwi两个函数,利用蒙特卡洛方法求积分的时候还可以加上可见性的判断。注意需要对模型的每一个顶点都需要计算一组基向量fi,如果加上自遮挡的判断会耗费大量时间,例如下面的dragon用了6000s来处理所有的模型顶点。
--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
贴一下代码
1. SHLightingMesh类 负责模型顶点的处理并打包绑定到gpu.
class SHLightingMesh : public glMeshWrapper { GLuint* vbos_4_shcoefs_gpu; float9* buff_4_shcoefs_cpu; public: SHLightingMesh(const redips::Triangles* model, const SamplesOnUnitSphere& samples, redips::ShaderSource shaderSource = redips::ShaderSource()) : glMeshWrapper(model, shaderSource, WrapOption::_genMipmap_){ bindVaoAttribData(0, 1, 2); if (shaderSource.sourceType == redips::ShaderSource::SourceType::_null_) { char strbuf[512]; sprintf_s(strbuf, "%s/OpenglWrappers/DemoMeshes/SHLighting", _REDIPS_ROOT_PATH_); useShader(strbuf); } //allocate memory vbos_4_shcoefs_gpu = new GLuint[mesh->groups.size()]; { glGenBuffers(mesh->groups.size(), vbos_4_shcoefs_gpu); int maxFaceCnt = 0; { for (auto& g : mesh->groups) maxFaceCnt = std::max(maxFaceCnt, g.faceCnt); _RUNTIME_ASSERT_(maxFaceCnt > 0, "assert [maxFaceCnt>0] failed"); buff_4_shcoefs_cpu = new float9[maxFaceCnt * 3]; } } //build-tree, accelerte intersection calculation const_cast<Triangles*>(model_ptr())->buildTree(); //calculate sh coefficients const std::vector<float3>& normals = mesh->normals; const std::vector<float3>& vertices = mesh->vertices; for (int i = 0; i < mesh->groups.size(); ++i) { auto& g = mesh->groups[i]; if (g.faceCnt <= 0) continue; _RUNTIME_ASSERT_(g.faceType >= GROUP_FACE_TYPE::_withnormal_, "assert [g.faceType>= GROUP_FACE_TYPE::_withnormal_] failed"); //fill cpu array float9* ptr_f9 = buff_4_shcoefs_cpu; int nid = g.nsid, fid = g.fsid, faceCnt = g.faceCnt; for (int neid = nid + faceCnt; nid < neid; ++nid,++fid) { const int3& indices_n = mesh->faces_vn[nid]; const int3& indices_v = mesh->faces_v [fid]; CalculateSHCoefficients(vertices[indices_v.x], normals[indices_n.x].unit(), samples, *ptr_f9++); CalculateSHCoefficients(vertices[indices_v.y], normals[indices_n.y].unit(), samples, *ptr_f9++); CalculateSHCoefficients(vertices[indices_v.z], normals[indices_n.z].unit(), samples, *ptr_f9++); } //copy to gpu and bind to vertex glBindVertexArray(vaos[i]); glBindBuffer(GL_ARRAY_BUFFER, vbos_4_shcoefs_gpu[i]); glBufferData(GL_ARRAY_BUFFER, sizeof(float) * 9 * faceCnt * 3, buff_4_shcoefs_cpu, GL_STATIC_DRAW); glVertexAttribPointer(3, 3, GL_FLOAT, GL_FALSE, sizeof(float) * 9, NULL); glEnableVertexAttribArray(3); glVertexAttribPointer(4, 3, GL_FLOAT, GL_FALSE, sizeof(float) * 9, (const void*)(sizeof(float)*3)); glEnableVertexAttribArray(4); glVertexAttribPointer(5, 3, GL_FLOAT, GL_FALSE, sizeof(float) * 9, (const void*)(sizeof(float)*6)); glEnableVertexAttribArray(5); }
delete buff_4_shcoefs_cpu; }; SHLightingMesh(const glMeshWrapper& another, redips::ShaderSource shaderSource = redips::ShaderSource()) = delete; ~SHLightingMesh() { glDeleteBuffers(mesh->groups.size(), vbos_4_shcoefs_gpu); } void draw() { if (!m_shader) { std::cerr << "shader error" << std::endl; return; }; m_shader->Use(); for (int i = 0; i < meshCnt; i++) { unsigned int flags = 0; if (meshMtls[i]->texture_ka != NULL && meshFaceTypes[i] == redips::GROUP_FACE_TYPE::_withtex_) { glActiveTexture(GL_TEXTURE0); glBindTexture(GL_TEXTURE_2D, mtlTextureHandle[meshMtls[i]->texture_ka]); glUniform1i(glGetUniformLocation(m_shader->Program, shaderAmbientTextureUniformStr), shaderAmbientTextureLocation); flags |= 1u; } { redips::float3 color = meshMtls[i]->ambient; glUniform3f(glGetUniformLocation(m_shader->Program, shaderAmbientColorUniformStr), color.x, color.y, color.z); } if (meshMtls[i]->texture_kd != NULL && meshFaceTypes[i] == redips::GROUP_FACE_TYPE::_withtex_) { glActiveTexture(GL_TEXTURE1); glBindTexture(GL_TEXTURE_2D, mtlTextureHandle[meshMtls[i]->texture_kd]); glUniform1i(glGetUniformLocation(m_shader->Program, shaderDiffuseTextureUniformStr), shaderDiffuseTextureLocation); flags |= 2u; } { redips::float3 color = meshMtls[i]->diffuse; glUniform3f(glGetUniformLocation(m_shader->Program, shaderDiffuseColorUniformStr), color.x, color.y, color.z); } glUniform1ui(glGetUniformLocation(m_shader->Program, shaderSurfaceTypeUniformStr), flags); glBindVertexArray(vaos[i]); glDrawArrays(GL_TRIANGLES, 0, 3 * meshFaceCnt[i]); } } private: void CalculateSHCoefficients(const float3& vertice, const float3& normal, const SamplesOnUnitSphere& samples, float9& result) { result = 0; redips::HitRecord hitRecord; for (int i = 0; i < samples.size(); ++i) { auto fs_i = std::max(samples[i].dot(normal.unit())*-1, 0.0f); if (fs_i > 0.0) { Ray ray(vertice, normal); if (!const_cast<Triangles*>(model_ptr())->intersect(ray, hitRecord)) { result = result + samples.samples_vec9[i] * fs_i; } } } result = result * (4 * PI / samples.size()); } const GLchar* shaderSurfaceTypeUniformStr = "material.flags"; const GLchar* shaderAmbientColorUniformStr = "material.ambient"; const GLchar* shaderDiffuseColorUniformStr = "material.diffuse"; const GLchar* shaderAmbientTextureUniformStr = "material.ambientTexture"; const GLchar* shaderDiffuseTextureUniformStr = "material.diffuseTexture"; const GLuint shaderAmbientTextureLocation = 0; const GLuint shaderDiffuseTextureLocation = 1; };
2. SamplesOnUnitSphere 类 负责在单位球面上随机生成均匀的采样点
class SamplesOnUnitSphere { public: std::vector<float9> samples_vec9; std::vector<redips::float3> samples_vec3; int size() const{ return samples_vec3.size(); } SamplesOnUnitSphere(int sampleCnt) { std::default_random_engine generator; std::uniform_real_distribution<float> normal_floats(0, 1); samples_vec3.resize(sampleCnt); samples_vec9.resize(sampleCnt); for (int i = 0; i < sampleCnt; ++i) { auto x = normal_floats(generator); auto y = normal_floats(generator); auto theta = acos(sqrt(x)) * 2 - PI * 0.5; auto phi = (2 * y - 1)*PI; samples_vec3[i] = redips::float3(cos(theta)*sin(phi), sin(theta), cos(theta)*cos(phi)); samples_vec9[i] = float9(samples_vec3[i]); } } redips::float3& operator[](size_t idx) { _RUNTIME_ASSERT_(idx < samples_vec3.size(), "assert [idx<samples_vec3.size()] failed"); return samples_vec3[idx]; } const redips::float3& operator[](size_t idx) const { _RUNTIME_ASSERT_(idx < samples_vec3.size(), "assert [idx<samples_vec3.size()] failed"); return samples_vec3[idx]; } };
3. Cubemap类 里面的computeSH函数负责将来自六个面的环境光编码成球谐基向量
class Cubemap { redips::FImage* faces[6]; //a cube with side-length 2 centered at (0,0,0) redips::Triangles cube = redips::float3(2, 2, 2); public: Cubemap(const char* picnames[6]) { for (int i = 0; i < 6; ++i) faces[i] = new redips::FImage(picnames[i]); } ~Cubemap() { for (int i = 0; i < 6; ++i) delete faces[i]; } redips::float3 texture(const redips::float3& direction) { using namespace redips; //find the intersection redips::Ray ray({ 0,0,0 }, direction); redips::HitRecord recorder; cube.intersect(ray, recorder); auto hitp = ray.ori + ray.dir * recorder.distance; //fabs(one-dim-of-$hitp$) must equal 1, because $hitp$ is on one face of $cube$ int FaceId = -1; float2 texcoord; const float epsi = 1e-5; if (fabs(hitp.x - 1) < epsi) { FaceId = 0; texcoord = float2(hitp.z, hitp.y)*0.5 + float2{0.5f, 0.5f}; } else if (fabs(hitp.x + 1) < epsi) { FaceId = 1; texcoord = float2(hitp.z*-1, hitp.y)*0.5 + float2{ 0.5f, 0.5f }; } else if (fabs(hitp.y - 1) < epsi) { FaceId = 2; texcoord = float2(hitp.x, hitp.z)*0.5 + float2{ 0.5f, 0.5f }; } else if (fabs(hitp.y + 1) < epsi) { FaceId = 3; texcoord = float2(hitp.x, hitp.z*-1)*0.5 + float2{ 0.5f, 0.5f }; } else if (fabs(hitp.z - 1) < epsi) { FaceId = 4; texcoord = float2(hitp.x*-1, hitp.y)*0.5 + float2{ 0.5f, 0.5f }; } else if (fabs(hitp.z + 1) < epsi) { FaceId = 5; texcoord = float2(hitp.x, hitp.y)*0.5 + float2{ 0.5f, 0.5f }; } if (FaceId < 0) { std::cerr << "[Cubemap] : Ops, something wrong in Cubemap.texture" << std::endl; return redips::float3(); } return faces[FaceId]->tex2d(texcoord.x, texcoord.y); } void computeSH(SamplesOnUnitSphere& samples, float9& channel_r, float9& channel_g, float9& channel_b) { channel_r = 0, channel_g = 0, channel_b = 0; for (int i = 0; i < samples.size(); ++i) { auto Es_i = texture(samples[i]*-1); auto& Ys = samples.samples_vec9[i]; channel_r = channel_r + Ys*Es_i.x; channel_g = channel_g + Ys*Es_i.y; channel_b = channel_b + Ys*Es_i.z; } auto scale = (4 * PI / samples.size()); { channel_r = channel_r * scale; channel_g = channel_g * scale; channel_b = channel_b * scale; } } };
4. main.cpp负责组装
#include <common/glfwApp.h> #include <OpenglWrappers/glHelper.h> #include <OpenglWrappers/DemoMeshes/Skybox.h> #include <OpenglWrappers/DemoMeshes/SHLightingMesh.h> //#define MODEL_PATH "E:/Documents/models/spheres/sphere0.obj" #define MODEL_PATH "E:/Documents/models/dragon.obj" //window size const redips::int2 Winsize{ 960,960 }; //glfw setup auto application = redips::glfw::getInstance(Winsize.x, Winsize.y); //screen capture auto screenCapture = redips::glScreenCapture::getInstance(Winsize); //skybox const char* skybox_faces[] = { "D:/Documents/resources/skyboxes/pure/black.jpg", "D:/Documents/resources/skyboxes/pure/blue.jpg", "D:/Documents/resources/skyboxes/pure/black.jpg", "D:/Documents/resources/skyboxes/pure/red.jpg", "D:/Documents/resources/skyboxes/pure/black.jpg", "D:/Documents/resources/skyboxes/pure/green.jpg" }; auto skybox = redips::SkyBox(skybox_faces); //standard pinhole camera redips::PhC phc(60, 1.0f, 0.1f, 10000); //load a obj and then wrap into a glMesh。 redips::SphericalHarmonic::SHLightingMesh * shmesh = nullptr; //mesh center redips::float3 heart = redips::float3{ 0,0,0 }; //keyboard event void movement() { if (application->keydown(GLFW_KEY_M)) { screenCapture->capture("capture.bmp", 0); } } //need register a display-callback function to tell glfwApp what to render void display() { using namespace redips; shmesh->uniformMat44f("model", shmesh->model_ptr()->Transform().transpose().ptr()); shmesh->uniformMat44f("projection_view", phc.glProjectionView().ptr()); shmesh->draw(); skybox.render(&phc); movement(); } void initialize() { using namespace redips; //setup camera phc.lookAt(heart + float3(0, 0, 20), heart, float3(0, 1, 0)); phc.setResolution(Winsize.x, Winsize.y); //generate samples SphericalHarmonic::SamplesOnUnitSphere samples(500); //calculate light SphericalHarmonic::Cubemap cubemap(skybox_faces); SphericalHarmonic::float9 light_r, light_g, light_b; cubemap.computeSH(samples, light_r, light_g, light_b); //warp mesh shmesh = new redips::SphericalHarmonic::SHLightingMesh(new Triangles(MODEL_PATH), samples); shmesh ->littedBy(light_r, light_g, light_b);
//transform mesh const_cast<Triangles*>(shmesh->model_ptr())->setTransform( Mat44f::scale(float3(60,60,60))* //Mat44f::rotatey(RAD(90))* Mat44f::translation(shmesh->model_ptr()->aabb_R().heart()*-1) ); } int main() { initialize(); application->registerDisplayCallback(display); application->bindCamera(&phc); application->loop(); return 0; }
redips 是之前封装的基础操作的类库