diff --git a/CG/GAMES101/README.md b/CG/GAMES101/README.md new file mode 100644 index 0000000..aaa98cf --- /dev/null +++ b/CG/GAMES101/README.md @@ -0,0 +1,11 @@ +这里会整理GAMES101的各个assignments的代码以及相关的知识 + +* [Assignment1](assignment1.md) +* [Assignment2](assignment2.md) +* [Assignment3](assignment3.md) +* [Assignment4](assignment4.md) +* [Assignment5](assignment5.md) +* [Assignment6](assignment6.md) +* [Assignment7](assignment7.md) + +games101的assignment更新到此为止,assignment8是物理模拟相关的,暂时先不更了~ \ No newline at end of file diff --git a/CG/GAMES101/assignment1.md b/CG/GAMES101/assignment1.md new file mode 100644 index 0000000..8e67d99 --- /dev/null +++ b/CG/GAMES101/assignment1.md @@ -0,0 +1,149 @@ +--- +title: GAMES101-assignment1-画三角形 +date: 2021-10-24 +tags: [GAMES101] +--- +# Assignment1 + +## 运行结果 + +关于三角形颠倒的解决方法,请见 [结果出现上下颠倒,怎么办?只加一个负号即可 – 计算机图形学与混合现实研讨会 (games-cn.org)](http://games-cn.org/forums/topic/结果出现上下颠倒,怎么办?只加一个负号即可/) + +image-20211024134106868 + +## 实现细节 + +详细代码请见:[详细代码](https://github.com/LJHG/GAMES101-assignments) + +实现了main.cpp中的**get_model_matrix**, **get_rotation**, **get_projection_matrix**函数。 + +其中每一个的实现基本都能在下面的相关知识找到对应,无非就是把矩阵用代码写出来了而已。 + +这里记录一个需要注意的点,就是在**get_projection_matrix**函数中,在设置**n**和**f**(近平面和远平面)时,需要加一个负号: + +```cpp +float n = -zNear; +float f = -zFar; +``` + + + + +# 相关知识 + + + +## 关于三维的变换 + +### 一般的旋转变换 + +[旋转变换(一)旋转矩阵_Frank的专栏-CSDN博客_旋转变换矩阵](https://blog.csdn.net/csxiaoshui/article/details/65446125) + +### 绕任意轴的旋转变换 + +image-20211017163952680 + +## MVP + +### M: Model transformation + +模型变化,**应该**就是指模型的变化。包括平移、旋转、放缩等等操作。 + +在本次的作业中指的就是绕z轴旋转,这对应一个变换矩阵。 + + + +### V: View transformation + +视图变换,指的就是从相机出发观测物体,那么就需要将相机移动到原点,并且面朝-z方向,向上方向为y轴正方向。 + +将相机移动后,我们需要也对物体做一个移动,来保持相对位置不变。所以最后的结果就是,我们直接假定相机在原点,面朝-z方向,向上方向为y轴正方向,然后**只对物体做这个变换**就行了。 + +image-20211017141934274 + +该变换分为两步: + +这里我们把视图变换的矩阵令为 $$M_{view}$$ ,平移矩阵令为 $$T_{view}$$ ,旋转矩阵令为 $$R_{view}$$ + +1. 将相机移动到原点 + + 平移矩阵为 $$T_{view}$$ + +image-20211017135851689 + +2. 将相机旋转到面朝-z方向,且上方为y正方向。 + + 直接对$$R_{view}$$ 进行求解比较困难,所以我们转而去求它的逆变换,即将处于原点,面向-z,上朝y的相机旋转到当前的相机位置,用$$R_{view}^{-1}$$表示,那么有: + +image-20211017142628026 + +​ 因为旋转矩阵是一个正交矩阵,所以矩阵的逆等于矩阵的转置,因此我们可以求得$$R_{view}$$: + +image-20211017142711735 + +因为模型变换和视图变化都是在物体上进行变化,所以一般我们把他们一起叫做**模型视图变换**。 + + + +### P:Perspective projection + +#### 正交投影(Orthographic projection) + +正交投影就是在空间中确定一个立方体 (用来表示相机的可视区域),对应六个面在不同坐标轴下的距离为right, left, top, bottom, near, far,要做的事情就是把这个立方体移动到原点,并且缩放到 (-1, 1)。 + +image-20211017144539084 + +对应的矩阵$$M_{ortho}$$ 为: + +image-20211017144634725 + +#### 透视投影(Perspective projection) + +##### 矩阵推导 + +透视投影就是将正交投影里面的一个立方体变成了一个forstum,我们要做的就是先将这个forstum "挤" 成一个立方体,然后后续操作就和正交投影一样了。 + + + +其中,透视投影变换矩阵为 $$M_{persp}$$ ,我们只需要去求透视投影向正交投影的变换矩阵 $$M_{persp->ortho}$$ 即可 + +image-20211017145049085 + +$$M_{persp->ortho}$$ 的求法比较发杂,下面贴一些图: + +* 对于当$$M_{persp->ortho}$$ 作用到(x,y,z,1)后,对于frostum里的每一个点,可以得到: + +image-20211017145920777 + +* 根据这个,我们已经可以得出矩阵中的大部分元素 + +image-20211017145933930 + +* 再根据下面标红的两个发现(近平面的点不会变,远平面的z不会变)列出公式,求解最后剩余的? + +image-20211017145948598 + +image-20211017150723015 + +image-20211017150746852 + +image-20211017150800394 + +于是就求解得到$$M_{persp->ortho}$$ 了。 + + + +##### fov与aspect ratio + +fov是(field of view),即竖直方向的可视角度 + +aspect ratio是width/height + +只需要使用fov和aspect ratio就可以用来表示一个透视投影(n和f是给定的,所以我们再把fov和aspect ratio转换为l, r, b, t即可) + + + +image-20211017155422637 + +image-20211017155555652 + diff --git a/CG/GAMES101/assignment2.md b/CG/GAMES101/assignment2.md new file mode 100644 index 0000000..02bbd7f --- /dev/null +++ b/CG/GAMES101/assignment2.md @@ -0,0 +1,86 @@ +--- +title: GAMES101-assignment2-ZBuffuer +date: 2021-11-11 +tags: [GAMES101] +--- +# Assignment2 + +## 运行结果 + +image-20220101175340666 + + + + + +## 实现细节 + +[详细代码](https://github.com/LJHG/GAMES101-assignments) +我认为这次的assignment的重点有三个: + +1. 三角形内的判断 +2. 三角形边界的确定 +3. 深度buffer + + + +### 1. 三角形内的判断 + +判断一个点是否在三角形内是通过叉乘来判断的,假设一个点为**P**,我们与三角形的三个顶点相连,并分别令为**PA**,**PB**,**PC**,那么需要计算**PAPB**, **PBPC**, **PCPA**的叉乘,只要三个结果同号,那么该点就在三角形内。这里需要注意一下叉乘的顺序,倘若把**PCPA**写为**PAPC**,结果就会大不相同。 + +```cpp +static bool insideTriangle(int x, int y, const Vector3f* _v) +{ + //TODO : Implement this function to check if the point (x, y) is inside the triangle + + //suppose the coordinate of p is (x,y) + + x = x + 0.5f; + y = y + 0.5f; + std::pair pa = std::make_pair(_v[0][0]-x, _v[0][1]-y); + std::pair pb = std::make_pair(_v[1][0]-x, _v[1][1]-y); + std::pair pc = std::make_pair(_v[2][0]-x, _v[2][1]-y); + + float papb = cross_product(pa,pb); + float pbpc = cross_product(pb,pc); + float pcpa = cross_product(pc,pa); + + if( (papb > 0 && pbpc > 0 && pcpa >0) || (papb < 0 && pbpc < 0 && pcpa < 0) ) + return true; + return false; +} +``` + + + +### 2. 三角形边界的确定 + +image-20211025231838941 + +通过对三角形的边界进行确认,可以加速光栅化的过程,具体代码如下: + +```cpp +//确定bounding box的坐标 +int left = std::min(v[0][0],std::min(v[1][0],v[2][0])); +int bottom = std::min(v[0][1],std::min(v[1][1],v[2][1])); +int right = std::max(v[0][0],std::max(v[1][0],v[2][0])) + 1 ; //向上取整 +int top = std::max(v[0][1],std::max(v[1][1],v[2][1])) + 1; +``` + + + +### 3. 深度buffer + +代码中维护了一个**depth_buf**,如果当前要绘制的像素离摄像机近,那么才绘制,否则不绘制。具体的相关代码也很简单: + +```cpp +int index = get_index(x,y); +if(z_interpolated < depth_buf[index] ){ // if near + Eigen::Vector3f point; + point << x,y,z_interpolated; + Eigen::Vector3f color; + set_pixel(point,t.getColor()); + depth_buf[index] = z_interpolated; +} +``` + diff --git a/CG/GAMES101/assignment3.md b/CG/GAMES101/assignment3.md new file mode 100644 index 0000000..bbc0c81 --- /dev/null +++ b/CG/GAMES101/assignment3.md @@ -0,0 +1,150 @@ +--- +title: GAMES101-assignment3-Blinn-Phongf以及Shader +date: 2022-01-01 +tags: [GAMES101] +--- +# Assignment3 + +两个月没更新了,来填个坑 + +## 1. 运行结果 + +这里就贴一个texture shader运行的结果: + +image-20220101165837476 + +## 2. 实现细节 + +[详细代码](https://github.com/LJHG/GAMES101-assignments) + +个人认为这次作业的主要有三个: + +1. 重心插值算法 +2. Blinn-Phong 光照模型 +3. shader以及texture的相关设计 + + + +### 2.1 重心插值算法 + +这个其实没啥好说的,代码是现成的,就是在 `rasterizer.cpp` 里的 `computeBarycentric2D` 函数,最后返回的是三个值,对应三个顶点的权重。 + +```cpp +static std::tuple computeBarycentric2D(float x, float y, const Vector4f* v){ + float c1 = (x*(v[1].y() - v[2].y()) + (v[2].x() - v[1].x())*y + v[1].x()*v[2].y() - v[2].x()*v[1].y()) / (v[0].x()*(v[1].y() - v[2].y()) + (v[2].x() - v[1].x())*v[0].y() + v[1].x()*v[2].y() - v[2].x()*v[1].y()); + float c2 = (x*(v[2].y() - v[0].y()) + (v[0].x() - v[2].x())*y + v[2].x()*v[0].y() - v[0].x()*v[2].y()) / (v[1].x()*(v[2].y() - v[0].y()) + (v[0].x() - v[2].x())*v[1].y() + v[2].x()*v[0].y() - v[0].x()*v[2].y()); + float c3 = (x*(v[0].y() - v[1].y()) + (v[1].x() - v[0].x())*y + v[0].x()*v[1].y() - v[1].x()*v[0].y()) / (v[2].x()*(v[0].y() - v[1].y()) + (v[1].x() - v[0].x())*v[2].y() + v[0].x()*v[1].y() - v[1].x()*v[0].y()); + return {c1,c2,c3}; +} +``` + +然后我们就可以对三角形内部的点进行颜色、法线、以及各种坐标的插值。 + + + +### 2.2 Blinn-Phong 光照模型 + +根据课程中给出的公式,我们可以计算 $$ L_d $$ , $$ L_s $$, $$ L_a $$ + + +$$ +L_d = k_d(I/r^2)max(0,n\cdot l) +$$ + +$$ +L_s = k_s(I/r^2)max(0,n\cdot h)^p +$$ + +$$ +L_a = k_aI_a +$$ + +对应代码如下: + +```cpp +for (auto& light : lights) + { + // TODO: For each light source in the code, calculate what the *ambient*, *diffuse*, and *specular* + // components are. Then, accumulate that result on the *result_color* object. + float r = calDistance_p2p(point,light.position); //计算光源到点的距离 + Eigen::Vector3f l = (light.position - point).normalized(); + Eigen::Vector3f v = (eye_pos - point).normalized(); + Eigen::Vector3f h = (l + v).normalized(); //计算半程向量 + L_d += kd.cwiseProduct(light.intensity)*MAX(0,normal.dot(l))/(r*r); + L_s += ks.cwiseProduct(light.intensity)*pow(MAX(0,normal.dot(h)),8)/(r*r); //p设置为8 + } + L_a = ka.cwiseProduct(amb_light_intensity); + result_color = L_d + L_s + L_a; +``` + + + +### 2.3 shader以及texture的相关设计 + +这次作业涉及到shader的设计,尤其是涉及到贴图时,个人感觉代码有点难看懂,所以在这里梳理一下。 + +首先,从 `main.cpp` 出发,我们要初始化一个 `rasterizer`: + +```cpp +rst::rasterizer r(700, 700); +``` + +我们需要为这个 `rasterizer` 指定相关的贴图: + +```cpp +auto texture_path = "hmap.jpg"; +r.set_texture(Texture(obj_path + texture_path)); +``` + +同时为它指定相关的 `shader`, 此处的`active_shader`是一个函数,输入是`fragment_shader_payload` ,输出是`Eigen::Vector3f`, 也就是颜色值。我们可以实现不同的 `active_shader` 来达到不同的效果。 + +```cpp +r.set_fragment_shader(active_shader); +``` + +传入了对应的`shader函数`后,我们会在 `rasterizer.cpp` 中,绘制三角形时取颜色值时用到: + +```cpp +auto interpolated_color = interpolate(alpha,beta,gamma,t.color[0],t.color[1],t.color[2],1); +auto interpolated_normal = interpolate(alpha,beta,gamma,t.normal[0],t.normal[1],t.normal[2],1); +auto interpolated_texcoords = interpolate(alpha,beta,gamma,t.tex_coords[0],t.tex_coords[1],t.tex_coords[2],1); +auto interpolated_viewpos = interpolate(alpha,beta,gamma,view_pos[0],view_pos[1],view_pos[2],1); +fragment_shader_payload payload( interpolated_color, interpolated_normal.normalized(), interpolated_texcoords, texture ? &*texture : nullptr); +payload.view_pos = interpolated_viewpos; +auto pixel_color = fragment_shader(payload); +Eigen::Vector2i point; +point << x,y; +set_pixel(point,pixel_color); +depth_buf[index] = z_interpolated; +``` + +我们只要传入一个对应的 `payload` ,就可以通过`shader函数`来获取一个对应的颜色值。 + +关于texture,这里在初始化 `payload` 时将`rasterizer` 的成员变量中的 `texture` 传给了 `payload` 。 + + + +最后,展示一下所有不同shader实现的结果~ + +* normal + + image-20220101174322948 + +* Blinn-Phong + + image-20220101174038116 + +* texture + + image-20220101174123089 + +* bump + + image-20220101174204254 + +* displacement + + image-20220101174246060 + + + diff --git a/CG/GAMES101/assignment4.md b/CG/GAMES101/assignment4.md new file mode 100644 index 0000000..2cb0c1e --- /dev/null +++ b/CG/GAMES101/assignment4.md @@ -0,0 +1,48 @@ +--- +title: GAMES101-assignment4-贝塞尔曲线 +date: 2022-02-28 +tags: [GAMES101] +--- +# Assignment4 + +## 1. 运行结果 + +运行结果: + +image-20220228220709563 + +## 2. 实现细节 + +[详细代码](https://github.com/LJHG/GAMES101-assignments) + +这次作业比较简单,写一个递归程序就能实现功能,具体如下所示: + +```cpp +cv::Point2f recursive_bezier(const std::vector &control_points, float t) +{ + // TODO: Implement de Casteljau's algorithm + //当点数量为1时,返回 + if(control_points.size() == 1){ + return control_points[0]; + } + std::vector new_control_points; new_control_points.clear(); + cv::Point2f prev_point = control_points[0]; + for(int i=1;i &control_points, cv::Mat &window) +{ + // TODO: Iterate through all t = 0 to t = 1 with small steps, and call de Casteljau's + // recursive Bezier algorithm. + for(float t = 0; t < 1; t += 0.001){ + cv::Point2f point = recursive_bezier(control_points, t); + window.at(point.y, point.x)[1] = 255; + } +} +``` + + diff --git a/CG/GAMES101/assignment5.md b/CG/GAMES101/assignment5.md new file mode 100644 index 0000000..fe15a89 --- /dev/null +++ b/CG/GAMES101/assignment5.md @@ -0,0 +1,149 @@ +--- +title: GAMES101-assignment5-光线追踪 +date: 2022-03-06 +tags: [GAMES101] +--- +# Assignment5 + +光线追踪 + +## 1. 运行结果 + +运行结果: + +img + +## 2. 实现细节 + +[详细代码](https://github.com/LJHG/GAMES101-assignments) + +这次作业主要是为了实现光线追踪里的与三角形相交,主要为了实现了两个功能: + +1. 实现遍历每个像素,从每个像素发射出光线。 +2. 实现对于光线与三角形相交的判断。 + + + +## 2.1 遍历像素发射光线 + +这里要实现的目标就是遍历每一个像素,对于每一个像素计算出一个 `dir` 向量,来表示一个当前像素的方向向量。 + +具体的代码如下所示: + +```cpp +float scale = std::tan(deg2rad(scene.fov * 0.5f)); +float imageAspectRatio = scene.width / (float)scene.height; //宽高比 + +// Use this variable as the eye position to start your rays. +Vector3f eye_pos(0); +int m = 0; +for (int j = 0; j < scene.height; ++j) +{ + for (int i = 0; i < scene.width; ++i) + { + // generate primary ray direction + // TODO: Find the x and y positions of the current pixel to get the direction + // vector that passes through it. + // Also, don't forget to multiply both of them with the variable *scale*, and + // x (horizontal) variable with the *imageAspectRatio* + float x = (((i + 0.5)/scene.width)*2 - 1) * scale * imageAspectRatio; + float y = (1 - ((j + 0.5)/scene.height)*2) * scale; + Vector3f dir = Vector3f(x, y, -1); // Don't forget to normalize this direction! + framebuffer[m++] = castRay(eye_pos, dir, scene, 0); + } + UpdateProgress(j / (float)scene.height); +} +``` + +这里需要注意一下**坐标的空间变换**(详细的描述可以见[计算机图形学学习笔记——Whitted-Style Ray Tracing(GAMES101作业5讲解)_dong89801033的博客-CSDN博客](https://blog.csdn.net/dong89801033/article/details/114834898?ops_request_misc=%7B%22request%5Fid%22%3A%22162216944616780357298394%22%2C%22scm%22%3A%2220140713.130102334.pc%5Fall.%22%7D&request_id=162216944616780357298394&biz_id=0&utm_medium=distribute.pc_search_result.none-task-blog-2~all~first_rank_v2~rank_v29-2-114834898.pc_search_result_cache&utm_term=games101作业5&spm=1018.2226.3001.4187)),这里一共经历了四次变换,分别是: + +1. Raster space -> NDC space,把坐标归一化到(0,1) + + $$ + x = (i+0.5)/width + $$ + + $$ + y = (j+0.5)/height + $$ + + 其中,加的0.5是为了保持在像素中心。 + +2. NDC space -> Screen space, 把坐标映射到(-1,1) + + $$ + x = 2*x - 1 + $$ + + $$ + y = 1 - 2*y + $$ + + 这里需要注意一下,y的计算是 1-2*y,原因和坐标的正负有关,因为在原本的NDC space中,y的坐标是向下 为正,然而在Screen Space里,y的坐标是向上为正,所以需要取个符号,不然图像会颠倒。 + +3. 图像缩放 + + 因为x和y在进行归一化时,是各自被归一化到了(0,1),因此比例会变得不一样。为了保持图像比例正确,需要为横向乘上一个宽高比。 + + $$ + x = x * width / height + $$ + +4. 关于可视角度 + + 一般来说,摄像机(eye pos)和成像平面的距离为1,所以我们可以通过控制一个可视角度α来决定能够看到的东西的多少。例如,当α为90°时,tan(2/α) = 1, 也就是说可视角度为[-1,1]。 + + 所以,这里我们需要为每一个轴来乘上一个scale。 + + $$ + x = x*tan(\alpha/2) + $$ + + $$ + y = x*tan(\alpha/2) + $$ + + + +最后得到的坐标变化就是: + +```cpp +float x = (((i + 0.5)/scene.width)*2 - 1) * scale * imageAspectRatio; +float y = (1 - ((j + 0.5)/scene.height)*2) * scale; +``` + + + +## 2.2 光线与三角形相交 + +这一块的话直接套Möller-Trumbore算法的公式即可: + +image-20220306133045747 + +```cpp +bool rayTriangleIntersect(const Vector3f& v0, const Vector3f& v1, const Vector3f& v2, const Vector3f& orig, const Vector3f& dir, float& tnear, float& u, float& v) +{ + // TODO: Implement this function that tests whether the triangle + // that's specified bt v0, v1 and v2 intersects with the ray (whose + // origin is *orig* and direction is *dir*) + // Also don't forget to update tnear, u and v. + Vector3f e1 = v1 - v0; + Vector3f e2 = v2 - v0; + Vector3f s = orig - v0; + Vector3f s1 = crossProduct(dir,e2); + Vector3f s2 = crossProduct(s,e1); + + tnear = dotProduct(s2,e2)/dotProduct(s1,e1); + u = dotProduct(s1,s)/dotProduct(s1,e1); + v = dotProduct(s2,dir)/dotProduct(s1,e1); + if(tnear > 0 && u >=0 && v >=0 && (1-u-v) >=0){ + return true; + } + return false; +} +``` + +由于外边在循环里对于每一束光线,要和每一个三角形求交,所以,这里把u,v 以及tnear传出去是为了让外边记录最小的tnear以及u,v来记录和哪一个三角形以及三角形的哪一个地方求交了。 + + + diff --git a/CG/GAMES101/assignment6.md b/CG/GAMES101/assignment6.md new file mode 100644 index 0000000..a590340 --- /dev/null +++ b/CG/GAMES101/assignment6.md @@ -0,0 +1,117 @@ +--- +title: GAMES101-assignment6-光线追踪的加速 +date: 2022-03-07 +tags: [GAMES101] +--- +# Assignment6 + +光线追踪之BVH加速 + +## 1. 运行结果 + +运行结果: + +image-20220306213926220 + +## 2. 实现细节 + +[详细代码](https://github.com/LJHG/GAMES101-assignments) + +这次作业使用了BVH来对光线追踪实现了加速,其中有两个地方比较重要: + +1. 如何实现光线和**A**xis-**A**ligned **B**ounding **B**ox(AABB)求交。 +2. 如何实现光线和BVH的一个node求交。 + + + +## 2.1 光线与AABB求交 + +这个的具体做法在GAMES101课程里已经讲的很清楚了: + +使用六个平面来定义出一个box: + +image-20220307092826919 + +然后计算光线与这六个平面相交的时间,因为是axis-aligned的,所以可以采用正视图/侧视图/俯视图的方式来简易计算求交的时间。 + +image-20220307092721269 + +对于每一个axis,可以计算出两个平面的求交时间分别为t_min和t_max,于是最后我们可以得到三组t_min, t_max。 + +最后取一个交集,就能够计算出光线进入这个box的t_min, t_max: + +image-20220307093155285 + +具体的代码如下所示: + +```python +inline bool Bounds3::IntersectP(const Ray& ray, const Vector3f& invDir, + const std::array& dirIsNeg) const +{ + // invDir: ray direction(x,y,z), invDir=(1.0/x,1.0/y,1.0/z), use this because Multiply is faster that Division + // dirIsNeg: ray direction(x,y,z), dirIsNeg=[int(x>0),int(y>0),int(z>0)], use this to simplify your logic + // TODO test if ray bound intersects + + /** + * 包围盒的元素:Vector3f pMin, pMax; + * 由 location = origin + time * direction --> time = (location - origin)/direction + */ + float t1,t2 = 0; + //这里不能在最后直接乘一个+-1来改变time的正负,不然要出事 + t1 = (pMin.x - ray.origin.x) * invDir.x; + t2 = (pMax.x - ray.origin.x) * invDir.x; + float x_tmin = (dirIsNeg[0] > 0) ? t1 : t2; + float x_tmax = (dirIsNeg[0] > 0) ? t2 : t1; + t1 = (pMin.y - ray.origin.y) * invDir.y; + t2 = (pMax.y - ray.origin.y) * invDir.y; + float y_tmin = (dirIsNeg[1] > 0) ? t1 : t2; + float y_tmax = (dirIsNeg[1] > 0) ? t2 : t1; + t1 = (pMin.z - ray.origin.z) * invDir.z; + t2 = (pMax.z - ray.origin.z) * invDir.z; + float z_tmin = (dirIsNeg[2] > 0) ? t1 : t2; + float z_tmax = (dirIsNeg[2] > 0) ? t2 : t1; + + float t_enter = std::max(x_tmin,std::max(y_tmin,z_tmin)); + float t_exit = std::min(x_tmax,std::min(y_tmax,z_tmax)); + + if(t_enter < t_exit && t_exit >= 0){ + return true; + } + return false; +} +``` + +注意这里是通过判断光线的方向是否为正/负来判断为每一个轴的t_min和t_max来赋什么值,这里不能直接给t直接乘一个正负号,因为后面要判断t_exit >=0 ,如果正负号改变了,那么后面应该也要改,就比较复杂了。 + + + +## 2.2 光线与BVH的一个node求交 + +这就是一个简单的递归: + +1. 首先让光线与node对应的bounding box求交,如果没有交点,那么返回空。 +2. 如果光线与node的bounding box 有交点而且node是叶子节点,那么让光线与叶子节点里的物体求交,并且返回 intersection。(根据这个项目其他部分的代码可以看出,每一个叶子只有一个物体,不会包含多个物体) +3. 如果光线和node的bounding box 有交点而且node不是叶子节点,那么与左右儿子分别相交获取交点,并且返回比较近的那个。 + +```cpp +Intersection BVHAccel::getIntersection(BVHBuildNode* node, const Ray& ray) const +{ + Intersection _intersection; + // TODO Traverse the BVH to find intersection + Vector3f invDir = ray.direction_inv; + std::array dirIsNeg = {int(ray.direction.x > 0), int(ray.direction.y > 0), int(ray.direction.z) > 0}; + + if(!node->bounds.IntersectP(ray,invDir,dirIsNeg) || node == nullptr){ + return _intersection; + } + if(node->left == nullptr && node->right == nullptr){ + //如果是叶节点,与叶节点的objects相交,并且返回intersection + return node->object->getIntersection(ray); + } + Intersection hit1 = getIntersection(node->left,ray); + Intersection hit2 = getIntersection(node->right,ray); + + return (hit1.distance < hit2.distance)?hit1:hit2; +} +``` + diff --git a/CG/GAMES101/assignment7.md b/CG/GAMES101/assignment7.md new file mode 100644 index 0000000..7837008 --- /dev/null +++ b/CG/GAMES101/assignment7.md @@ -0,0 +1,166 @@ +--- +title: GAMES101-assignment7-路径追踪 +date: 2022-03-13 +tags: [GAMES101] +--- +# Assignment7 + +光线追踪之路径追踪 + +## 1. 运行结果 + +当把spp设置为32时,结果如图所示: + +image-20220313154640560 + + + +## 2. 实现细节 + +[详细代码](https://github.com/LJHG/GAMES101-assignments) + +这次作业主要实现了路径追踪,其中涉及到的坑点比较多,主要是以下几个: + +1. 相较于assignment6有几个地方需要修改一下,不然会出现**很错误的结果**。 +2. 直接光照的实现 +3. 间接光照的实现 +4. 还可以拓展的地方 + + + +## 2.1 在assignment6的基础上修改(注意) + +首先,需要对 `BVH.cpp` 的 `BVHAccel::getIntersection` 函数进行修改: + +```cpp +std::array dirIsNeg = {int(ray.direction.x > 0), int(ray.direction.y > 0), int(ray.direction.z) > 0}; +``` + +这里需要将 > 修改为 >= ,不然的话光线和物体求交全是空,结果就是全是黑的。。。 + +```cpp +std::array dirIsNeg = {int(ray.direction.x >= 0), int(ray.direction.y >= 0), int(ray.direction.z) >= 0}; +``` + + + +然后,需要对 `Bounds3.hpp` 的 `Bounds3::IntersectP` 函数进行修改: + +```cpp +if(t_enter < t_exit && t_exit >= 0){ + return true; + } +``` + +这里在判断光线与包围盒相交时要把 < 改为 <= 才对(表示非常不理解),不然的话场景的一大半都会是黑色的。。。 + +```cpp +if(t_enter <= t_exit && t_exit >= 0){ + return true; + } +``` + + + +## 2.2 直接光照的实现 + +直接光照基本就是参照着课程来实现的,说几个需要注意的地方吧。 + +1. 第一个,由于直接光照是直接对光源来进行采样的,所以渲染方程是需要修改的,修改后的代码如下所示: + + ```cpp + float distance = dotProduct(hit_obj_to_light, hit_obj_to_light); + float cosine1 = dotProduct(hit_obj.normal, hit_obj_to_light_dir); + float cosine2 = dotProduct(light_intersection.normal, -hit_obj_to_light_dir); + if(light_pdf > pdf_epsilon) + Lo_direct = light_intersection.emit * f_r * cosine1*cosine2 / light_pdf / distance; + ``` + +2. 作业中存在着大量的向量,并且很多时候不同的向量就是方向取了一个反(对应到代码里就是取了一个负号),我这里的作业里面的一些向量的方向是我自己认为是对的,不保证正确(但是经过测试好像大多数的向量的方向正负并没有什么影响,很迷)。 + +3. 在直接光照的实现中,从物体向光源打出的线需要判断中间是否有物体遮挡,这里不能直接用 `Intersection.happened` 来判断是否有交点,因为要保证交点在物体和光线之间,所以代码应该这样写: + + ```cpp + if(t.distance - hit_obj_to_light.norm() > -epsilon){ + //光线没有被遮挡 + xxxx + } + ``` + + + +## 2.3 间接光照的实现 + +同样,间接光照也指出几个需要注意的地方: + +1. 刚刚直接光照是对光源进行了采样,这里间接光照也需要再进行一次采样,不过就是普通的采样。 + + ```cpp + Vector3f obj1_to_obj2_dir = hit_obj.m->sample(-wo, hit_obj.normal).normalized(); + float obj2_pdf = hit_obj.m->pdf(-wo,obj1_to_obj2_dir,hit_obj.normal); + ``` + +2. 从物体向采样的方向射出一条光线,打中一个物体,这里还需要判断打中的物体是否会发光,不然会出现很多的白色噪点: + + ```cpp + if(t_2.happened && !t_2.m->hasEmission()){ //这里需要确定交点处不会发光,来确保这是间接光照而不是直接射到了光源 + //光线和物体相交 + Vector3f f_r = hit_obj.m->eval(-obj1_to_obj2_dir, wo, hit_obj.normal); + float cosine = dotProduct(hit_obj.normal, obj1_to_obj2_dir); + if(obj2_pdf > pdf_epsilon) + Lo_indirect = shade(t_2, -obj1_to_obj2_dir) * f_r * cosine / obj2_pdf / RussianRoulette; + } + ``` + +3. 最后补充一下,需要单独判断物体是否会发光,如果会发光直接返回,不然的话结果里是没有灯的: + + ```cpp + //如果物体会发光,直接返回光照 + if (hit_obj.m->hasEmission()) + { + return hit_obj.m->getEmission(); + } + ``` + + + +## 2.4 拓展 + +这个作业还有一些可以补充的地方,比如说多线程加速,比如说看到论坛里有人说随机数的生成很大的影响了性能,不过这些都先留着好了。。。 + + + +# 3. 结果展示 + +最后展示一下将spp分别设置为1,16,32的结果: + +* spp = 1 + + image-20220313154827341 + +* spp = 16 + + image-20220313155000558 + +* spp = 32 + + image-20220313155028295 + + + + + + + + + + + + + + + + + + + diff --git a/CG/GAMES202/README.md b/CG/GAMES202/README.md new file mode 100644 index 0000000..3e7a4d3 --- /dev/null +++ b/CG/GAMES202/README.md @@ -0,0 +1,3 @@ +这里会整理GAMES202的各个assignments的代码以及相关的知识 + +* [Assignment1](assignment1.md) diff --git a/CG/GAMES202/assignment1.md b/CG/GAMES202/assignment1.md new file mode 100644 index 0000000..8164000 --- /dev/null +++ b/CG/GAMES202/assignment1.md @@ -0,0 +1,185 @@ +--- +title: GAMES202-Assinment1-实时阴影 +date: 2022-04-09 +tags: [GAMES202] +--- +> 日期:2022.4.11 + +# Assignment1-实时阴影 + +开始看games202了,第一次的作业是实时阴影,参考了一些GitHub的代码然后自己实现了一下,虽然感觉还有不少的问题... + +## 运行结果 + +这里贴一个pcss的运行结果,感觉看起来不太对。。。 + +image-20220411200837460 + +## 实现细节 + +### 0. 坐标转换 + +首先,在写shader之前,需要缕清楚坐标的变换,比如说这里,在真正开始进行深度测试之前,你需要得到shadowmap上真正对应的坐标,也就是**从光源的角度看过去的一个转化为了ndc的坐标**。 + +那么就需要做两步: + +1. 坐标变换,把世界坐标系下(?)转换到从光源角度看过去的坐标。 +2. ndc变化 + +具体的实现就是: + +```glsl +// 坐标变换,在vertex shader里 +vPositionFromLight = uLightMVP * vec4(aVertexPosition, 1.0); + +// ndc变换,在fragment shader里 +vec3 shadowCoord = (vPositionFromLight.xyz / vPositionFromLight.w) * 0.5 + 0.5; //为啥是 *0.5 + 0.5 ? ... +``` + +这里的`*0.5 + 0.5` 我比较迷惑(我抄的github),根据我的理解,前面的`(vPositionFromLight.xyz / vPositionFromLight.w)`是做一个归一化,那么后面的`*0.5 + 0.5`应该就是一个平移的操作,但是我觉得具体乘多少加多少要根据数据的范围来定。 + +> 我猜测进行了`(vPositionFromLight.xyz / vPositionFromLight.w)`的数据范围在 [-1,1],这样的话[-1\*0.5 + 0.5, 1\*0.5+0.5] 刚好就被shift到了[0,1]。。。 + + + +### 1. 硬阴影 + +既然得到了以光源为原点的坐标系的坐标,那么做硬阴影就很简单了,代码如下: + +```glsl +// 直接使用shadowMap产生硬阴影 +float useShadowMap(sampler2D shadowMap, vec4 shadowCoord){ + // 使用shadowCoord在shadowMap采样获得深度值,然后与实际的z进行比较 + vec4 depthPack = texture2D(shadowMap, shadowCoord.xy); //因为从深度图中查询到的实际上是一个颜色 + float depth = unpack(depthPack); //所以在这里需要把深度图查询到的颜色转换为一个浮点数 + if(depth < shadowCoord.z){ + return 0.0; + } + return 1.0; +} +``` + +思路就是在shadowMap里去查那个点的深度,然后与实际的深度比,如果实际的深度更大,那么就没挡住,否则就挡住了。 + +需要稍微注意的是`unpack`操作,因为深度图实际上是一张图片,所以需要将那个点的rgba像素转换为浮点数后才能进行比较。 + +使用硬阴影的效果图如下: + +image-20220411200710483 + +### 2. PCF + +PCF思想也很简单,就是在自己周围做一个采样,然后看有多少物体被挡住(1),有多少没被挡住(0),然后做一个平均来表示自身的visibility。 + +作业给出的框架是使用了均匀采样和泊松采样,比如说像这样: + +```glsl +void poissonDiskSamples( const in vec2 randomSeed ) { + + float ANGLE_STEP = PI2 * float( NUM_RINGS ) / float( NUM_SAMPLES ); + float INV_NUM_SAMPLES = 1.0 / float( NUM_SAMPLES ); + + float angle = rand_2to1( randomSeed ) * PI2; + float radius = INV_NUM_SAMPLES; + float radiusStep = radius; + + for( int i = 0; i < NUM_SAMPLES; i ++ ) { + poissonDisk[i] = vec2( cos( angle ), sin( angle ) ) * pow( radius, 0.75 ); + radius += radiusStep; + angle += ANGLE_STEP; + } +} +``` + +我们只需要给它穿入一个seed,然后就可以通过访问`poissonDisk[i]`去得到一个采样后的值,不过我有点没太搞懂这个得到的采样后的值的一个大致范围,我只知道它是均匀的,所以后面我就开始瞎搞了。 + +PCF的代码大致如下: + +```glsl +float PCF(sampler2D shadowMap, vec4 coords, float filterSize) { + uniformDiskSamples(coords.xy); + // poissonDiskSamples(coords.xy); + float ans = 0.0; + // 其实这个值应该是有一个具体含义,然后计算出来的,不过不太清楚poissonDisk里的值的一个具体范围,所以就随便设了。 + for(int i=0;i coords.z - 0.01) ? 1.0:0.0; // 必须加一个 -0.01来做一个边界的过滤,不然整个屏幕都是糊的 + } + return ans/float(NUM_SAMPLES); +} +``` + +可以看到,这里我乘上了一个`filterSize`,用来表示一个采样的范围,这个应该是可以通过一些推导得出来一个比较合理的值的,不过我直接就把它设成0.01了,效果还不错。。大概长这样: + +image-20220411200749844 + + + +### 3. PCSS + +PCSS就是估计一个平均的blocker depth后去求一个半影长度,后面就和PCF一模一样了。 + +首先是求平均的blocker depth,大概长这样: + +```glsl +float findBlocker( sampler2D shadowMap, vec2 uv, float zReceiver ) { + uniformDiskSamples(uv); //采样 + float allDepth = 0.0; + int blockerCnt = 0; + for(int i=0;i zReceiver - 0.01) ? 1:0; + } + return float(allDepth) / float(blockerCnt); +} +``` + +在求平均的blocker depth时,采样范围时多少又是一个问题。。。我记得在课上好像讲了一下使用光源和blocker以及shading point的一个距离来取一个采样范围,大致就是blocker离光源近,采样范围就大,blocker离光源远,采样范围就小。 + +图好像是这样的: + +image-20220411200135154 + +不过我没有用这种方法(抱歉),我又直接设了一个 `BLOCKER_SIZE = 0.01`,可能就是这个原因让后面的结果看起来很怪。 + +求出了 blocker depth 后,再算一个半影长度(把半影长度当成一个filter size),最后就可以直接套PCF了,代码大致如下: + +```glsl +float PCSS(sampler2D shadowMap, vec4 coords){ + + // STEP 1: avgblocker depth + float averageDepth = findBlocker(shadowMap, coords.xy, coords.z); + + // STEP 2: penumbra size + // 使用相似三角形计算半影 + float dReceiver = coords.z; + float dBlocker = averageDepth; + float penumbraSize = float(L_WIDTH) * (dReceiver - dBlocker)/ dBlocker * 0.01; + // STEP 3: filtering + return PCF(shadowMap, coords, penumbraSize); +} +``` + +关于使用相似三角形计算半影的图长这样: + +image-20220411195928455 + +至于为什么我在算半影时后面又乘了一个0.01,是因为不乘0.01范围又大了。。就把这个0.01当成filter_size吧。 + +最后PCSS的结果如下: + +> 这肯定不对,可能还是算平均blocker depth那里有问题🤔,这里的图是 BLOCKER_DPETH 为0.01的结果,随着把BLOCKER_DPETH变大,实的影子的长度会越来越短,0.05~0.08时效果还行,但是都是野路子。。。 + +image-20220411200837460 + + + +## 总结 + +感觉这次做的太糙了,但是也不想做了,暂时先这样,写shader太坐牢了,连print都不能打,调试地狱。所以我还是不知道到底采样后的数组里的数据范围是多少,这也是我在使用万能的0.01的原因。 + + + diff --git a/CG/OpenGL/README.md b/CG/OpenGL/README.md new file mode 100644 index 0000000..99406c1 --- /dev/null +++ b/CG/OpenGL/README.md @@ -0,0 +1,2 @@ +记录在opengl下踩的坑 +* [如何使用assimp读取没有贴图的obj模型](assimp_material.md) \ No newline at end of file diff --git a/CG/OpenGL/assimp_material.md b/CG/OpenGL/assimp_material.md new file mode 100644 index 0000000..2528973 --- /dev/null +++ b/CG/OpenGL/assimp_material.md @@ -0,0 +1,112 @@ +--- +title: 如何使用assimp读取obj文件的材质 +date: 2022-03-25 +tags: [OpenGL] +--- +> 日期: 2022年3月25日 + + + +# 如何使用assimp读取obj文件的材质 + +之前在看LearnOpenGL的[这一节](https://learnopengl-cn.github.io/03 Model Loading/03 Model/),主要讲的就是使用assimp来加载模型,但是给的参考代码只针对了模型有材质贴图的情况,也就是通过去获取diffuse贴图以及specular贴图的方式获得模型的材质。 + +比如这是官方给的示例的文件,其中包含了很多的材质贴图: + +image-20220325135545963 + +在相关的mtl文件里也链接到了具体的贴图文件: + +image-20220325135723346 + +但是我们在使用blender等很多软件建模时是没有材质贴图的,一般来说只有一个obj文件和一个mtl文件: + +image-20220325135810722 + +image-20220325135828015 + +所以只要想办法获取到mtl文件里的Ka, Kd, Ks这几项,然后想办法把它传进shader渲染就行了。 + + + +首先,使用assimp来获取这几个值的核心代码在这里: + +```cpp +Material loadMaterialWithoutTextures(aiMaterial *mat){ + Material result; + aiColor3D color; + //读取mtl文件顶点数据 + mat->Get(AI_MATKEY_COLOR_AMBIENT, color); + result.ambient = glm::vec3(color.r, color.g, color.b); + mat->Get(AI_MATKEY_COLOR_DIFFUSE, color); + result.diffuse = glm::vec3(color.r, color.g, color.b); + mat->Get(AI_MATKEY_COLOR_SPECULAR, color); + result.specular = glm::vec3(color.r, color.g, color.b); + return result; + } +``` + +需要进行一些判断来这个mesh是否是拥有材质贴图的mesh,如果没有材质贴图,就去读取Ka, Kd, Ks,否则就直接去读取贴图(原本的实现方法)。 + +```cpp +if(material->GetTextureCount(aiTextureType_DIFFUSE) == 0 && material->GetTextureCount(aiTextureType_SPECULAR) == 0){ + //没有贴图,之前读取模型的material的值 + Material colorMaterial = loadMaterialWithoutTextures(material); + return Mesh(vertices, indices, colorMaterial); + } + else{ + vector diffuseMaps = loadMaterialTextures(material, + aiTextureType_DIFFUSE, "texture_diffuse"); + textures.insert(textures.end(), diffuseMaps.begin(), diffuseMaps.end()); + vector specularMaps = loadMaterialTextures(material, + aiTextureType_SPECULAR, "texture_specular"); + textures.insert(textures.end(), specularMaps.begin(), specularMaps.end()); + return Mesh(vertices, indices, textures); + } +``` + + + +这里仅仅列举了一部分需要修改的代码,实际上还需要修改很多地方,比如说Mesh相关的构造函数,Draw函数等等,比如说,我这里还修改了相关的Shader,通过在shader代码里创建两个不同的材质结构体,以及一个bool量,来实现决定渲染时,应该采用哪一种方式: + +```cpp +//带有贴图的Material +struct Material { + sampler2D diffuse; + sampler2D specular; + float shininess; +}; +uniform Material material; + +//不带有贴图的Material +struct MaterialWithoutTexture { + vec3 ambient; + vec3 diffuse; + vec3 specular; + float shininess; +}; +uniform MaterialWithoutTexture materialWithoutTexture; + +uniform bool hasTexture; +``` + +```cpp +vec3 ambient; + vec3 diffuse; + vec3 specular; + if(hasTexture){ + ambient = light.ambient * vec3(texture(material.diffuse, TexCoords)); + diffuse = light.diffuse * diff * vec3(texture(material.diffuse, TexCoords)); + specular = light.specular * spec * vec3(texture(material.specular, TexCoords)); + }else{ + ambient = light.ambient * materialWithoutTexture.diffuse ; //虽然material里有ambient,但是不使用ambient + diffuse = light.diffuse * diff * materialWithoutTexture.diffuse; + specular = light.specular * spec * materialWithoutTexture.specular; + } +``` + + + +最后就可以直接对没有贴图的obj来进行渲染了~ + +image-20220325151732301 diff --git a/CG/README.md b/CG/README.md new file mode 100644 index 0000000..ce08f91 --- /dev/null +++ b/CG/README.md @@ -0,0 +1,2 @@ +这里会整理一些关于计算机图形学的知识。 +* [GAMES101 assignments](GAMES101/README.md) \ No newline at end of file diff --git a/CG/SoftRasterizer/3d.md b/CG/SoftRasterizer/3d.md new file mode 100644 index 0000000..b50cca5 --- /dev/null +++ b/CG/SoftRasterizer/3d.md @@ -0,0 +1,232 @@ +--- +title: 软光栅器MicroRenderer(二) 进入三维 +date: 2022-05-27 23:38:00 +tags: [MicroRenderer] +--- +# 软光栅器MicroRenderer(二) 进入三维 +代码见:https://github.com/LJHG/MicroRenderer + +本次主要讨论的内容是 mvp变换,三角形的裁切以及相机控制。 + +## 1. mvp变换 + +具体的mvp矩阵推导这里就不细讲了,关于mvp变换我基本上是参照了games101课程里的内容,可以看一下这里:[Assignment1 · GitBook (ljhblog.top)](https://www.ljhblog.top/CG/GAMES101/assignment1.html) 。其中我的view matrix 和 projection matrix 是完全按照课程上的描述实现的,不过最后的我的视图矩阵和透视矩阵好像实现出来好像都有点问题。。。所以目前我的这两个矩阵都是抄的[别人的](https://yangwc.com/2019/05/01/SoftRenderer-Math/),暂时还没啥问题。 + +一旦你会计算mvp矩阵了,剩下的事情就很简单了,只需要去继承Shader然后实现一个新的类,并且把vertex shader重新实现为三维版本就可以了,在三维空间的shader里,vertex shader的职责就是把顶点进行mvp变换。 + +```cpp +// 3d shader, implement mvp transformation for vertex shader, fragment shader just pass color +class ThreeDShader: public Shader{ +public: + ThreeDShader(glm::mat4 _modelMatrix, glm::mat4 _viewMatrix, glm::mat4 _projectionMatrix); + VertexOutData vertexShader(VertexData &v) override; + glm::vec4 fragmentShader(VertexOutData &v) override; +}; +``` + +```cpp +VertexOutData ThreeDShader::vertexShader(VertexData &v) { + VertexOutData vod; + vod.position = glm::vec4(v.position,1.0f); + vod.position = projectionMatrix * viewMatrix * modelMatrix * vod.position; //mvp transformation + vod.color= v.color; + vod.normal = v.normal; + return vod; +} + +glm::vec4 ThreeDShader::fragmentShader(VertexOutData &v) { + glm::vec4 result; + result = v.color; + return result; +} +``` + +因为不涉及具体的着色,所以这里我们的 fragment shader仍然传递颜色即可。 + + + +## 2. 三角形的裁切 + +> 究极野路子,谨慎参考 + +关于三角形的裁切这个东西,一般来说应该是放在光栅化之前来做:[渲染管线中,背面剔除和齐次空间的裁剪哪个先进行?为什么? - 知乎 (zhihu.com)](https://www.zhihu.com/question/469259481/answer/1973677886)。裁切有很多很复杂的算法。一般来说,大家会把三角形裁切和背面剔除放在一起来讨论。 + +不过我这里就比较拉了,为了实现最简单的方法,我把裁切放在了光栅化的时候来做。具体的原理就是,通过对经过了视口变换的坐标进行判断,也就是说如果x坐标不在 `[0,width]` 内,y坐标不在 `[0,height]`内,z坐标不在`[-1,1]`内,那么就不去画它。。。 + +具体的代码就是这样: + +```cpp +void ShadingPipeline::fillRasterization(VertexOutData &v1, VertexOutData &v2, VertexOutData &v3) { + // I don't know name for this algorithm + // ref: GAMES101 assignment2: https://www.ljhblog.top/CG/GAMES101/assignment2.html + float x1 = v1.position[0]; float y1 = v1.position[1]; float z1 = v1.position[2]; + float x2 = v2.position[0]; float y2 = v2.position[1]; float z2 = v2.position[2]; + float x3 = v3.position[0]; float y3 = v3.position[1]; float z3 = v3.position[2]; + + /** + * 关于clip是否应该在光栅化时做的问题 + * 一般来说,clip好像应该在做光栅化前来做,但是如果直接在光栅化时做会十分方便... + * 这里我在光栅化时作了两个clip: + * 1. 一个判断z,如果z不在[-1,1]内,那么就说明不再视线范围内,直接不画 + * 2. 另一个对left right bottom top 做一个裁剪,不去管屏幕空间外的像素了 + */ + // clip1 + if(z1<-1.0f || z1 >1.0f || z2<-1.0f || z2 > 1.0f || z3<-1.0f || z3>1.0f){ + //out of depth range, clip + return; + } + + // get bounding box of the triangle + int left = static_cast(std::min(std::min(x1,x2),x3)); + int bottom = static_cast(std::min(std::min(y1,y2),y3)); + int right = static_cast(std::max(std::max(x1,x2),x3)) + 1; //ceil + int top = static_cast(std::max(std::max(y1,y2),y3)) + 1; //ceil + + // clip2 + left = std::max(left,0); + bottom = std::max(bottom,0); + right = std::min(right,width-1); // need to minus 1, or will draw at the other edge... + top = std::min(top,height-1); // need to minus 1, or will draw at the other edge... + + for(int x = left;x<=right;x++){ + for(int y=bottom;y<=top;y++){ + if(MathUtils::insideTriangle(static_cast(x),static_cast(y),x1,y1,x2,y2,x3,y3)){ + VertexOutData tmp = MathUtils::barycentricLerp(x,y,v1,v2,v3); + int index = y*width+x; + if(tmp.position[2] < zBuffer[index]){ + zBuffer[index] = tmp.position[2]; + //fragment shader + glm::vec4 color = shader->fragmentShader(tmp); + int colorIndex = index*3; //multiply channel + image[colorIndex +0] = static_cast(color[0]); + image[colorIndex +1] = static_cast(color[1]); + image[colorIndex +2] = static_cast(color[2]); + } + } + + } + } +} +``` + +同时要注意的是,裁切非常的重要,不仅仅是性能的问题,还有显示的问题。如果你不对z方向进行裁切,会发生非常诡异的结果(明明不在视线内的东西还会画出来,所以其实最开始我是没发现要裁切的,直到我发现了这个诡异的bug)。 + + + +## 3. 相机控制 + +关于相机的控制,主要的实现流程就是控制相机,得到相机变换后的 位置,看的地方,然后直接算视图矩阵就行了。 + +这里要稍微注意一下,那就是在计算视图矩阵时: + +```cpp +glm::mat4 MathUtils::calViewMatrix(glm::vec3 cameraPos, glm::vec3 target, glm::vec3 worldUp) +``` + +这个 `target` 实际上是一个位置。一般来说,在创建相机时,我们需要一个变量来保存相机的朝向,例如: + +```cpp +class Camera { +public: + Camera(float _width, float _height); + glm::vec3 cameraPos; // camera position + glm::vec3 target; // camera focus position, target = cameraPos + cameraFront + glm::vec3 worldUp; // (0,1,0) + float fov; + float aspectRatio; + float zNear; + float zFar; + float width; + float height; + float cameraMoveSpeed; + float cameraTurnSpeed; + glm::vec3 cameraFront; // camera direction +}; +``` + +由于相机看的地方肯定是在自己的朝向的前面,所以就可以使用 `target = cameraPos + cameraFront` 来计算相机看的地方了。 + +这部分具体的东西就是一些SDL2相关的操作了,贴下代码吧: + +```cpp +if(event.type == SDL_KEYDOWN){ + switch (event.key.keysym.sym) { + case 'w': + { + camera.cameraPos += camera.cameraFront * camera.cameraMoveSpeed; + break; + } + case 'a': + { + camera.cameraPos += glm::normalize(glm::cross(camera.cameraFront, camera.worldUp)) * camera.cameraMoveSpeed; + break; + } + case 's': + { + camera.cameraPos -= camera.cameraFront * camera.cameraMoveSpeed; + break; + } + case 'd': + { + camera.cameraPos -= glm::normalize(glm::cross(camera.cameraFront, camera.worldUp)) * camera.cameraMoveSpeed; + break; + } + case SDLK_UP: + { + camera.cameraFront = glm::normalize(camera.cameraFront + camera.worldUp * camera.cameraTurnSpeed); + break; + } + case SDLK_DOWN: + { + camera.cameraFront = glm::normalize(camera.cameraFront - camera.worldUp * camera.cameraTurnSpeed); + break; + } + case SDLK_LEFT: + { + camera.cameraFront = glm::normalize(camera.cameraFront + glm::cross(camera.cameraFront, camera.worldUp) * camera.cameraTurnSpeed); + break; + } + case SDLK_RIGHT: + { + camera.cameraFront = glm::normalize(camera.cameraFront - glm::cross(camera.cameraFront, camera.worldUp) * camera.cameraTurnSpeed); + break; + } + } + camera.target = camera.cameraPos + camera.cameraFront; +} +if(event.type == SDL_MOUSEBUTTONDOWN){ + mouseClick = true; + lastMouseX = event.button.x; + lastMouseY = event.button.y; + +} +if(event.type == SDL_MOUSEBUTTONUP){ + mouseClick = false; +} +if(mouseClick){ + if(event.type == SDL_MOUSEMOTION){ + glm::vec3 front = getCameraFront(camera,event.button.x,event.button.y); + camera.cameraFront = glm::normalize(camera.cameraFront + front * camera.cameraTurnSpeed * mouseDragSpeed); + camera.target = camera.cameraPos + camera.cameraFront; + } +} +``` + + + +## 4. 效果展示 + +最后就可以得到这样的效果: + +![camera_mvp](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202205272336655.gif) + + + + + + + + + + + diff --git a/CG/SoftRasterizer/README.md b/CG/SoftRasterizer/README.md new file mode 100644 index 0000000..8d2f821 --- /dev/null +++ b/CG/SoftRasterizer/README.md @@ -0,0 +1 @@ +写个软光栅器 \ No newline at end of file diff --git a/CG/SoftRasterizer/pipeline.md b/CG/SoftRasterizer/pipeline.md new file mode 100644 index 0000000..9fd2895 --- /dev/null +++ b/CG/SoftRasterizer/pipeline.md @@ -0,0 +1,195 @@ +--- +title: 软光栅器MicroRenderer(一) 渲染管线与三角形 +date: 2022-05-24 16:52:48 +tags: [MicroRenderer] +--- +# 软光栅器MicroRenderer(一) 渲染管线与三角形 +代码见:https://github.com/LJHG/MicroRenderer + +最近仿照了这位大佬的博客[软渲染器Soft Renderer:3D数学篇 | YangWC's Blog](https://yangwc.com/2019/05/01/SoftRenderer-Math/)来实现一个软光栅器,其中的代码结构基本上参照了它的代码,一边参考一遍自己写把渲染管线实现了一下,在这里稍微记录一下。 + +## 1. 代码结构 + +目前基本的代码结构是这样的: + +image-20220524170700703 + + + +这里大致说一下每一个文件的作用(虽然后续有重构的可能性): + +1. CommonUtils.cpp。用来存放一些通用函数,比如说读取图片之类的。 +2. MathUtils.cpp。用来存放数学相关的函数,比如说计算MVP矩阵,插值等。 +3. Renderer.cpp。渲染相关代码,用来管理渲染的物体以及调用渲染管线渲染等。 +4. Shader.cpp。存放shader代码。 +5. ShadingPipeline.cpp。渲染管线代码,生成渲染结果。 +6. Structure.cpp。定义vetex, frament 以及 mesh的数据类型。 +7. WindowApp.cpp。SDL2相关函数。 + +> 关于Renderer 和ShadingPipeline 以及 Shader 的组织问题: +> +> 关于这三个类的所属关系大概是: +> +> ShadingPipeline 是 Renderer的一个成员,同时 Shader 是 Shadingpipeline的一个成员。 +> +> Renderer主要负责管理场景里的Mesh,同时调用Shadingpipeline去渲染。 +> +> ShadingPipeline主要负责渲染管线的一整套流程,比如说根据vertices indices 构建三角形,比如说光栅化。 +> +> Shader主要负责实现自己的vertexShader以及FragmentShader函数。 + +## 2. 具体细节 + +### 2.1 像素展示 + +基本上可以把渲染的整个过程抽离为 生成像素buffer -> 显示像素buffer。显示像素buffer这个部分就交给了SDL2来做,只需要在每一帧的时候把一个像素数组传递给SDL2就可以展示了,所以我们的渲染器所需要的就是在每一帧生成出这个像素数组就可以了。这里的选取是非常随意的,在GAMES101里使用的opencv来把结果存成图片,[这位大佬的博客](https://yangwc.com/2019/05/01/SoftRenderer-Rasterization/)使用了qt来展示结果,后来他后续把代码重构了,也是使用的SDL2。 + +比如说其实核心就是 `updateCanvas`这个函数: + +```cpp +void WindowApp::updateCanvas(uint8_t* pixels, int width, int height, int channel) { + SDL_LockSurface(windowSurface); + Uint32* surfacePixels = (Uint32*)windowSurface->pixels; //获取当前屏幕的像素指针 + for(int i=0;iformat, + pixels[index * channel + 0], + pixels[index * channel + 1], + pixels[index * channel + 2]); + surfacePixels[index] = _color; + } + } + SDL_UnlockSurface(windowSurface); + SDL_UpdateWindowSurface(window); +} +``` + + 通过每次传入一个pixels矩阵就可以实现画面的更新了。 + +### 2.2 shader设计 + +shader的设计大概是这样的: + +```cpp +class Shader { +public: + virtual VertexOutData vertexShader(VertexData& v) {VertexOutData _v; return _v;}; //meaningless implementation + virtual glm::vec4 fragmentShader(VertexOutData& v) {glm::vec4 r; return r;}; //meaningless implementation + void setModelMatrix(glm::mat4 matrix); + void setViewMatrix(glm::mat4 matrix); + void setProjectionMatrix(glm::mat4 matrix); +private: + glm::mat4 modelMatrix; + glm::mat4 viewMatrix; + glm::mat4 projectionMatrix; +}; + +class SimpleShader: public Shader{ +public: + VertexOutData vertexShader(VertexData &v) override; + glm::vec4 fragmentShader(VertexOutData &v) override; +}; +``` + +如果要实现一个shader,比如说phong shader, 那么就去继承基类Shader,然后实现它的虚函数`vertexShader`以及 `fragmentShader`。一般来说,`vertexShader`里要做的事就是实现mvp变换, `fragmentShader`要做的事就是对于一个像素,给出它的颜色(这里也不完全是,但是对于这里的最简单的SIMPLE_SHADER是这样的)。 + +### 2.3 渲染管线 + +在这里重温一下渲染管线: + +![OpenGL Rendering Pipeline | Download Scientific Diagram](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202205241725833.ppm) + +首先传入vertex data, 经过vertex shader后会得到一个新的数据,里面会包括经过坐标变换后的坐标,再经过了视图变换(我不是很清楚视图变换是否应该写在vertex shader里)过后就可以进行光栅化了。在进行光栅化的过程中,我们需要决定每一个像素应该给什么颜色,所以就需要使用插值来插值出那个像素点的世界坐标,纹理坐标等数据,然后把它传给fragment shader,最后得到颜色。 + +这里要说的是,我这里写的代码fragment shader的输入数据是VertextOut,这个数据类型是vextex shader输出的数据类型。通过对不同VertextOut进行差值,我们可以得到一个新的VertextOut来作为结果给fragment shader,所以好像也没毛病。。。 + +这里展示一下渲染管线的代码: + +```cpp +void ShadingPipeline::shade(const std::vector& _vertices, + const std::vector& _indices, + int rasterizingMode) { + //according to indices, every 3 indices organize as a triangle, len(indices) could be greater than len(_vertices) + + VertexData v1,v2,v3; + VertexOutData v1o,v2o,v3o; + for(int i=0;i<_indices.size()/3;i++){ + v1 = _vertices[_indices[i*3+0]]; + v2 = _vertices[_indices[i*3+1]]; + v3 = _vertices[_indices[i*3+2]]; + //vertex shader + v1o = shader->vertexShader(v1); + v2o = shader->vertexShader(v2); + v3o = shader->vertexShader(v3); + + //view port transformation + v1o.position = viewPortMatrix * v1o.position; + v2o.position = viewPortMatrix * v2o.position; + v3o.position = viewPortMatrix * v3o.position; + + //rasterization + // the triangle will appear upside down because it goes like ➡️ x ⬇️ y, but never mind... + if(rasterizingMode == LINE){ + // BresenHam line drawing algorithm + bresenhamLineRasterization(v1o,v2o); + bresenhamLineRasterization(v1o,v3o); + bresenhamLineRasterization(v3o,v2o); + }else if(rasterizingMode == FILL){ + // bounding box inside triangle fill algorithm -> games101 assignment2 and assignment3 + fillRasterization(v1o,v2o,v3o); + } + + // double buffer + swapBuffer(); + } + } +``` + +因为这里是直接画三角形,我的vertex shader里直接传递数据,同时fragment shader也是直接传递颜色。所以我传的vertex data就是直接的ndc坐标,经过视图变化过后得到像素坐标,然后就可以进行光栅化了。 + +其中,填充的光栅化算法长这样: + +```cpp +void ShadingPipeline::fillRasterization(VertexOutData &v1, VertexOutData &v2, VertexOutData &v3) { + // I don't know name for this algorithm + // ref: GAMES101 assignment2: https://www.ljhblog.top/CG/GAMES101/assignment2.html + float x1 = v1.position[0]; float y1 = v1.position[1]; float z1 = v1.position[2]; + float x2 = v2.position[0]; float y2 = v2.position[1]; float z2 = v2.position[2]; + float x3 = v3.position[0]; float y3 = v3.position[1]; float z3 = v3.position[2]; + + // get bounding box of the triangle + int left = static_cast(std::min(std::min(x1,x2),x3)); + int bottom = static_cast(std::min(std::min(y1,y2),y3)); + int right = static_cast(std::max(std::max(x1,x2),x3)) + 1; //ceil + int top = static_cast(std::max(std::max(y1,y2),y3)) + 1; //ceil + + for(int x = left;x<=right;x++){ + for(int y=bottom;y<=top;y++){ + if(MathUtils::insideTriangle(static_cast(x),static_cast(y),x1,y1,x2,y2,x3,y3)){ + VertexOutData tmp = MathUtils::barycentricLerp(x,y,v1,v2,v3); + int index = y*width+x; + if(tmp.position[2] < zBuffer[index]){ + zBuffer[index] = tmp.position[2]; + //fragment shader + glm::vec4 color = shader->fragmentShader(tmp); + int colorIndex = index*3; //multiply channel + imageSwap[colorIndex +0] = static_cast(color[0]); + imageSwap[colorIndex +1] = static_cast(color[1]); + imageSwap[colorIndex +2] = static_cast(color[2]); + } + } + + } + } +} +``` + +最后得到的结果大概长这样: + +image-20220524174029635 + +image-20220524174119417 + + diff --git a/CG/SoftRasterizer/shading.md b/CG/SoftRasterizer/shading.md new file mode 100644 index 0000000..6f67b7f --- /dev/null +++ b/CG/SoftRasterizer/shading.md @@ -0,0 +1,517 @@ +--- +title: 软光栅器MicroRenderer(三) 着色器(完结) +date: 2022-06-01 12:06:00 +tags: [MicroRenderer] +--- +# 软光栅器MicroRenderer(三) 着色器 +代码见:https://github.com/LJHG/MicroRenderer + +本节会对软渲染器剩下的内容进行说明,主要包括: +1. 模型加载 +2. 贴图uv颜色采样 +3. 更加复杂的shader的设计,本项目中主要实现了 Gouraud shader以及Phong shader。 + + + +## 1. 模型加载 + +首先来看一下Mesh这个类都包含哪些成员变量: + +```cpp +class Mesh{ +public: + Mesh(); + + /** getter and setter for vertices and indices **/ + std::vector& getVertices(); + std::vector& getIndices(); + Material getMaterial(); + std::string getTextureUrl(); + void setVertices(const std::vector& _vertices); + void setIndices(const std::vector& _indices); + void setMaterial(Material _material); + void addVertex(VertexData v); + void addIndex(unsigned int idx); + void setTextureUrl(std::string _textureUrl); + + /** set model matrix **/ + void setModelMatrix(const glm::mat4& m); + glm::mat4 getModelMatrix(); + + + /** quick create mesh functions **/ + void asTriangle(VertexData v1, VertexData v2, VertexData v3); // read 3 vertices and as a triangle + void asCube(glm::vec3 center, float size, Material _material); // use fot generate point light + + +private: + std::vector vertices; + std::vector indices; + Material material; + std::string textureUrl; + glm::mat4 modelMatrix; +}; +``` + +可以看到,首先是最重要的顶点和下标,也就是vertices 和 indices,这存储了顶点的信息以及顶点之间的连接关系。同时mesh还包含了材质mateiral,包含了ka, kd 以及ks,一般来说我们认为一个mesh的材质是固定的。mesh还包含了textureUrl用来保存贴图的路径,以及一个modelMatrix,来处理有的mesh需要移动的情况。 + +同时,VertexData的成员主要如下所示: + +```cpp +struct VertexData{ + glm::vec3 position; + glm::vec3 normal; + glm::vec4 color; //如果是Phong shader, 那么就不使用color这个属性,而是去查询纹理或者直接用材质 + glm::vec2 textureCoord; // 纹理坐标 +}; +``` + +对于color这个属性,这里要加以说明一下。之前我们在使用 THREE_D_SHADER 绘制三角形时,直接给顶点赋了颜色,同时vertex shader 和 fragment shader 也是直接传递颜色,那个时候就是直接使用color这个属性来进行绘制。 + +在使用 Gouraud shader 或者 Phong shader时,顶点的颜色就不会被指定了。所以VertexData这个结构体的color属性就不会被用到了。这两个Shader具体怎么着色会在后续进行说明。 + + + +上面对于一个Mesh的成员变量已经阐述完毕了,所以接下来我们就要在加载obj模型时把这些成员变量给填满就行了。 + +在这里模型加载我使用了[assimp](https://github.com/assimp/assimp)这个来进行实现。只要给定obj的文件路径,就可以从**obj里面解析出若干个mesh**。这里要稍微注意一下,一个obj文件里是可以包括多个mesh的,所以不能把加载模型的这个函数作为Mesh类的一个成员函数。因此,我单独在外边写了一个加载模型的函数,并且返回的是一个Mesh的vector。 + +```cpp +static std::vector loadObjModel(std::string filepath) +{ + std::vector meshes; + Assimp::Importer import; + const aiScene* scene = import.ReadFile(filepath, aiProcess_Triangulate | aiProcess_FlipUVs); + // 异常处理 + if (!scene || scene->mFlags == AI_SCENE_FLAGS_INCOMPLETE || !scene->mRootNode) + { + std::cout << "读取模型出现错误: " << import.GetErrorString() << std::endl; + exit(-1); + } + // 模型文件相对路径 + std::string rootPath = filepath.substr(0, filepath.find_last_of('/')); + + // 循环生成 mesh + for (int i = 0; i < scene->mNumMeshes; i++) + { + // 创建mesh + Mesh mesh; + + // 获取 assimp 的读取到的 aiMesh 对象 + aiMesh* aiMesh = scene->mMeshes[i]; + + VertexData vertexData; + // 我们将数据传递给我们自定义的mesh + for (int j = 0; j < aiMesh->mNumVertices; j++) + { + // 顶点 + glm::vec3 vvv; + vvv.x = aiMesh->mVertices[j].x; + vvv.y = aiMesh->mVertices[j].y; + vvv.z = aiMesh->mVertices[j].z; + vertexData.position = vvv; + + // 法线 + vvv.x = aiMesh->mNormals[j].x; + vvv.y = aiMesh->mNormals[j].y; + vvv.z = aiMesh->mNormals[j].z; + vertexData.normal = glm::normalize(vvv); // 法线记得要归一化 + + // 纹理坐标: 如果存在则加入。assimp 默认可以有多个纹理坐标 我们取第一个(0)即可 + glm::vec2 vv(0, 0); + if (aiMesh->mTextureCoords[0]) + { + vv.x = aiMesh->mTextureCoords[0][j].x; + vv.y = aiMesh->mTextureCoords[0][j].y; + } + vertexData.textureCoord = vv; + vertexData.color = glm::vec4(128,128,128,0); //设置一个默认颜色。。 + mesh.addVertex(vertexData); + } + + // 传递面片索引 + for (int j = 0; j < aiMesh->mNumFaces; j++) + { + aiFace face = aiMesh->mFaces[j]; + for (int k = 0; k < face.mNumIndices; k++) + { + mesh.addIndex(face.mIndices[k]); + } + } + + // 读取材质和贴图 + if(aiMesh->mMaterialIndex >= 0){ + aiMaterial *material = scene->mMaterials[aiMesh->mMaterialIndex]; + // 这里只考虑漫反射贴图了,暂时先不去管镜面反射贴图 + if(material->GetTextureCount(aiTextureType_DIFFUSE) == 0){ + //没有漫反射贴图,直接读取 ka kd ks + std::cout<<"do not has texture, read material"<Get(AI_MATKEY_COLOR_AMBIENT, color); + glm::vec3 ka = glm::vec3(color.r,color.g, color.b); + material->Get(AI_MATKEY_COLOR_DIFFUSE, color); + glm::vec3 kd = glm::vec3(color.r,color.g, color.b); + material->Get(AI_MATKEY_COLOR_SPECULAR, color); + glm::vec3 ks = glm::vec3(color.r,color.g, color.b); + mesh.setMaterial(Material(ka,kd,ks,16.0f)); // set shininess at 16.0f... + } + else{ + //有漫反射贴图,记录 ka ,ks,并且把贴图路径记录下来 + std::cout<<"has texture"<Get(AI_MATKEY_COLOR_AMBIENT, color); + glm::vec3 ka = glm::vec3(color.r,color.g, color.b); + glm::vec3 kd = glm::vec3(0.0f,0.0f, 0.0f); + material->Get(AI_MATKEY_COLOR_SPECULAR, color); + glm::vec3 ks = glm::vec3(color.r,color.g, color.b); + mesh.setMaterial(Material(ka,kd,ks,16.0f)); // set shininess at 16.0f... + // 贴图路径 + aiString _textureUrl; + material->GetTexture(aiTextureType_DIFFUSE,0,&_textureUrl); + mesh.setTextureUrl(rootPath + "/" + _textureUrl.C_Str()); + } + } + + + // add mesh + meshes.push_back(mesh); + } + + return meshes; +} +``` + +这里要稍微说一下关于材质和贴图的处理。有的obj文件是没有贴图的,对于这种情况就直接读取模型的材质ka kd ks即可,模型的材质路径不需要去单独设置了,不去设置会默认设置为none。如果obj文件有贴图,那么就把kd这一项给设置为全0,并且为mesh设置材质的路径,后续读取到贴图的颜色后再把颜色赋值给kd。 + +> 这里之所以这样做是因为我们这里只考虑了漫反射贴图,其实有的时候还会有镜面反射贴图的,不过这里先不管这些了。在只考虑漫反射贴图时,我们可以使用在贴图上颜色去替代材质的kd这一项。 + +至此,模型的加载已经结束。 + +## 2. 贴图uv颜色采样 + +刚刚在读取模型时我们已经在Mesh里面拿到了贴图的文件路径。根据顶点的纹理坐标,我们可以在贴图上去拿到颜色值,然后赋给kd,进行着色。 + +在这里我们使用 [stb_image.h](https://github.com/nothings/stb/blob/master/stb_image.h) 这个头文件来进行图片的读取。 + +```cpp +Image CommonUtils::loadImage(std::string url) { + int width,height,channel; + uint8_t* rgb_image = stbi_load(url.c_str(), &width, &height, &channel, 3); + Image image(rgb_image,width,height,channel); + return image; +} +``` + +这里读取来的图片会被保存在一个 char* 的数组里。图片的保存格式是RGBRGBRGB...左上角对应index 0,右下角对应 index 最大,同时颜色按行存储。比如说第一行的index对应 0,1,2,3,4,第二行对应 5,6,7,8,9。这非常的重要,这影响到后面颜色采样是否对的问题。 + +后续在着色时,我们会再次把mesh的贴图路径传递给shader,并且在shader里把贴图读取到,然后就可以开始进行颜色的采样了。 + +我这里时直接把颜色采样写在了shader里,没有单独再去封装函数了,大概长这样: + +```cpp +if(textureUrl != "none"){ + int xIdx = v.textureCoord[0] * texture.width; + int yIdx = (1-v.textureCoord[1]) * texture.height; + int idx = (yIdx * texture.width + xIdx)*3; + kd = glm::vec3(texture.pixels[idx+0]/255.0f, + texture.pixels[idx+1]/255.0f, + texture.pixels[idx+2]/255.0f); +} +``` + +因为纹理坐标是在(0,1)的,所以需要分别乘上一个宽度和高度,再经过把float -> int,直接向下取整得到纹理坐标。 + +> 关于这里为什么是 (1-v.textureCoord[1]),是因为uv坐标,也就是纹理坐标的开始位置(0,0)在左下角,而我们读取的图片的开始位置(0,0)在左上角,所以在y分量的处理上要这样做。 + +同时,在计算对应像素的index的时候记得要乘上一个3,也就是通道数。 + +>tips:要实现更精细的采样纹理坐标,还可以做双线性插值或者三线性插值,这里先不弄了。 + + + +## 3. shader编写 + +根据不同的着色频率,我们可以编写不同的shader。 + +在Blinn-Phong着色模型中,根据不同的着色频率,可以分为三个shader: + +1. Flat shader,计算三角形的法线,为三角形进行着色。 +2. Gouraud shader, 为每个顶点计算法线,计算出每个顶点的颜色然后插值出每个像素的颜色。 +3. Phong shader,为每个像素插值出法线,直接为每个像素计算着色。 + +这里我们仅仅实现下面两个,也就是Gouraud shader以及Phong shader。 + +### 写在之前 + +关于各种类的设计是非常多的,我这里对与shader这个类的设计可能不是非常优雅,但是基本可以实现功能。不过在这里还是稍微阐述一下shader的设计思路。 + +首先渲染是交给渲染管线来做的: + +```cpp +class ShadingPipeline { +public: + ShadingPipeline(int _width,int _height,Shader* _shader); + void shade(const std::vector& _vertices, + const std::vector& _indices, + int rasterizingMode); + void clearBuffer(); + uint8_t* getResult(); + Shader* shader; +private: + int width; + int height; + uint8_t* image; + float* zBuffer; + glm::mat4 viewPortMatrix; + void bresenhamLineRasterization(VertexOutData& from, VertexOutData& to); + void fillRasterization(VertexOutData& v1, VertexOutData& v2, VertexOutData& v3); +}; +``` + +每一个渲染管线有一个Shader变量,Shader变量负责提供vertexShader以及fragmentShader在渲染管线里进行调用。 + +这里其实设计得不是很好,因为我这里给渲染管线的Shader变量的类型是Shader的基类,也就是说,如果后续我写的一个Shader类继承了这个基类,并且定义了一些变量以及函数,那么在渲染管线里也还是拿不到的。 + +所以没有办法,我就只好在Shader的基类里把所有可能用到的变量以及函数先定义好,然后在Shader基类的继承子类里只需要去实现 vertexShader 以及 fragmentShader这两个基函数即可。 + +```cpp +class Shader { +public: + Shader() = default; + virtual VertexOutData vertexShader(VertexData& v) {VertexOutData _v; return _v;}; //meaningless implementation + virtual glm::vec4 fragmentShader(VertexOutData& v) {glm::vec4 r; return r;}; //meaningless implementation + void setModelMatrix(glm::mat4 matrix); + void setViewMatrix(glm::mat4 matrix); + void setProjectionMatrix(glm::mat4 matrix); + + //Gouraud and phong shader needed + void setMaterial(Material _material); + void setEyePos(glm::vec3 _eyePos); + void addDirectionLight(DirectionLight* light); + void addPointLight(PointLight* light); + void setTexture(std::string _textureUrl); + +protected: + glm::mat4 modelMatrix; + glm::mat4 viewMatrix; + glm::mat4 projectionMatrix; + + //Gouraud and phong shader needed + std::string textureUrl; + Image texture; + std::vector directionLights; // multiple lights + std::vector pointLights; + Material material; + glm::vec3 eyePos; +}; +``` + +比如这里定义的很多变量,例如 material, eyepos什么的,在 THREE_D_SHADER 里是完全用不到的。 + + + +### 3.1 Gouraud shader + +首先要为shader添加光源,从上面可以看到,shader里有两个vector,一个用来存放平行光,一个用来存放点光源。通过调用 `addDirectionLight` 以及 `addPointLight`函数,我们就可以为shader添加光源。 + +> 在使用Blinn-Phong着色模型时,必须要添加光源,因为颜色其实就是光与材质作用的结果。 + +简单展示一下光源的代码: + +```cpp +class Light{ +public: + glm::vec3 ambient; + glm::vec3 diffuse; + glm::vec3 specular; + +}; + +class DirectionLight : public Light{ +public: + DirectionLight(glm::vec3 _ambient, glm::vec3 _diffuse, glm::vec3 _specular, glm::vec3 _direction); + glm::vec3 direction; +}; + +class PointLight : public Light{ +public: + PointLight(glm::vec3 _ambient, glm::vec3 _diffuse, glm::vec3 _specular, glm::vec3 _position); + void setLightPos(glm::vec3 _position); + glm::vec3 position; +}; +``` + +然后就是 Gouraud shader 对于自己的vertex shader以及fragment shader的具体实现了,具体如下所示: + +```cpp +VertexOutData GouraudShader::vertexShader(VertexData &v) { + VertexOutData vod; + vod.position = glm::vec4(v.position,1.0f); + vod.position = projectionMatrix * viewMatrix * modelMatrix * vod.position; //mvp transformation + vod.normal = v.normal; + + // 计算着色 + glm::vec3 ka = material.ka; + glm::vec3 kd = material.kd; + if(textureUrl != "none"){ + int xIdx = v.textureCoord[0] * texture.width; + int yIdx = (1-v.textureCoord[1]) * texture.height; + int idx = (yIdx * texture.width + xIdx)*3; + kd = glm::vec3(texture.pixels[idx+0]/255.0f, + texture.pixels[idx+1]/255.0f, + texture.pixels[idx+2]/255.0f); + } + glm::vec3 ks = material.ks; + + glm::vec3 color(0.0f,0.0f,0.0f); + + //directional lights + for(const auto& light:directionLights){ + glm::vec3 lightDir = glm::normalize(-light->direction); // pos -> light + glm::vec3 viewDir = glm::normalize(eyePos - v.position); // pos -> view + glm::vec3 half = glm::normalize(lightDir + viewDir); // 半程向量 + glm::vec3 L_a = light->ambient * kd; + glm::vec3 L_d = light->diffuse * kd * std::max(0.0f,glm::dot(v.normal,lightDir)); + glm::vec3 L_s = light->specular * ks * pow(std::max(0.0f,glm::dot(v.normal,half)),material.shininess); + + + color += (L_a + L_d + L_s)*255.0f; + } + + //point lights + for(const auto& light:pointLights){ + glm::vec3 lightDir = glm::normalize(light->position - v.position); // pos -> light + glm::vec3 viewDir = glm::normalize(eyePos - v.position); // pos -> view + glm::vec3 half = glm::normalize(lightDir + viewDir); // 半程向量 + + float distance = MathUtils::calPoint2PointSquareDistance(v.position,light->position); + + glm::vec3 L_d = light->diffuse * kd * std::max(0.0f,glm::dot(v.normal,lightDir)); + glm::vec3 L_s = light->specular * ks * pow(std::max(0.0f,glm::dot(v.normal,half)),material.shininess); + + + color += (L_d/distance + L_s/distance)*255.0f; + } + + + color = glm::vec3(std::min(color[0],255.0f), std::min(color[1],255.0f),std::min(color[2],255.0f)); + + vod.color = glm::vec4(color,1.0f); + + return vod; +} + +glm::vec4 GouraudShader::fragmentShader(VertexOutData &v) { + return v.color; +} +``` + +具体没什么好说的,关于Blinn-Phong着色模型的原理可以看一下[这里](https://www.ljhblog.top/CG/GAMES101/assignment3.html)。需要注意的是 Gouraud shader 的主要代码部分是在 vertexShader,而fragmentShader就直接传递颜色即可。 + +最后实现出来的结果如下所示: + +![1](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202206011450942.gif) + +可以看到,地面没有被照亮,因为点光源里地板的四个顶点太远了,点光源对于四个顶点的颜色贡献不大。 + +![2](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202206011452620.gif) + +加载纹理的情况下,Gouraud shader的效果也不是很好,因为等于说纹理坐标都没有进行插值,而是直接拿纹理坐标去查询,效果肯定不会特别好。 + + + +### 3.2 Phong shader + +对于Phong shader,其实基本上就是把 Gouraud shader 的vertex shader搬到了fragment shader中去。 + +> 稍微提一下,在把世界坐标,纹理坐标等传递给fragment shader之前是需要进行插值的,这里我还是选择了使用了二维的中心坐标插值,可能会有一点问题。 + +shader的代码主要如下所示: + +```cpp +VertexOutData PhongShader::vertexShader(VertexData &v) { + VertexOutData vod; + vod.position = glm::vec4(v.position,1.0f); + vod.position = projectionMatrix * viewMatrix * modelMatrix * vod.position; //mvp transformation + vod.worldPos = modelMatrix * glm::vec4(v.position,1.0f); // vec3 = vec4... + vod.normal = v.normal; + vod.textureCoord = v.textureCoord; + + return vod; +} + +glm::vec4 PhongShader::fragmentShader(VertexOutData &v) { + // 计算着色 + glm::vec3 ka = material.ka; + glm::vec3 kd = material.kd; + if(textureUrl != "none"){ + int xIdx = v.textureCoord[0] * texture.width; + int yIdx = (1-v.textureCoord[1]) * texture.height; + int idx = (yIdx * texture.width + xIdx)*3; + kd = glm::vec3(texture.pixels[idx+0]/255.0f, + texture.pixels[idx+1]/255.0f, + texture.pixels[idx+2]/255.0f); + } + glm::vec3 ks = material.ks; + + glm::vec3 color(0.0f,0.0f,0.0f); + + //directional lights + for(const auto& light:directionLights){ + glm::vec3 lightDir = glm::normalize(-light->direction); // pos -> light + glm::vec3 viewDir = glm::normalize(eyePos - v.worldPos); // pos -> view + glm::vec3 half = glm::normalize(lightDir + viewDir); // 半程向量 + glm::vec3 L_a = light->ambient * kd; + glm::vec3 L_d = light->diffuse * kd * std::max(0.0f,glm::dot(v.normal,lightDir)); + glm::vec3 L_s = light->specular * ks * pow(std::max(0.0f,glm::dot(v.normal,half)),material.shininess); + + + color += (L_a + L_d + L_s)*255.0f; + } + + //point lights + for(const auto& light:pointLights){ + glm::vec3 lightDir = glm::normalize(light->position - v.worldPos); // pos -> light + glm::vec3 viewDir = glm::normalize(eyePos - v.worldPos); // pos -> view + glm::vec3 half = glm::normalize(lightDir + viewDir); // 半程向量 + + float distance = MathUtils::calPoint2PointSquareDistance(v.worldPos,light->position); + + glm::vec3 L_d = light->diffuse * kd * std::max(0.0f,glm::dot(v.normal,lightDir)); + glm::vec3 L_s = light->specular * ks * pow(std::max(0.0f,glm::dot(v.normal,half)),material.shininess); + + + color += (L_d/distance + L_s/distance)*255.0f; + } + + + color = glm::vec3(std::min(color[0],255.0f), std::min(color[1],255.0f),std::min(color[2],255.0f)); + + + + return glm::vec4(color,1.0f); +} +``` + +最后实现出来的结果如下所示: + +![phong_shading](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202206011457405.gif) + +可以看到,地面有明显的光照反射了,因为现在是逐像素着色了。 + +![rock](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202206011459984.gif) + +同时,石头的纹理也更加清楚了,因为现在是把纹理坐标插值后再去查询颜色。 + + + +## 4. 最后 + +目前这个软渲染器差不多就实现这么多了,不会再添加例如PCF,PCSS这些新功能了,基本上触及到我的知识盲区了。同时现在跑一个Phong shader已经非常卡了,感觉一秒只有5帧吧,当然这和我写的不好也有问题,期间基本上没考虑性能的问题,越写到后面越放飞自我,然后结构也开始比较混乱了。。。 + +不过一路把一个软渲染器从头到尾实现收获还是非常大,让我更加清楚地理解了渲染管线的运作流程,比如说什么地方该做mvp变换,什么地方应该做插值等等。。。虽然有的地方我的代码架构和实际的硬件实现区别会很大,不过我感觉基本上应该也差不大是这么一个流程吧。。。 + +总之暂时先这样吧。。。🤓 + + + diff --git a/CG/imgui_giveup.md b/CG/imgui_giveup.md new file mode 100644 index 0000000..2c3c789 --- /dev/null +++ b/CG/imgui_giveup.md @@ -0,0 +1,55 @@ +--- +title: 我不想用imgui了 +date: 2022-09-15 20:34:29 +tags: [imgui] +--- +# 我不想用imgui了 +弄这个imgui弄了好多天,从最开始跑通它官网上的demo,到琢磨应该怎么把它整合到我的软渲染器里,最后 -> 放弃。 + +弄了加起来可能还是有好几个小时,为了不让这段时间白费,写点东西记录一下。 + +原因主要有三点吧,当然技术上的原因比较多。 + +首先是我觉得自己的这个软渲染器不需要那么多的按钮啥的来交互,我能想到的就是imgui能够有一个很好的显示fps的窗口,然后就是我能够在ui里选择渲染的后端(光栅器/光追),当然现在我连光追一点都没写,感觉ui也不是特别有必要了。 + +其次就是感觉这个项目没有必要导入太多的其他库,不然就太庞大了,也失去了软光栅器那种徒手造轮子的感觉了。现在我已经导入了sdl, assimp这些庞大的库了,已经很难受了。再来一个glew和imgui,库就太多了。 + +最后,当然也是最重要的原因就是我菜,opengl太麻烦了,看得我头皮发麻。 + +在这里稍微整理一下目前遇到的坑吧,如果以后我会了,可能会来填。 + +目前观察到的现象就是,我的软渲染器是通过`SDL_UpdateWindowSurface(window)` 这个函数来实现窗口的像素更新的,我还在stackoverflow上[问了一波](https://stackoverflow.com/questions/73726756/sdl-updatewindowsurface-cover-the-result-of-sdl-gl-swapwindow),别人的说法是: + +image-20220915204747082 + +意思就是你用windowSurface时就不能再去调其它的3D图形的API了 + +但是Imgui调用的`SDL_GL_SwapWindow(window)`刚好就是 + +```cpp +/** + * Update a window with OpenGL rendering. + * + * This is used with double-buffered OpenGL contexts, which are the default. + * + * On macOS, make sure you bind 0 to the draw framebuffer before swapping the + * window, otherwise nothing will happen. If you aren't using + * glBindFramebuffer(), this is the default and you won't have to do anything + * extra. + * + * \param window the window to change + * + * \since This function is available since SDL 2.0.0. + */ +extern DECLSPEC void SDLCALL SDL_GL_SwapWindow(SDL_Window * window); +``` + +所以这两个函数应该是不能同时调的,我的实验结果就是`SDL_UpdateWindowSurface(window)`会覆盖掉`SDL_GL_SwapWindow(window)`,不论先后。 + +中途我还尝试了一下用opengl的`glDrawPixels`搞了半天,但是没搞定,我也不知道咋回事,opengl真不熟。 + +但是解决方案其实也是有的,比如说把要展示的图片弄成一个opengl的2D贴图,然后传给`Imgui::Image()`,但是,又是opengl... 我顶不住了orz。 + +有人也给出了相关的实现[【Dear imgui】有关Image的实现_Gary的面包屑小道的博客-CSDN博客](https://blog.csdn.net/DY_1024/article/details/105690739),但是最后好像他也没能够把图片弄在背景里,而是搞到一个小窗口里了。 + +总之先这样吧,反正我个人感觉软渲染器加上Imgui会有点臃肿而且**麻烦**,如果有一天我熟练掌握了opengl,或者vulkan...? 到时候这可能也不是啥问题了。🤔 \ No newline at end of file diff --git a/CG/imgui_run.md b/CG/imgui_run.md new file mode 100644 index 0000000..dbcef4d --- /dev/null +++ b/CG/imgui_run.md @@ -0,0 +1,52 @@ +--- +title: 如何将imgui跑起来 +date: 2022-09-05 15:58:31 +tags: [imgui] +--- + +# 如何将imgui跑起来 + +不得不说imgui的新手指引真的太劝退了,居然没有一个傻瓜式的直接用cmake编译的例子,这对于cmake小白(我)来说太不友好了。加之imgui本身支持很多后端(opengl, vulkan, metal...),第一眼看过去根本不知所措,我连要include哪些文件都不知道。 + +于是,在经过一天多的摸索后,我终于把imgui跑起来了。 + +image-20220905202805000 + +首先,如果要将imgui导入你的项目,你需要将整个imgui项目导入你的项目中,比如这样: + +image-20220905203020499 + +然后你需要编写相关的CmakeLists... + +```cmake +cmake_minimum_required(VERSION 3.21) +project(imgui_demo) + +set(CMAKE_CXX_STANDARD 14) + +find_package(glfw3 REQUIRED) # 添加glfw +find_package(GLEW REQUIRED) # 添加glew + +include_directories(imgui) +include_directories(imgui/backends) + +add_executable(imgui_demo + main.cpp + imgui/imgui.cpp + imgui/imgui_draw.cpp + imgui/imgui_tables.cpp + imgui/imgui_demo.cpp + imgui/imgui_widgets.cpp + imgui/backends/imgui_impl_glfw.cpp + imgui/backends/imgui_impl_opengl3.cpp + ) + +target_link_libraries(imgui_demo glfw GLEW::GLEW) +``` + +你需要确定自己要用什么后端来跑,这里我们就选择最简单的opengl好了,那么你就需要导入 glfw+glew(这里我是这么选的,其实也可以glfw + glad ?),总之你需要确定自己需要确定自己要导入什么库。 + +接着include相关的文件夹,再在add_executable里确定相关的文件就可以跑了。。。 + +这里的main.cpp我选择了imgui官方提供的example_glfw_opengl3里的main.cpp来跑。 + diff --git a/CPP/README.md b/CPP/README.md new file mode 100644 index 0000000..154b140 --- /dev/null +++ b/CPP/README.md @@ -0,0 +1 @@ +整理一些CPP相关的知识 \ No newline at end of file diff --git a/CPP/compiler_speedUp.md b/CPP/compiler_speedUp.md new file mode 100644 index 0000000..2fcf3f1 --- /dev/null +++ b/CPP/compiler_speedUp.md @@ -0,0 +1,232 @@ +--- +title: 一些CPP简单的提速方法 +date: 2022-06-18 18:00:00 +tags: [CPP] +--- +# 一些CPP简单的提速方法 + +最近看了小彭老师的c++课[第四讲](https://www.bilibili.com/video/BV12S4y1K721/?spm_id_from=333.788&vd_source=c0c1ccbf42eada4efb166a6acf39141b), 讲的是c++的编译器优化,打开了新世界的大门。 + +在这里简单记录一下作业里面我用到的优化方法。 + +> 在windows wsl2平台进行测试,CPU为i7-6700HQ, 从上至下迭代优化 +> +> > 为什么不在macOS测试捏?macOS SOA居然比AOS慢,完全搞不懂为什么,放弃。 + +### 不做任何优化:1570ms + +首先贴一下官方的源程序: + +```cpp +#include +#include +#include +#include +#include + +float frand() { + return (float)rand() / RAND_MAX * 2 - 1; +} + +struct Star { + float px, py, pz; + float vx, vy, vz; + float mass; +}; + +std::vector stars; + +void init() { + for (int i = 0; i < 48; i++) { + stars.push_back({ + frand(), frand(), frand(), + frand(), frand(), frand(), + frand() + 1, + }); + } +} + +float G = 0.001; +float eps = 0.001; +float dt = 0.01; + +void step() { + for (auto &star: stars) { + for (auto &other: stars) { + float dx = other.px - star.px; + float dy = other.py - star.py; + float dz = other.pz - star.pz; + float d2 = dx * dx + dy * dy + dz * dz + eps * eps; + d2 *= sqrt(d2); + star.vx += dx * other.mass * G * dt / d2; + star.vy += dy * other.mass * G * dt / d2; + star.vz += dz * other.mass * G * dt / d2; + } + } + for (auto &star: stars) { + star.px += star.vx * dt; + star.py += star.vy * dt; + star.pz += star.vz * dt; + } +} + +float calc() { + float energy = 0; + for (auto &star: stars) { + float v2 = star.vx * star.vx + star.vy * star.vy + star.vz * star.vz; + energy += star.mass * v2 / 2; + for (auto &other: stars) { + float dx = other.px - star.px; + float dy = other.py - star.py; + float dz = other.pz - star.pz; + float d2 = dx * dx + dy * dy + dz * dz + eps * eps; + energy -= other.mass * star.mass * G / sqrt(d2) / 2; + } + } + return energy; +} + +template +long benchmark(Func const &func) { + auto t0 = std::chrono::steady_clock::now(); + func(); + auto t1 = std::chrono::steady_clock::now(); + auto dt = std::chrono::duration_cast(t1 - t0); + return dt.count(); +} + +int main() { + init(); + printf("Initial energy: %f\n", calc()); + auto dt = benchmark([&] { + for (int i = 0; i < 100000; i++) + step(); + }); + printf("Final energy: %f\n", calc()); + printf("Time elapsed: %ld ms\n", dt); + return 0; +} +``` + +### 加-O3: 没啥用 + +是的,加O3优化没起作用,我也不知道为啥。。。 + +```Cmake +set(CMAKE_CXX_FLAGS "-O3") +``` + +### 加-ffast-math: 1190ms + +增加`-ffast-math`指令可以让GCC更大胆地尝试浮点运算的优化,如果能够保证程序中不会出现无穷大和NAN。 + +```cmake +set(CMAKE_CXX_FLAGS "-O3 -ffast-math") +``` + +### 加-march=native: 970ms + +通过 `-march=native` 来让编译器检测电脑是否支持AVX,来实现更快的SIMD。 + +```cpp +set(CMAKE_CXX_FLAGS "-O3 -ffast-math -march=native") +``` + +### AOS改为SOA:240ms + +AOS(Array of Struct): 传统的一般的面向对象的编程方式,例如上面代码中的: + +```cpp +struct Star { + float px, py, pz; + float vx, vy, vz; + float mass; +}; + +std::vector stars; +``` + +因为内存排布不友好,所以很难做SIMD。 + +SOA(Struct of Array): 高性能计算中常见的编程方式,可以将上面的写法更改为: + +```cpp +template +struct Star { + float px[N]; + float py[N]; + float pz[N]; + float vx[N]; + float vy[N]; + float vz[N]; + float mass[N]; +}; +``` + +这样的内存排布友好,所以方便做SIMD。 + +同样的,循环访问变量时的代码也要进行修改,例如: + +```cpp +for(size_t i=0;i<48;i++){ + for(size_t j=0;j<48;j++){ + float dx = stars.px[j] - stars.px[i]; + float dy = stars.py[j] - stars.py[i]; + float dz = stars.pz[j] - stars.pz[i]; + float d2 = dx * dx + dy * dy + dz * dz + eps * eps; + d2 *= std::sqrt(d2); + stars.vx[i] += dx * stars.mass[j] * G * dt / d2; + stars.vy[i] += dx * stars.mass[j] * G * dt / d2; + stars.vz[i] += dx * stars.mass[j] * G * dt / d2; + } + } +``` + +可以看到优化效果是非常明显的(但是在mac上不是。。至少在M1的mac上不是,实测还变慢了,不知道是clang/LLVM的锅还是M1的锅。。。) + +### sqrt改为std::sqrt: 170ms + +传统的一些数学函数是C语言的遗产,不要用了,比如说sqrt函数只支持double类型,自然就会很慢。 + +所以这里修改为std::sqrt,这是C++的数学函数。 + +```cpp +d2 *= std::sqrt(d2); +``` + + + +### 将有些变量改为constexpr,有的函数改为static: 没啥用 + +```cpp +static float frand() { + return (float)rand() / RAND_MAX * 2 - 1; +} + +constexpr float G = 0.001; +constexpr float eps = 0.001; +constexpr float dt = 0.01; +``` + +constexpr是可以把一些变量直接进行编译器求值,static应该也差不多吧(我也不是很确定)。 + +但是不知道为什么加了没用。。。又或者说有一点用但是不明显吧。 + + + +## 总结 + +从1570ms最后优化到了170ms,接近10倍,牛逼疯了。 + +我感觉我以前写的代码可能全是垃圾代码。。。 + + + + + + + + + + + diff --git a/CPP/handcrafted_vector.md b/CPP/handcrafted_vector.md new file mode 100644 index 0000000..4a84d48 --- /dev/null +++ b/CPP/handcrafted_vector.md @@ -0,0 +1,198 @@ +--- +title: 手写一个vector +date: 2023-09-19 22:23:00 +tags: [vector, cpp] +--- +# 手写一个vector +> 日期: 2023-09-19 + +前两天在b站刷到一个手写vector的视频,于是就照着写了一下,新的c++知识增加了。。。 + +[手写一个 std::vector 可以有多复杂?_哔哩哔哩_bilibili](https://www.bilibili.com/video/BV1iX4y1w7x4/?spm_id_from=333.999.0.0&vd_source=c0c1ccbf42eada4efb166a6acf39141b) + +直接上代码(似乎太长就没有了高亮。。。?将就看吧): + +```cpp +#include +using namespace std; + +template +class _vector{ +public: + using value_type = T; + using reference_type = T&; + using const_reference_type = const T&; + using iterator = T*; //其实还可以有 const_iterator,但是这里不多实现了 +private: + T* m_data; + size_t m_size; + size_t m_capacity; +public: + // 默认构造函数 + constexpr _vector() noexcept : m_data(), m_size(), m_capacity() {} + // 析构函数 + ~_vector(){ + // 先去调用每一个元素的析构函数,然后再删除这一段内存 + for(size_t i=0; i(operator new(rhs.m_capacity * sizeof(T))); + // 对于每一个元素调用它的拷贝构造函数 + for(size_t i = 0; iemplace_back(rhs); + } + void push_back(T&& rhs){ + cout<<"push_back func2"<emplace_back(std::move(rhs)); + } + + template + T& emplace_back(ArgsT&&... args){ + if(m_size < m_capacity){ + m_size += 1; + }else{ + //扩容 + // 先申请一段新的内存 + m_capacity = ceil((m_capacity+1) * 1.0f * 1.5); + T* m_data_new = static_cast(operator new(m_capacity * sizeof(T))); //记得static cast + // 把原本内存的元素移动到新的地方 + for(int i=0;i(args)...); + return *cur; + } +}; + +void test1(){ + _vector _v; + _v.push_back(1); + _v.push_back(2); + _v.push_back(3); + _v.push_back(4); + _v.push_back(5); + _v.push_back(6); + _v.push_back(7); + cout<<_v.size()<<" "<<_v.capacity()< _v2(_v); + cout<<_v2.size()<<" "<<_v2.capacity()< 首先分配和传入对象大小相同的一段内存(**怎么分配一段内存**),然后对于内存的每一个位置,对于每一个元素调用拷贝构造函数来新建对象(**怎么在一个指定地址新建对象**)。 +2. 析构函数怎么实现 -> 先调用每一个元素的析构(**怎么调用析构**),然后delete这段内存。 +3. 因为push_back就是去调emplace_back,所以重点在于emplace_back应该怎么实现。不过还需要注意的是,对于接受右值引用为参数的push_back在继续向下传递时,要使用std::move,不然的话传进去就变左值了。 +4. emplace_back怎么实现 -> 其实就是扩容后,需要重新分配更大的内存,然后把原来的元素通过移动构造的方式在新的内存新建。分配完后在末尾原地构造新元素,因为emplace_back可以接收各种参数,所以这里如何构造取决于外部传入的参数(如果是一堆参数那么就是普通构造,如果是左值引用那么就是拷贝构造,如果是右值引用那么就是移动构造。。。) + + + +一些新学习的语法: + +* `using value_type = T;` 其实和typedef 差不多 +* `m_data[i].~T();` 调用析构 +* `m_data = static_cast(operator new(rhs.m_capacity * sizeof(T)));` 使用operator new 分配一段内存 +* `new(&m_data[i]) T(rhs.m_data[i]);` 在指定地址新建对象 + diff --git "a/Debug/cuda\345\222\214cpu\347\232\204\350\256\241\347\256\227\350\257\257\345\267\256.md" "b/Debug/cuda\345\222\214cpu\347\232\204\350\256\241\347\256\227\350\257\257\345\267\256.md" new file mode 100644 index 0000000..4fe2f61 --- /dev/null +++ "b/Debug/cuda\345\222\214cpu\347\232\204\350\256\241\347\256\227\350\257\257\345\267\256.md" @@ -0,0 +1,38 @@ +--- +title: 记录cuda和cpu的计算误差带来的坑 +date: 2022-10-28 17:07:05 +tags: [cuda] +--- +# 记录cuda和cpu的计算误差带来的坑 +> 日期: 2022-10-28 + +最近在用torchac糊demo,之前光速写完了然后一直发现编码和解码的点云对应不上,然后这个bug我de了好多天。 + +最开始以为是库的精度问题带来的误差,然后我把神经网络的输入输出都打出来看了,发现编码端和解码端的神经网络的输入输出都是一样的,这让我非常迷惑。 + +由于实在不知道啥情况了,我打算把输出输出用pickle存入文件,然后单独写一个py文件来debug。结果发现编码解码端的神经网络输出是一样的,而且在单独的文件里做算术编解码时,输出的结果居然都是对的。 + +这让我开始怀疑是不是编码的bytestream出了问题,于是我将实际编码以及单独文件里编码的bytestream打出来看了一下,发现居然不一样。 + +这让我非常的震惊,明明两个函数的输入都是一样的,唯一的区别在于一个是直接跑的,另外一个是被我存进了文件,然后再读出来再跑的。 + +最后发现坑在于.cuda()和.cpu()。原本我有一个变量是存在gpu上的,但是在保存为文件后,就默认保存为.cpu()格式了,然后再读取时也没有将其转为.cuda()。 + +具体的问题呢应该是出在这里,我试了一下,输入的pdf分别是.cuda()和.cpu()时,结果是不一样的。 + +行吧。。。 +```python +def pdf2cdf(pdf): + # if symbols nums is N, need cdf of length of N+1 + # pdf -> cdf + # B * N * 256 B * N * 257 + # https://blog.csdn.net/OrdinaryMatthew/article/details/125757080 + # pdf = pdf/torch.sum(pdf,-1,keepdim=True) # normalize, 如果本身输入的pdf之和为1,则不需要这行 + pdf = F.softmax(pdf,dim=-1) + # print(pdf) + cdf = torch.cumsum(pdf, -1) + cdf = torch.cat([torch.zeros_like(cdf[:,:,:1]),cdf], -1) # add zero in the front + cdf = torch.cat([cdf[:,:,:-1],torch.ones_like(cdf[:,:,:1])], -1) # add one in the back (其实最后是有0的,但是一般来说会超过1一点点,所以这里直接把最后的1.000几几直接替换为1) + return cdf + +``` \ No newline at end of file diff --git "a/Debug/pybind\345\206\205\345\255\230\346\263\204\346\274\217.md" "b/Debug/pybind\345\206\205\345\255\230\346\263\204\346\274\217.md" new file mode 100644 index 0000000..43b1140 --- /dev/null +++ "b/Debug/pybind\345\206\205\345\255\230\346\263\204\346\274\217.md" @@ -0,0 +1,68 @@ +--- +title: 记录一次使用pybind内存泄漏的debug过程 +date: 2022-03-30 +tags: [pybind] +--- +> 日期:2022.3.30 + +# 记录一次使用pybind内存泄漏的debug过程 + +最初发现bug是因为使用nohup命令炼丹时,炼着炼着任务就自己停了,没有任何的报错信息,当时十分迷惑,于是又跑了一次,发现任务又在差不多的位置停了。当时就断定应该是内存泄漏,操作系统把任务给我杀了,于是就开始查看为啥会内存泄漏。 + +最开始以为是pytorch的问题,比如说一些梯度的累加啥的,当时就在网上搜到说 loss 的记录与累加要写成: + +``` +totalLoss += loss.item() +``` + +的形式,然后因为我之前写的是: + +``` +totalLoss += loss.data +``` + +于是我就改了,事实发现没什么卵用。 + +于是开始在网上搜一些内存泄漏相关的debug工具,于是搜到了 `memory_profiler`。 + +这个真的是个神器,只需要在函数前面加上 `@profile` 修饰符,就可以在函数执行结束时查看这个函数执行过程中内存的变化。 + +最开始查看当然从 `train` 入手,于是操作一波发现每一次的内存增加都来自于 `dataLoader` ,当时我就在猜测应该是 `pybind` 的问题了,毕竟自己用 `pybind` 魔改了一波`dataLoader`,很难不保证不出什么幺蛾子。 + +> pybind是一个可以在python调用c++的一个库,比numba等直接在python里面写的库更加自由,因为pybind能够直接编写c++的代码(埋下伏笔),但是相对来说也更难。 + +> 当初我的需求是使用c++来编写相关的代码来加速点云转换为八叉树,因为在python里直接转八叉树实在是太慢了,所以当初就想当然地在技术选型时使用了pybind,事实证明,这是一个还不错的决定,如果你的c++写得很好,而不是像我一样瞎用指针。 + +于是再次深入,去 `dataLoader` 里使用`@profile` 修饰符查看内存的变化情况,最后发现每一次的内存增加都来自于一个 我用pybing写的函数,这时候就基本确定了bug来自于何处。 + +仔细一看代码: + +```cpp +process_info* bfs_process_octree(shared_ptr root) { + + map nodedict; nodedict.clear(); + map layerIndexs; layerIndexs.clear(); + queue> octantQueue; + .... +} +``` + +好家伙,到处都是指针,当时我就在猜测应该是内存释放出了问题,所以我开始尝试一些例如 `del`, `gc.collecte()` 的一些操作,试图在python里把内存释放掉,然而没什么用。 + +当初我写指针也是想着既然C++和python数据通信,地址一定会更快一些,于是在面对一些数据量比较大的变量时都是用的指针,结果没想到因为内存释放出了大问题。 + +所以最后没办法,我就把所有的指针又改回了对象...(有点low) + +```cpp +process_info bfs_process_octree(shared_ptr root) { + + map nodedict; nodedict.clear(); + map layerIndexs; layerIndexs.clear(); + queue> octantQueue; + .... +} +``` + +不过问题解决了,而且好像速度并没有变得更慢🤔? + +总之以后要慎用指针了... \ No newline at end of file diff --git a/Debug/readme.md b/Debug/readme.md new file mode 100644 index 0000000..65f7ebd --- /dev/null +++ b/Debug/readme.md @@ -0,0 +1 @@ +记录一些Debug/踩坑日常 \ No newline at end of file diff --git a/DistributedSys/README.md b/DistributedSys/README.md new file mode 100644 index 0000000..80bcb89 --- /dev/null +++ b/DistributedSys/README.md @@ -0,0 +1 @@ +[MIT6.824课程记录](./mit6.824/readme.md) \ No newline at end of file diff --git a/DistributedSys/mit6.824/lab1.md b/DistributedSys/mit6.824/lab1.md new file mode 100644 index 0000000..bf2217c --- /dev/null +++ b/DistributedSys/mit6.824/lab1.md @@ -0,0 +1,95 @@ +--- +title: MIT6.824-lab1-MapReduce +date: 2022-04-20 +tags: [MIT6.824] +--- +> 日期: 2022.4.20 + +# mit6.824 lab1 + +lab1是要求实现mapreduce,下面就简单记录一下做这个lab的一些注意点。 + +## 1. map function & reduce function + +首先要搞明白map reduce 里的map function 以及reduce function都在干什么。 + +### map + +拿这里的wordcount的为例,map函数要实现的就是对与每一个输入的filename,生成 形式为 <单词,1>的键值对,并且要将其保存为一个中间形式(intermediate)。在这里,我们将中间形式保存为文件,并且将其命名为 mr-X-Y,X是map task的id,Y是reduce task的id,也就是说,当一个中间文件形成时,就已经决定了它应该被哪个reduce task来处理了。 + +需要注意的是,在产生中间文件时,mr-X-Y的Y的选择很重要。在lab中提供了一个 `ihash` 函数,这个函数的作用是输入一个字符串,返回一个int的数,然后你用这个int的数去mod你的reduce task的数目就可以决定mr-X-Y的Y。最开始我是直接拿了要处理的txt文件名去做hash,也就是说一个map任务的所有处理结果全部写入一个文件,**但其实这是错误的**。 + +实际上在map reduce中是有一个shuffle阶段的,这个阶段就是说对于map产生的结果做一个shuffle,来使得在做reduce task的时候能够让一个reduce task能够在不同的中间结果中处理相同的部分内容。就比如说,reduce task 1只去处理对于字符 a,b,c的数量统计,reduce task 2只去处理对于字符 e,f,g的数量统计, 不过都是针对所有的中间结果的。 + +如果使用文件名来做hash,可以预想到的结果是,每一个Y里的中间结果都是差不多的,这是我们不想看到的,因为这样的话reduce task就没办法去挑自己想处理的文件来处理了。怎么办呢?其实这里应该是拿word去做hash,也就是说,一个map任务的处理结果会存到多个中间文件里去,这样一想就很清楚了。对于0号的map task,它的结果可能会存到 mr-0-1, mr-0-2。。。然后mr-X-1里全是存的 <单词1,1> <单词2, 1>, 然后mr-X-2里全是存的<单词3, 1> <单词4, 1>这种。 + +### reduce + +reduce task 要实现的就是对于一堆输入 <单词,1>的键值对,将他们做一个求和,比如说把 <单词1,1><单词1,1><单词1,1>变成 <单词1,3>,就这么简单。 + + + +## 2. worker & master + +这里我写的有点糙,感觉有一堆bug,尤其是crash test,就偶尔能过,不过摆了。 + +### worker + +先来说worker的设计,相对来说简单一些。 + +worker主要就是通过rpc来向master要任务,如果说是map任务,那么就执行上面说的map task需要做的东西;如果是reduce任务,那么就执行上面说的reduce task需要做的东西;如果要等待,就啥也不干嗯等;不然就直接退出了。 + + + +### master + +接下来就是master了,感觉这就是一个加锁地狱,尤其是当调crash test时。 + +首先你需要给master维护一堆变量,在创建时我们需要创建: + +```go +type Master struct { + // Your definitions here. + mapFiles []string // txt文件名数组 + mapTasksTaken map[int]bool // map task 被taken了的记录 + mapOKTable map[int]bool // map task 做完了的记录 + mapTaskTimestamp map[int]int // map任务开始时间 + nReduce int // reduce task的数目 + reduceTasks map[int][]string // int -> reduce task的数目 []string -> reduce task 要处理的中间文件名 + reduceTasksTaken map[int]bool // reduce task 被taken了的记录 + reduceOKTable map[int]bool // reduce task 做完了的记录 + reduceTaskTimestamp map[int]int // 记录reduce任务开始时间 +} +``` + +基本上需要维护上面这些变量。 + +此外,你还需要创建很多锁: + +```go +var mapOkMutex = sync.Mutex{} +var reduceOkMutex = sync.Mutex{} +var reduceTakenMutex = sync.Mutex{} +var mapTakenMutex = sync.Mutex{} +var mapTimestampMutex = sync.Mutex{} +var reduceTimestampMutex = sync.Mutex{} +var reduceTasksMutex = sync.Mutex{} +var allLock = sync.Mutex{} +``` + +其实这一部分我还是不太清楚加这么多锁有没有用,因为最开始我是直接在`GetTask`函数的开头加了一个锁,在函数最后释放锁,就这么粗暴的操作可以通过除了最后一个crash test的test,然而当我加了这么多锁过后,crash test还是过不了,so... + +简单说一下master都干了什么,master的主要任务就是实现 `GetTask`函数,使得在worker向master要任务时,要给它一个合理的安排。 + +首先要知道,在mapreduce里,需要在map任务都做完了过后再去做reduce任务,所以master自己要知道map任务是否做完了,reduce任务是否做完了,于是就这用到了我们上面列出的 `mapOKTable` 和 `reduceOKTable`, 当worker的map或者reduce任务做完后,会通过rpc调用一个master的函数来让master知道任务做完了。其余的部分好像也没啥好说的了,无非是map task没做完,给worker分一个map task,reduce task没做完,给worker分一个reduce task。 + +这里还要注意一下,考虑到crash的情况(虽然到最后我也没有完全解决这个问题),我还维护了两个变量 `mapTasksTaken` 和 `reduceTasksTaken`,也就是说每一个任务会经历 `untaken -> taken -> done` 的过程,并不是说任务分配出去了就做完了。如何判断做完了?-> 调用了rpc来让master知道做完了就做完了。如何判断没做完? -> 这里还记录了时间戳 `mapTaskTimestamp` 和 `reduceTaskTimestamp`,超过了一定的时间没有做完那么任务就会又变为 untaken,然后就会分配给别人。ps: 不过这里还有bug,直到最后我还发现会出现都开始执行reduce了,还有map task完成的情况,那么是不是就是因为 master以为任务挂了分配给了另一个worker,然后其实它没挂,然后就开始鬼畜了。。。所以我觉得加时间戳这个方法不是特别好。。。总之很迷,整个crash test都很迷而且难调,所以。。。摆了。 + + + +## 3. 总结 + +所以这就是整个lab的总结(摆烂实录)了 :) + +![截屏2022-04-20 下午4.13.11](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/202204201710677.png) + diff --git a/DistributedSys/mit6.824/readme.md b/DistributedSys/mit6.824/readme.md new file mode 100644 index 0000000..1afed4f --- /dev/null +++ b/DistributedSys/mit6.824/readme.md @@ -0,0 +1,3 @@ +这里会整理MIT6.824的各个lab的记录 + +* [lab1-mapreduce](lab1.md) diff --git a/MLDLRL/DL/README.md b/MLDLRL/DL/README.md new file mode 100644 index 0000000..3bc3c61 --- /dev/null +++ b/MLDLRL/DL/README.md @@ -0,0 +1,4 @@ +这里记录的都是关于深度学习一些比较general的知识, 不会具体涉及到CV, NLP或者数据挖掘等具体应用。 + +* [batch与并行计算速度的讨论](batch-and-speed.md) +* [使用larq对BNN从训练到部署](bnndemo.md) \ No newline at end of file diff --git a/MLDLRL/DL/bnndemo.md b/MLDLRL/DL/bnndemo.md new file mode 100644 index 0000000..5e66f4c --- /dev/null +++ b/MLDLRL/DL/bnndemo.md @@ -0,0 +1,198 @@ +--- +title: 使用larq对BNN从训练到部署 +date: 2021-05-29 20:30:04 +tags: [二值神经网络] +--- + +由于一些机缘巧合,接触到了二值神经网络,于是它就成为了我毕业设计的选题。经过一番挣扎过后,也算是简单做了一点东西,在此记录一下。 + + + +## 1. 模型的搭建以及训练 + +目前的二值神经网络的搭建与部署一般有两种方法: + +1. 使用pytorch搭建模型,训练后转为onnx,然后使用**dabnn**或者**bolt**进行部署。这种方法的好处在于Pytorch框架本身比较好用,也比较自由。但是在部署时常常会出现dabnn或者bolt转换模型失败的情况,有时还是比较头疼。尤其是在dabnn不再维护,bolt的二值网络仅支持armv7.2指令集以上的手机时,这种方法让我直接走不下去了。 +2. 使用[**larq**](https://github.com/larq/larq)进行搭建并且使用LCE进行部署。这种方法的好处在于模型搭建简单以及模型转换以及部署的简单,这也是我最后使用的方法,唯一的缺点可能在于tensorflow相较于pytorch没那么好用,而且这边的模型搭建一般推荐使用Sequential的方式,总之最好不要用Model子类的方式来搭建模型,因为larq对那种方法支持不是很好。 + +在这里我使用了Bi-Real Net18来进行训练,使用了larq以及tensorflow来搭建模型。模型的搭建大致是参考larq_zoo里的代码,因为我是在CIFAR100数据集上进行训练,所以我对网络的开头做了一点修改,把卷积核以及步长改小了一点,去掉了池化层,这样效果会好一些。模型训练完毕后,将模型保存为.h5格式。 + + + +## 2. 模型转换 + +如果使用larq的话,模型转换这一步就比较简单,通过调用以下几行代码,直接将模型转换为tflite格式即可。 + +```python +import larq as lq +import larq_compute_engine as lce +import tensorflow as tf + +m = tf.keras.models.load_model("mobilev1.h5") +lq.models.summary(m) +with open("mobilenetV1.tflite", "wb") as flatbuffer_file: + flatbuffer_bytes = lce.convert_keras_model(m) + flatbuffer_file.write(flatbuffer_bytes) +``` + + + +## 3. 模型部署 + +模型部署需要使用到LCE,因为最后需要将模型在app中跑起来,所以在这一步中我们的目标是编写相关的C++代码实现推理,将代码使用JNI进行封装,并且最后将代码编译成为.so文件,即动态链接库,然后最后将这个.so文件放入Android Studio来进行调用。详细的流程如下所示。 + +### 3.1 安装LCE + +模型推理的过程,使用了[LCE](https://github.com/larq/compute-engine),直接去github页面看安装过程就好了,他们的文档写得很详细。 + +### 3.2 编写推理代码并编译 + +LCE安装好后,我们在LCE根目录下新建一个文件夹叫做jni_lce,里面用来存放我们即将编写的代码。在该文件夹下编写lce.cc文件。这里截取部分代码这个文件很简单,实现了两个函数,一个是loadModel,一个是predict。此处要注意JNI编程的函数命名规范。举个例子,Java_com_ljh_bnndemo_Net_loadModel这个函数对应的就是com.ljh.bnndemo包下的Net类的loadModel方法。此外,还要注意的是,我这里直接使用了一个全局变量来保存读取的模型(interpreter),这种写法其实很烂,但我也不知道怎么写才能更好了2333。 + +```c++ + +//use a interpreter as a global variable +std::unique_ptr interpreter; + +extern "C" JNIEXPORT jboolean JNICALL +Java_com_ljh_bnndemo_Net_loadModel( + JNIEnv* env, + jobject thiz, + jobject model_buffer) { + + char* buffer = static_cast(env->GetDirectBufferAddress(model_buffer)); + size_t buffer_size = static_cast(env->GetDirectBufferCapacity(model_buffer)); + + // Load model + std::unique_ptr model = + tflite::FlatBufferModel::BuildFromBuffer(buffer,buffer_size); + TFLITE_MINIMAL_CHECK(model != nullptr); + + // Build the interpreter + tflite::ops::builtin::BuiltinOpResolver resolver; + compute_engine::tflite::RegisterLCECustomOps(&resolver); + + InterpreterBuilder builder(*model, resolver); + builder(&interpreter); + TFLITE_MINIMAL_CHECK(interpreter != nullptr); + + // Allocate tensor buffers. + TFLITE_MINIMAL_CHECK(interpreter->AllocateTensors() == kTfLiteOk); + + LOGI("model load succeed!!!"); + + return true; +} + +extern "C" JNIEXPORT jfloatArray JNICALL +Java_com_ljh_bnndemo_Net_predict( + JNIEnv* env, + jobject thiz, + jfloatArray arr) { + + float *jInput; + jInput = env->GetFloatArrayElements(arr, 0); + const jint length = env->GetArrayLength(arr); + + LOGI(".................start to predict...................."); + // // Fill input buffers + // TODO(user): Insert code to fill input tensors + float* input = interpreter->typed_input_tensor(0); + + for(int i=0;i<1024;i++) + { + input[i*3 + 0] = jInput[i*3 + 0]; + input[i*3 + 1] = jInput[i*3 + 1]; + input[i*3 + 2] = jInput[i*3 + 2]; + } + + // Run inference + TFLITE_MINIMAL_CHECK(interpreter->Invoke() == kTfLiteOk); + + // Read output buffers + float* output = interpreter->typed_output_tensor(0); + + //输出 + //CIFAR100对应100分类 + float *log_mel = new float[100]; + for(int i=0;i<100;i++){ + log_mel[i] = output[i]; + } + jfloatArray result = env->NewFloatArray(100); + env -> SetFloatArrayRegion(result,0,100,log_mel); + + LOGI("predict over \n"); + + return result; +} + +``` + +代码编写完毕后,需要编写build文件来进行编译,LCE使用的是bazel来进行管理项目,所以build文件的编写如下 + +``` +load("@org_tensorflow//tensorflow/lite:build_def.bzl", "tflite_linkopts") + +package( + default_visibility = ["//visibility:public"], + licenses = ["notice"], # Apache 2.0 +) + +cc_binary( + name = "liblce.so", + srcs = [ + "lce.cc", + ], + linkopts = tflite_linkopts() + select({ + "@org_tensorflow//tensorflow:android": [ + "-pie", # Android 5.0 and later supports only PIE + "-lm", # some builtin ops, e.g., tanh, need -lm + ], + "//conditions:default": [], + }), + deps = [ + "//larq_compute_engine/tflite/kernels:lce_op_kernels", + "@org_tensorflow//tensorflow/lite:framework", + "@org_tensorflow//tensorflow/lite/kernels:builtin_ops", + ], + linkshared=True, +) +``` + +随后需要对项目进行编译: + +在LCE根目录下使用bazel对该项目进行编译。如下: + + ```shell + bazel build --config=android_arm64 //jni_lce:liblce.so + ``` + +编译后会生成一个LCE/bazel-bin文件夹。LCE/bazel-bin/jni_lce文件夹下找到**liblce.so**动态链接库文件,得到这个文件后,就可以将该文件添加到Android Studio中的工程项目中,并使用相关java进行调用了。 + + + +### 3.3 在Android Studio中调用 + +首先,需要将动态链接库添加到android studio中,此处需要在项目的**main**文件夹下创建**jniLibs**文件夹,随后在jniLibs文件夹下创建**arm64-v8a**文件夹,随后将liblce.so文件放在arm64-v8a文件夹内。 + +然后需要为创建相关的类以及方法来实现对于C++编写的函数的调用。在这里我们在com.ljh.bnndemo包下创建Net类,并且类中相关的代码编写如下: + + ```java + package com.ljh.bnndemo; + + import java.nio.ByteBuffer; + + public class Net { + static { + System.loadLibrary("lce"); + } + + public native boolean loadModel(ByteBuffer modelBuffer); + public native float[] predict(float[] input); + } + ``` + +接下来只需要对这两个方法进行调用就可以了。 + +至此就完成了从BNN的训练到部署啦。关于app的详细代码可以见https://github.com/LJHG/BNNDemo + diff --git a/MLDLRL/DL/pytorch-mps.md b/MLDLRL/DL/pytorch-mps.md new file mode 100644 index 0000000..f76df7e --- /dev/null +++ b/MLDLRL/DL/pytorch-mps.md @@ -0,0 +1,80 @@ +--- +title: mac pytorch gpu加速尝鲜( +date: 2023-12-22 16:00:22 +tags: [mac, pytorch] +--- +# mac pytorch gpu加速尝鲜( +> 日期: 2023-12-22 + +前段时间就听说mac好像支持gpu加速了,当时似乎是pytorch官方发了一个博客[Introducing Accelerated PyTorch Training on Mac | PyTorch](https://pytorch.org/blog/introducing-accelerated-pytorch-training-on-mac/) (卧槽原来已经一年了),我记得刚出的时候是只支持cpu来着,难道说一开始就支持mps了🤔?无所谓了,反正我一直没用过。 + +然后最近突然发现pytorch竟然已经更新到2.1了,pytorch给mac提供了一个叫做MPS (Metal Performance Shaders)的东西,来实现gpu加速(苹果真的得给这些社区磕头)。 + +于是就简单写了个脚本,测试一下(on 2021 mbp m1pro 16g)。。。 + +```python +import torch +import time + +if __name__ == '__main__': + if(torch.backends.mps.is_available()): + print("## enabling mps ##") + device = torch.device("mps") + elif(torch.cuda.is_available()): + print("## enabling cuda ##") + device = torch.device("cuda") + else: + print("## enabling cpu ##") + device = torch.device("cpu") + + layer = torch.nn.TransformerEncoderLayer(512, 8 ,2048) + model = torch.nn.TransformerEncoder(layer, 6).to(device) + start = time.time() + for i in range(10): + res = model(torch.zeros(1024,512).to(device)) + end = time.time() + print("time cost: ", end-start) + + + +``` + +简单写了个transoformer encoder,然后跑了一下 + +``` +# torch 2.1 mps +time cost: 0.658048152923584 + +# torch 1.12 mps(没错,那个时候已经有mps了) +time cost: 1.2606778144836426 + +# torch 2.1 cpu +time cost: 4.325138092041016 + +# torch 1.12 cpu +time cost: 4.158311128616333 + +# torch 1.13 GTX1060 +time cost: 1.8937978 + +# torch 1.13 i7-7700 +time cost: 5.250705 + + + + + +``` + +可以看出来,gpu加速还是会有4~7倍的提升,随着版本更新,cpu推理没啥变化,但是gpu变快了好多。。。可能之前我看的那个blog post就是说的gpu加速变快了,但是我也找不到了。 + +还拿实验室的祖传1060来测了一下,还跑不过m1 pro。。。 + +虽然测得非常不严谨,不过mac上的表现真的还不错了,而且苹果的统一内存也挺香的,配个好点的机器,在果子上炼丹不是梦哇。。 + + + + + + + diff --git a/MLDLRL/ML/README.md b/MLDLRL/ML/README.md new file mode 100644 index 0000000..d23e651 --- /dev/null +++ b/MLDLRL/ML/README.md @@ -0,0 +1,3 @@ +这里记录的都是关于传统机器学习的算法。 + +* [KNN计算的三种姿势](knn-distance-computing.md) \ No newline at end of file diff --git a/MLDLRL/ML/knn-distance-computing.md b/MLDLRL/ML/knn-distance-computing.md new file mode 100644 index 0000000..4be60a1 --- /dev/null +++ b/MLDLRL/ML/knn-distance-computing.md @@ -0,0 +1,176 @@ +--- +title: 计算欧式距离的三种姿势 +date: 2020-11-02 15:58:46 +tags: [机器学习, machine learning, numpy] +--- +最近在看cs231n,第一个作业里的第一个就是KNN,总所周知,KNN中距离矩阵的计算是非常重要的,这个作业里列举了三种距离的计算方式,当时就震惊了,这里来记录一下。 + + + +# 双重循环 + +emmm,最显而易见的方式,不多说,直接上代码。 +ps: 突然发现typora里面的插入代码方式hexo也可以认识,终于可以不用写那个啥奇怪的codeblock了。 + +```python +def compute_distances_two_loops(self, X): + """ + Compute the distance between each test point in X and each training point + in self.X_train using a nested loop over both the training data and the + test data. + Inputs: + - X: A numpy array of shape (num_test, D) containing test data. + - num_test是测试集里样本的数量,D是像素点的个数(这里是3072,3*32*32) + Returns: + - dists: A numpy array of shape (num_test, num_train) where dists[i, j] + is the Euclidean distance between the ith test point and the jth training + point. + """ + num_test = X.shape[0] + num_train = self.X_train.shape[0] + dists = np.zeros((num_test, num_train)) + for i in xrange(num_test): + for j in xrange(num_train): + ##################################################################### + # TODO: # + # Compute the l2 distance between the ith test point and the jth # + # training point, and store the result in dists[i, j]. You should # + # not use a loop over dimension. # + ##################################################################### + dists[i][j] = np.linalg.norm(X[i] - self.X_train[j]) #求欧式距离 + pass + ##################################################################### + # END OF YOUR CODE # + ##################################################################### + return dists +``` + +呃,这里要说一下,计算距离时不是一定要写成np.linalg.norm(X[i] - self.X_train[j]),其实写成和下面的单重循环那种形式一样也可以,没区别。 + +# 单重循环 +单重循环应该说是很好的利用了numpy的特性吧,等于就是让numpy来帮你做了循环,这里举个例子。 + +假设有代码如下: + +```python +a = np.array([[1,2,3]]) +b = np.array([1]) +print(a+b) +``` + +这种相加在numpy中是可行的,结果为: + +```python +[[2 3 4]] +``` + +既然可以这样操作,那么如果将其运用于距离计算中,便可以少写一重循环了。 + +代码如下: + +```python + def compute_distances_one_loop(self, X): + """ + Compute the distance between each test point in X and each training point + in self.X_train using a single loop over the test data. + Input / Output: Same as compute_distances_two_loops + """ + num_test = X.shape[0] + num_train = self.X_train.shape[0] + dists = np.zeros((num_test, num_train)) + for i in xrange(num_test): + ####################################################################### + # TODO: # + # Compute the l2 distance between the ith test point and all training # + # points, and store the result in dists[i, :]. # + ####################################################################### + dists[i] = np.sqrt(np.sum(np.square(self.X_train-X[i]),axis=1)) + pass + ####################################################################### + # END OF YOUR CODE # + ####################################################################### + return dists +``` + + + +# 不用循环 + +这个方法超级叼,将计算欧式距离完全转换为了矩阵运算 + +推导过程可以见[这一篇博文](https://blog.csdn.net/IT_forlearn/article/details/100022244),公式写的很详细。 + +具体代码如下: + +```python +def compute_distances_no_loops(self, X): + """ + Compute the distance between each test point in X and each training point + in self.X_train using no explicit loops. + Input / Output: Same as compute_distances_two_loops + """ + num_test = X.shape[0] + num_train = self.X_train.shape[0] + dists = np.zeros((num_test, num_train)) + ######################################################################### + # TODO: # + # Compute the l2 distance between all test points and all training # + # points without using any explicit loops, and store the result in # + # dists. # + # # + # You should implement this function using only basic array operations; # + # in particular you should not use functions from scipy. # + # # + # HINT: Try to formulate the l2 distance using matrix multiplication # + # and two broadcast sums. # + ######################################################################### + d1 = np.sum(np.square(X),axis=1,keepdims=True) #shape(num_test,1) + d2 = np.sum(np.square(self.X_train),axis=1) #shape(1,num_train) + d3 = np.multiply(np.dot(X,self.X_train.T),-2) #shape(num_test,num_train) + dists = np.sqrt(d1+d2+d3) + #numpy真的神奇,d1和d2居然可以加,一个是 num_test*1 一个是 1*num_train,两个一加变成了num_test*num_train 而每一行就是num_test里的每一个元素加上num_train的一行的那个元素 + pass + ######################################################################### + # END OF YOUR CODE # + ######################################################################### + return dists +``` + + + +# 三种方法的耗时比较 + +```python +# Let's compare how fast the implementations are +def time_function(f, *args): + """ + Call a function f with args and return the time (in seconds) that it took to execute. + """ + import time + tic = time.time() + f(*args) + toc = time.time() + return toc - tic + +two_loop_time = time_function(classifier.compute_distances_two_loops, X_test) +print('Two loop version took %f seconds' % two_loop_time) + +one_loop_time = time_function(classifier.compute_distances_one_loop, X_test) +print('One loop version took %f seconds' % one_loop_time) + +no_loop_time = time_function(classifier.compute_distances_no_loops, X_test) +print('No loop version took %f seconds' % no_loop_time) + +# you should see significantly faster performance with the fully vectorized implementation + +-------- +Two loop version took 37.711502 seconds +One loop version took 96.402259 seconds +No loop version took 0.381577 seconds +``` + + + + + +最后吐槽一下,hexo写博客插图片真的太蠢了,居然不能插相对路径,别人github都可以,搞得我在typora上写起来及其不友好,只好不插图片了唉。 \ No newline at end of file diff --git a/MLDLRL/ML/scatter_onehot.md b/MLDLRL/ML/scatter_onehot.md new file mode 100644 index 0000000..e9e9ea3 --- /dev/null +++ b/MLDLRL/ML/scatter_onehot.md @@ -0,0 +1,73 @@ +--- +title: 使用scatter函数来计算onehot向量 +date: 2022-04-06 +tags: [pytorch] +--- +> 日期:2022.4.6 + +## 循环太慢 + +最近在炼丹的过程中的测试阶段有一个需求,就是对于每一个输出的结果,需要对它做一个log,然后再把所有的结果加起来。因为本来的输入又是一个序列,加上batch_size,所以一共有三个维度,因此最开始我计算的时候是做了两重循环: + +```python +for i in range(batch_size): + for j in range(sequence_size): + ans += -torch.log2(pred[i][j][int(gt[i][j])]) +``` + +这里说明一下,pred,也就是预测输出的shape是`batch_size * sequence_size * classes`;gt,也就是实际label的shape是`batch_size * sequence_size`。 + +可以看到,这里的两重循环是非常辣眼睛的,经过测试,哪怕是将 `batch_size` 设置为512,`sequence_size` 设置128,一趟计算下来也要3秒钟,在数据量很大的时候,这是非常难以接受的。 + +## 矩阵运算替代 + +所以我就想到或许可以用矩阵运算的方式来代替掉这个循环,经过一阵推导,发现确实能够这样做,步骤大致如下: + +1. 将 gt 矩阵变为 one-hot 形式,也就是说把 `batch_size * sequence_size` 变为 `batch_size * sequence_size * classes`的形式,也就是把最后一个维度变为one-hot。 +2. 将pred矩阵与gt矩阵进行点乘。 +3. 将点乘后的结果中所有的0替换为1。 +4. 在整个矩阵上求 -log2 +5. 矩阵求和 + +### one-hot矩阵生成 + +所以可以看到,一阵操作里面最复杂的应该就是如何把一个矩阵的最后一维变为one-hot形式,经过一番搜索,我发现pytorch里提供了一个叫做 `scatter_` 的函数,其大致的描述是这样的: + +image-20220406203621393 + +用通俗易懂的话来说就是,给定一个source矩阵,和一个index矩阵,选定一个维度,根据index矩阵从source矩阵里**挑一些值**来作为结果(self)。 + +于是用下面两行代码就可以生产一个one-hot矩阵: + +```python +one_hot_gt = torch.zeros(batch_size,sequence_size,256) +one_hot_gt = one_hot_gt.scatter_(dim=2,index=gt.reshape(batch_size,sequence_size,1).data.long().cpu(),value=1) +``` + +解释一下: + +1. 首先,生成一个全零的`batch_size * sequence_size * classes`形式的矩阵。 +2. 对于`scatter_`函数的参数说明一下。`dim=2`是因为要在最后一维做one-hot。由于输入的index(gt)和self(one_hot_gt)需要维度一致,所以我们将它reshape一下,变为 `batch_size * sequence_size * 1`。`value=1`是因为要把所有的值都替换为1,这也就是source。 + +如此操作就可以达成以下的效果: + +```python +one_hot_gt[i][j][gt[i][j][0]] = 1 # gt[i][j][0]就是类别的编号 +``` + +### 后续操作 + +一旦生成了one-hot矩阵,后续的操作就很简单了,完整的代码如下: + +```python +# 首先将 gt 修改为one-hot形式 (batch_size, sequence_size, 256) +one_hot_gt = torch.zeros(batch_size,sequence_size,256) +one_hot_gt = one_hot_gt.scatter_(dim=2,index=gt.reshape(batch_size,sequence_size,1).data.long().cpu(),value=1) +# 将one_hot_gt和pred点乘 +ans = torch.mul(one_hot_gt.cuda(),pred.cuda()) +# 将所有的0替换为1 +ans[ans == 0] = 1 +ans = -torch.log2(ans).sum() +``` + +最后测试得到,原本要3s的运算只需要0.1s了~ \ No newline at end of file diff --git a/MLDLRL/README.md b/MLDLRL/README.md new file mode 100644 index 0000000..a22727d --- /dev/null +++ b/MLDLRL/README.md @@ -0,0 +1 @@ +相关的文章整理得不多,所以把ML/DL/RL这几个部分的文章整理在一起 \ No newline at end of file diff --git a/MLDLRL/RL/README.md b/MLDLRL/RL/README.md new file mode 100644 index 0000000..8cc3fa7 --- /dev/null +++ b/MLDLRL/RL/README.md @@ -0,0 +1,3 @@ +不会强化学习,但是也会学习一些相关的知识。 + +* [value&policy iteration](value&policy_iteration.md) \ No newline at end of file diff --git a/MLDLRL/RL/value&policy_iteration.md b/MLDLRL/RL/value&policy_iteration.md new file mode 100644 index 0000000..a796448 --- /dev/null +++ b/MLDLRL/RL/value&policy_iteration.md @@ -0,0 +1,89 @@ +--- +title: 值迭代和策略迭代 +date: 2021-11-27 +tags: [强化学习] +--- +最近的随机过程作业涉及到了value iteration和policy iteration,于是就去搜了相关的网课[CS229-lecture17](https://www.youtube.com/watch?v=d5gaWTo6kDM)来看,讲的很好,于是做了一些笔记。 +由于这些东西和强化学习能扯上一点关系,所以就分类到了强化学习。 + +### 关于π和$$V^{\pi}$$ + +π是optimal polocy,是一个state -> action 的映射 + +比如: + +![image-20211227183556201](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227183556201.png) + +$$V^{\pi}$$ 是 从某一个位置开始(当作初始位置) 所获得的一个 reward + +![image-20211227183739382](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227183739382.png) + + + + + +#### bellman equation + +![image-20211227184936151](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227184936151.png) + + + +#### 例子 + +比如说要求 (3,1) 这个位置的 $$V^{\pi}$$ + +![image-20211227185217516](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227185217516.png) + +那么求法就是: + +![image-20211227185300830](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227185300830.png) + +如果把每一个状态对应的$$V^{\pi}$$当作是一个未知数,由于已经给出了确定的policy $$\pi$$ (也就是说action是确定的),那么根据bellman equation,每一个$$V^{\pi}$$都可以写出一个方程。所以可以用一个linear solver来构建方程并且求解。 + + + + + +### 关于$$\pi^*$$ 和 $$V^*$$ + +$$\pi^*$$ 是 optimal policy + +$$ V^* $$是 optimal policy 对应的value function,公式如下: + +![image-20211227191645200](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227191645200.png) + +#### 对应的bellman equation + +![image-20211227192056237](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227192056237.png) + + + +### value iteration + +由于V(s)直接存在相互依赖关系,value iteration有两种,一种是 synchrounous,就是每一个V(s)同步更新,另一种就是asynchrounous,就是V(s)的更新不同步。不过都差不多。 + +![image-20211227193331244](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227193331244.png) + +经过多次迭代,V就会很快收敛到V*。最后收敛完毕后,再去计算对应的 π(s) 。 + + + +### policy iteration + +policy iteration的重点在于π,在循环中。 + +第一步:根据 π(s) 计算出对应的V。 如何计算? -> 如前面所述,使用linear solver来计算。 + +第二步:假设V是optimal value function,即V*,然后更新 π(s)。 + +![image-20211227195556271](https://raw.githubusercontent.com/ljhgpp/whatisthis/main/static/image-20211227195556271.png) + + + +### pros and cons of two methods + +由于policy iteration实际是一个 linear solver based 的 方法,所以当状态数很少时,求解速度很快,但是当状态数目很多时,求解速度就会变慢。所以在状态数目很多时,倾向于使用 value iteration。 + + + + diff --git a/README.md b/README.md new file mode 100644 index 0000000..83e0bad --- /dev/null +++ b/README.md @@ -0,0 +1,7 @@ +# jh的个人站 + +我的知识整理站点,什么都写一点。 + +> 从PaperMod白嫖过来的logo :D + +