球谐函数实现环境光照漫反射实践

news2024/10/6 4:06:03

该文章以及代码主要来自
图形学论文解析与复现:【论文复现】An Efficient Representation for Irradiance Environment Maps
作者:Monica的小甜甜

与原文的不同

  • 对一些有问题的地方进行了修改
  • 添加了注释
  • 对有疑问的地方添加了疑问点
  • 引入了其他一些Blog填充了原文中忽略的信息

1、预计算球面谐波函数系数

首先根据上一篇【球谐函数在环境光照中的使用原理】得到的最终公式:

在这里插入图片描述
我们需要预计算 L l m L_l^m Llm的值。计算公式为:
在这里插入图片描述

Ω \Omega Ω为球面积分,这里对应对天空盒逐像素积分。
积分代码为:

void Harmonics::Evaluate()//求值
{
	m_Coefs = vector<glm::vec3>(m_Degree, glm::vec3());
	
	//6张图
	for (int k = 0; k < 6; k++)
	{
		cv::Mat img = m_Images[k];
		int w = m_Images[k].cols;
		int h = m_Images[k].rows;
		//逐像素
		for (int j = 0; j < w; j++)
		{
			for (int i = 0; i < h; i++)
			{
				// 像素点位置
				float px = (float)i + 0.5;
				float py = (float)j + 0.5;
				// 像素点UV 【-1,1】:以摄像机正对位置的(0,0)
				float u = 2.0 * (px / (float)w) - 1.0;
				float v = 2.0 * (py / (float)h) - 1.0;
				// 像素间UV的一半的偏移量
				float d_x = 1.0 / (float)w;
				// (x0,y0)像素左下角 (x1,y1)像素右上角
				float x0 = u - d_x;
				float y0 = v - d_x;
				float x1 = u + d_x;
				float y1 = v + d_x;
				// 计算Cubemap的一个像素对应的立体角的大小
				float d_a = surfaceArea(x0, y0) - surfaceArea(x0, y1) - surfaceArea(x1, y0) + surfaceArea(x1, y1);
				// 纹理像素点 转化为 世界坐标点
				u = (float)j / (img.cols - 1);
				v = 1.0f - (float)i / (img.rows - 1);
				glm::vec3 p = CubeUV2XYZ({ k, u, v });
				// 获取当前像素颜色
				auto c = img.at<cv::Vec3f>(i, j);
				glm::vec3 color = {c[2], c[1], c[0]};
				// 得到基函数计算结果列表
				vector<float> Y = Basis(p);
				// 计算系数
				for (int i = 0; i < m_Degree; i++)
				{
					m_Coefs[i] = m_Coefs[i] + Y[i] * color * d_a;
				}
			}
		}
	}
}

其中 计算Cubemap的一个像素对应的立体角的大小原理可参照
Solid Angle of A Cubemap Texel - 计算Cubemap的一个像素对应的立体角的大小

我们将得到的积分结果保存在一个文件中【SHCoefficients.txt】,用于之后读取。

2、预计算BRDF的LUT图

LUT(Look up Table)图,预计算了任意一个天空盒下,已知法线和视口的夹角以及材质粗糙度,查找得到Frenel项

然而这个LUT图和IBL中的LUT有一些不同。
因为IBL中的LUT加入了 n ⋅ w n\cdot w nw 光照衰减项。
而在球谐函数中, n ⋅ w n\cdot w nw 作为 t l 参与运算 t_l参与运算 tl参与运算,因此在球谐函数的IBL中删除了 n ⋅ w n\cdot w nw

main函数计算

	for(int i = 0; i < N; i++){
		for (int j = 0; j < N; j++)
		{
			float NoV = (i + 0.5f) * (1.0f / N);
			float roughness = (j + 0.5f) * (1.0f / N);
			glm::vec2 eval = IntegrateBRDF(NoV, roughness);
			tex.store<glm::vec2>({ i, N - j - 1 }, 0, eval);
		}
	}

其他被调用函数

const float PI = 3.14159265358979323846264338327950288;

float RadicalInverse_VdC(unsigned int 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 float(bits) * 2.3283064365386963e-10;
}

glm::vec2 Hammersley(unsigned int i, unsigned int N)
{
	return glm::vec2(float(i) / float(N), RadicalInverse_VdC(i));
}

glm::vec3 ImportanceSampleGGX(glm::vec2 Xi, float roughness, glm::vec3 N)
{
	float a = roughness * roughness;

	float phi = 2.0 * PI * Xi.x;
	float cosTheta = sqrt((1.0 - Xi.y) / (1.0 + (a*a - 1.0) * Xi.y));
	float sinTheta = sqrt(1.0 - cosTheta * cosTheta);

	// from spherical coordinates to cartesian coordinates
	glm::vec3 H;
	H.x = cos(phi) * sinTheta;
	H.y = sin(phi) * sinTheta;
	H.z = cosTheta;

	// from tangent-space vector to world-space sample vector
	glm::vec3 up = abs(N.z) < 0.999 ? glm::vec3(0.0, 0.0, 1.0) : glm::vec3(1.0, 0.0, 0.0);
	glm::vec3 tangent = normalize(cross(up, N));
	glm::vec3 bitangent = cross(N, tangent);

	glm::vec3 sampleVec = tangent * H.x + bitangent * H.y + N * H.z;
	return normalize(sampleVec);
}

float GeometrySchlickGGX(float NdotV, float roughness)
{
	float a = roughness;
	float k = (a * a) / 2.0;

	float nom = NdotV;
	float denom = NdotV * (1.0 - k) + k;

	return nom / denom;
}

float GeometrySmith(float roughness, float NoV, float NoL)
{
	float ggx2 = GeometrySchlickGGX(NoV, roughness);
	float ggx1 = GeometrySchlickGGX(NoL, roughness);

	return ggx1 * ggx2;
}

glm::vec2 IntegrateBRDF(float NdotV, float roughness, unsigned int samples = 1024)
{
	glm::vec3 V;
	V.x = sqrt(1.0 - NdotV * NdotV);
	V.y = 0.0;
	V.z = NdotV;

	float A = 0.0;
	float B = 0.0;

	glm::vec3 N = glm::vec3(0.0, 0.0, 1.0);

	for (unsigned int i = 0u; i < samples; ++i)
	{
		glm::vec2 Xi = Hammersley(i, samples);
		glm::vec3 H = ImportanceSampleGGX(Xi, roughness, N);
		glm::vec3 L = normalize(2.0f * dot(V, H) * H - V);

		float NoL = glm::max(L.z, 0.0f);
		float NoH = glm::max(H.z, 0.0f);
		float VoH = glm::max(dot(V, H), 0.0f);
		float NoV = glm::max(dot(N, V), 0.0f);

		if (NoL > 0.0)
		{
			float G = GeometrySmith(roughness, NoV, NoL);
			float G_Vis = (G * VoH) / (NoH * NoV) / NoL;
			float Fc = pow(1.0 - VoH, 5.0);

			A += (1.0 - Fc) * G_Vis;
			B += Fc * G_Vis;
		}
	}

	return glm::vec2(A / float(samples), B / float(samples));
}

在这里插入图片描述

3、将计算数据传入Shader

  • 传入BRDFLUT纹理
  • 传入球谐函数系数列表
void CShadingPass::initV()
{

	auto m_LUTTexture = std::make_shared<ElayGraphics::STexture>();
	loadTextureFromFile("../Textures/BRDFLUT/BRDFLut.dds", m_LUTTexture);
	getCoefs();
	ElayGraphics::Camera::setMainCameraFarPlane(100);
	ElayGraphics::Camera::setMainCameraPos({ -1.57278, 0.244948, 0.367194 });
	ElayGraphics::Camera::setMainCameraFront({ 0.967832, -0.112856, -0.224865 });
	ElayGraphics::Camera::setMainCameraMoveSpeed(0.5);
	m_pShader = std::make_shared<CShader>("Sponza_VS.glsl", "Sponza_FS.glsl");
	m_pSponza = std::dynamic_pointer_cast<CSponza>(ElayGraphics::ResourceManager::getGameObjectByName("Sponza"));
	m_pShader->activeShader();
	m_pShader->setTextureUniformValue("u_BRDFLut", m_LUTTexture);
	m_pShader->setMat4UniformValue("u_ModelMatrix", glm::value_ptr(m_pSponza->getModelMatrix()));
	for (int i = 0; i < m_Coefs.size(); i++)
	{
		m_pShader->setFloatUniformValue("u_Coef[" + std::to_string(i) + "]", m_Coefs[i].x, m_Coefs[i].y, m_Coefs[i].z);
	}
	m_pSponza->initModel(*m_pShader);
}

4、 Draw

#version 430 core

in  vec3 v2f_FragPosInViewSpace;
in  vec2 v2f_TexCoords;
in  vec3 v2f_ViewSpaceNormal;
in  vec3 v2f_WorldSpaceNormal;

layout (location = 0) out vec4 Albedo_;

const float PI = 3.1415926535897932384626433832795;
uniform vec3 u_Coef[16];
uniform vec3 u_DiffuseColor;
uniform sampler2D u_BRDFLut;


vec3 FresnelSchlickRoughness(float cosTheta, vec3 F0, float roughness)
{
    return F0 + (max(vec3(1.0 - roughness), F0) - F0) * pow(max(1.0 - cosTheta, 0.0), 5.0);
}  

void main()
{	
	if((abs(v2f_ViewSpaceNormal.x) < 0.0001f) && (abs(v2f_ViewSpaceNormal.y) < 0.0001f) && (abs(v2f_ViewSpaceNormal.z) < 0.0001f))
	{
		Albedo_ = vec4(0, 0, 0, 1);
		return;
	}

	float Basis[9];
	float x = v2f_WorldSpaceNormal.x;
	float y = v2f_WorldSpaceNormal.y;
	float z = v2f_WorldSpaceNormal.z;
	float x2 = x * x;
	float y2 = y * y;
	float z2 = z * z;
    
	//这里所有系数应该为乘PI------------------个人认为
    Basis[0] = 1.f / 2.f * sqrt(1.f / PI);
    Basis[1] = 2.0 / 3.0 * sqrt(3.f / 4.f * PI) * z;
    Basis[2] = 2.0 / 3.0 * sqrt(3.f / 4.f * PI) * y;
    Basis[3] = 2.0 / 3.0 * sqrt(3.f / 4.f * PI) * x;
	
    Basis[4] = 1.0 / 4.0 * 1.f / 2.f * sqrt(15.f * PI) * x * z;
    Basis[5] = 1.0 / 4.0 * 1.f / 2.f * sqrt(15.f * PI) * z * y;
    Basis[6] = 1.0 / 4.0 * 1.f / 4.f * sqrt(5.f * PI) * (-x2 - z2 + 2 * y2);
    Basis[7] = 1.0 / 4.0 * 1.f / 2.f * sqrt(15.f * PI) * y * x;
    Basis[8] = 1.0 / 4.0 * 1.f / 4.f * sqrt(15.f * PI) * (x2 - z2);

	vec3 Diffuse = vec3(0,0,0);

	vec3 F0 = vec3(0.2,0.2,0.2);
	float Roughness = 0.5;
	vec3 N = normalize(vec4(v2f_ViewSpaceNormal,1.0f)).xyz;//viewMatrix * 
	vec3 V = -normalize(v2f_FragPosInViewSpace);
	//vec3 R = reflect(-V, N); 
	F0        = FresnelSchlickRoughness(max(dot(N, V), 0.0), F0, Roughness);
	vec2 EnvBRDF  = texture(u_BRDFLut, vec2(max(dot(N, V), 0.0), Roughness)).rg;
	vec3 LUT = (F0 * EnvBRDF.x + EnvBRDF.y);

	for (int i = 0; i < 9; i++)
		Diffuse += u_Coef[i] * Basis[i] * (1-LUT);
	Albedo_ = vec4(Diffuse);
}

结果展示

在这里插入图片描述
在这里插入图片描述
只有漫反射的效果:
在这里插入图片描述
只有镜面反射的效果:
在这里插入图片描述

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/989609.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

【笔试强训选择题】Day37.习题(错题)解析

作者简介&#xff1a;大家好&#xff0c;我是未央&#xff1b; 博客首页&#xff1a;未央.303 系列专栏&#xff1a;笔试强训选择题 每日一句&#xff1a;人的一生&#xff0c;可以有所作为的时机只有一次&#xff0c;那就是现在&#xff01;&#xff01; 文章目录 前言一、Day…

camx camera initial

qnx 平台中的camera hal 接口 HAL3Module&#xff1a;chi_hal_override_entry 在android 的中使用Camx 打开com.qti.chi.override.so进行注册hal ops 操作接口 camhal3module.cpp 中的构造函数HAL3Module中 CHIHALOverrideEntry funcCHIHALOverrideEntry reinterpret_cast…

即拼七人拼团系统开发模式的奖励机制都有哪些?

即拼七人拼团是市场上非常火爆的商业模式之一&#xff0c;它通过团购机制、互动社交和抽奖等方式&#xff0c;有效解决了电商平台的复购难题。不仅可以降低商品价格&#xff0c;还能够增加用户的参与感和购物乐趣&#xff0c;提升平台的用户粘性和产品销量。下面就来具体说一下…

驱动轴相机参数设置Web前端界面开发

一、基于Django的Web应用界面的开发&#xff1a; 在Realtimeresults.html上添加一个按钮组件&#xff0c;获取检测到的轴型和车轮信息&#xff0c;点击后可以获取package.json里存放的json数据&#xff0c;效果如下&#xff1a; 实现逻辑&#xff1a;需要从URL设置、视图函数、…

YOLOV8从零搭建一套目标检测系统(修改model结构必看)附一份工业缺陷检测数据集

目录 1.YOLOV8介绍 2.YOLOV8安装 2.1环境配置 3.数据集准备 1.YOLOV8介绍 Yolov8结构图&#xff1a; YoloV8相对于YoloV5的改进点&#xff1a; Replace the C3 module with the C2f module. Replace the first 6x6 Conv with 3x3 Conv in the Backbone. Delete two Convs …

Mysql--事务

事务 开始之前&#xff0c;让我们先想一个场景&#xff0c;有的时候&#xff0c;为了完成某个工作&#xff0c;需要完成多种sql操作 比如转账 再比如下单 第一步 我的账户余额减少 第二步 商品的库存要减少 第三步 订单表中要新增一项 事务的本质&#xff0c;就是为了把多个操…

Excel数学、工程和科学计算插件:FORMULADESK Studio

如果 Excel 是您的武器 - 让我们磨砺您的剑&#xff01;为整天使用 Excel 的人们提供创新的 Excel 加载项&#xff0c;你需要这个 FORMULADESK Studio。。。 Excel 插件为任何使用 Excel 执行数学、工程和科学计算的人提供了必备工具。 * 将公式视为真正的数学方程 * 为您的公…

vue3中如何实现通过点击不同的按钮切换不同的页面

完成以上需求&#xff0c;我们可以使用vue中的component标签来实现。 component是Vue.js中一个特殊的标签&#xff0c;用于动态地绑定其它组件。它可以与v-bind:is指令一起使用&#xff0c;来决定要渲染哪个组件。下面是示例代码 <template><div class"app-conte…

OpenCV 07(图像滤波器)

一、卷积 什么是图片卷积? 图像卷积就是卷积核在图像上按行滑动遍历像素时不断的相乘求和的过程 步长 步长就是卷积核在图像上移动的步幅. 上面例子中卷积核每次移动一个像素步长的结果, 如果将这个步长修改为2, 结果会如何? 为了充分扫描图片, 步长一般设为1. padding …

【操作系统】电脑上没有IIS怎么办

文章目录 前言一、查看二、解决 前言 有的新机刚开始在计算机-管理-服务下没有IIS网络服务怎么办。 一、查看 桌面计算机/此电脑 鼠标右键&#xff1a;管理 服务和应用 发现没有IIS 二、解决 控制面板 程序和功能 启动或关闭Windows功能 IIS相关的所有功能选中&#xff…

【JavaScript】JS语法入门到实战

文章目录 一、初识JavaScript1. 什么是JavaScript&#xff1f;2. JavaScript 和 HTML 和 CSS 之间的关系3. JavaScript的运行过程4. JavaScript的组成 二、JavaScript的书写形式三、变量1. 输入输出2. 变量的使用3. 数据类型 四、运算符五、分支和循环语句1. 分支语句2. 循环语…

将PyCharm中的终端运行前面的PS修改成当前环境

最近使用Pycharm中的Terminal来pip安装一些pakage&#xff0c;发现Terminal运行前面的显示的是PS&#xff0c;然后输入安装指令报错。“python无法将“pip”项识别为 cmdlet、函数、脚本文件或可运行程序的名称。” 解决方法&#xff1a; 只需要在pycharm的设置中修改一些termi…

Java缓存理解

CPU占用&#xff1a;如果你有某些应用需要消耗大量的cpu去计算&#xff0c;比如正则表达式&#xff0c;如果你使用正则表达式比较频繁&#xff0c;而其又占用了很多CPU的话&#xff0c;那你就应该使用缓存将正则表达式的结果给缓存下来。 数据库IO性能&#xff1a;如果发现有大…

基于YOLOv8和WiderFace数据集的人脸目标检测系统(PyTorch+Pyside6+YOLOv8模型)

摘要&#xff1a;基于YOLOv8和WiderFace数据集的人脸目标检测系统可用于日常生活中检测与定位人脸目标&#xff0c;利用深度学习算法可实现图片、视频、摄像头等方式的目标检测&#xff0c;另外本系统还支持图片、视频等格式的结果可视化与结果导出。本系统采用YOLOv8目标检测算…

数据结构(C语言版)概念、数据类型、线性表

数据结构&#xff08;C语言&#xff09;基本概念 数据的基本单位 数据的基本单位是位&#xff08;bit&#xff09;和字节&#xff08;byte&#xff09;。位是最小的存储单位&#xff0c;它可以表示一个二进制的0或1。字节由8个位组成&#xff0c;用于表示一个字符或数字。在计…

STM32 Nucleo-144开发板开箱bring-up

文章目录 1. 开篇2. 开发环境搭建2.1 下载官方例程2.2 ST-Link安装 3. STM32F446ZE demo工程3.1 STM32F446ZE简介3.2 跑个demo试一试 1. 开篇 最近做项目&#xff0c;用到STM32F446ZET6这款MCU&#xff0c;为了赶进度&#xff0c;前期软件需要提前开发&#xff0c;于是在某宝买…

基于FPGA的RGB图像转Ycbcr实现,包括tb测试文件以及MATLAB辅助验证

目录 1.算法运行效果图预览 2.算法运行软件版本 3.部分核心程序 4.算法理论概述 5.算法完整程序工程 1.算法运行效果图预览 将FPGA的数据导入到matlab进行显示 2.算法运行软件版本 Vivado2019.2 matlab2022a 3.部分核心程序 timescale 1ns / 1ps // // Company: // E…

企业微信后台管理-关联小程序、H5/web

企业微信后台管理-小程序、web 应用-关联小程序应用-绑定h5(web页面的app)企业微信-工作台 应用-关联小程序 企业微信后台管理地址&#xff0c;管理员扫码登录后台管理&#xff0c;找到应用管理-自建-创建应用/小程序。 填写项目相关信息之后&#xff0c;如下图。 点击关联小…

Cpp/Qtday030908cpp基础

目录 目录 自行封装一个栈的类&#xff0c;包含私有成员属性&#xff1a;栈的数组、记录栈顶的变量 成员函数完成&#xff1a;构造函数、析构函数、拷贝构造函数、入栈、出栈、清空栈、判空、判满、获取栈顶元素、求栈的大小 头文件&#xff1a;stack.h 源文件: stack.cp…

数据结构之队列的实现(附源码)

目录 一、队列的概念及结构 二、队列的实现 拓展&#xff1a;循环队列 三、初学的队列以及栈和队列结合的练习题 一、队列的概念及结构 队列&#xff1a;只允许在一端进行插入数据操作&#xff0c;在另一端进行删除数据操作的特殊线性表&#xff0c;队列具有先进先出FIFO(Fi…