图形流水线中光栅化原理与实现

光栅化主要解决的问题
光栅化原理
- 判断像素在三角形内部or外部
- 顶点属性插值
- - 重心坐标系
  - 插值深度
  - - 直接插值Z值的问题
    - 正确的深度插值公式推导
  - 透视校正
代码实现

光栅化主要解决的问题

在传统图形学流水线，技术难点可以分为两大类：

可见性(visibility)
着色(shading)
文章图形流水线中坐标变换详解：模型矩阵、视角矩阵、投影矩阵中的坐标系变换和本文介绍的光栅化就是解决 可见性(visibility) 的关键技术。

着色(shading) 涉及到后续局部光照模型以及全局光照模型。这类知识点在以后我学明白后再写文章介绍给大家。

光栅化原理

那什么是可见性问题呢？简单来说就是处在后面的物体应该被前面的物体所遮挡，从而在渲染结果上表现为不可见（对于不透明物体来说是这样，对于透明物体来说是颜色混合）。
结合文章（图形流水线中坐标变换详解：模型矩阵、视角矩阵、投影矩阵），所有三角形点在经过模型矩阵变换、视角矩阵变换、投影矩阵变化以及透视除法后，坐标都变换到NDC坐标系下（x,y,z∈【-1， 1】）。而在知道输出屏幕大小的情况下，通过窗口变换可将x/y变换到窗口大小下（x∈【0，width】 y∈【0，height】z不变）。至此我们即将所有三角形投射到raster_space中。
光栅化的作用是判断哪些像素点与三角形“干涉”(在三角形内部)，并插值出三角形内部点的属性值（Z值、颜色、法向、纹理坐标等）。

光栅化解决可见性(visibility) 思想：借助Z-Buffer, FrameBuffer
- Z-Buffer是一个与光栅图像大小一致的二维数组，记录着每个像素点的深度值。
- Z-Buffer中所有值初始化为inifity，当三角形上点P的Z值小于Z-Buffer中记录的值时，Z-Buffer中相应位置深度值更新为P点的Z值。否则不进行更新。
- FrameBuffer也是一个与光栅图像大小一致的二维数据，记录着每个像素点的颜色值
- 当Z-Buffer深度值更新时，FrameBuffer对应处颜色值也进行更新。
光栅化流程：

//遍历所有三角形
FOR 每个已经转换到窗口坐标系(raster space)下的三角形：
//内部两重循环遍历所有像素点
FOR row in raster space：
FOR column in raster space：
pixelCoord = vec2(row, column)
IF pixelCoord 在三角形内部：
插值pixelCoord 的Z值
用Z值与Z-Buffer做深度测试
IF 通过深度测试：
更新Z-Buffer对应处的Z值
插值pixelCoord的其他顶点属性(color, normal, uv)
将像素点处的颜色值写入FrameBuffer

根据光栅化步骤，我们可以提出两个关键点：

如何判断某个像素在三角形内部
如何正确插值顶点属性

判断像素在三角形内部or外部

前提条件：三角形都已投影到raster space空间中

该问题可以简化为，在二维平面中，如何判断一个点在三角形内部。如果你学过计算几何中凸包求解问题。你可能记得里面有一个TO-LEFT测试：一条向量无限延长后可以将平面分成两部分，分别称为其左侧和右侧

给定三角形旋向，三角形三条边即可以确定三条向量，如果对于一个点P，对这三条向量来说，它都处于同一侧(都在左侧(或者边上)、都在右侧(或者边上))，那个我们可以称点P它在三角形内部。具体情况如下图所示：
用向量叉乘判断点P在向量的左侧或右侧或在直线上

因此只需要判断三次叉乘结果是否同号(或者=0)即可判断点是否在上三角形上
判断像素在三角形内部or外部伪代码

//edge Funciton
bool edgeFunction(vec2 p, vec2 a, vec2 b, vec2 c){
//三角形三个点abc构成三个向量ab bc ca
vec2 ab = b - a;
vec2 bc = c - b;
vec2 ca = a - c;

//待判断点p分别于abc构成三个向量 ap, bp, cp
vec2 ap = p - a;
vec2 bp = p - b;
vec2 cp = p - c;
//得到三次叉乘结果
result1 = ab * ap;
result2 = bc * bp;
result3 = ca * cp;
//判断是否同号
float threshold = 1e-5;
if(result1 > -threshold && result2 > -threshold && result3 > -threshold)//都为非负数
return true;
if(result1 < threshold && result2 < threshold && result3 ><threshold)//都为非正数
return true;
return false;
}

顶点属性插值

流水线中我们只会给三角形顶点赋值某些顶点属性，光栅化需要根据三个顶点信息插值三角形内部点的顶点属性。

重心坐标系

对于三角形内部点P来说，它可以由顶点V0/V1/V2唯一表示：
p = λ0 * V0 + λ1 * V1 + λ2 * V2 且 λ0 >=0 λ1 >=0 λ2 >=0; λ0 + λ1 + λ2 = 1
(λ0, λ1, λ2)即V0V1V2组成的三角形内部点P的重心坐标

同理我们可以用重心坐标来插值顶点属性
P_attribute = λ0 * V0_attribute + λ1 * V1_attribute + λ2 * V2_attribute

既然可以用重心坐标系来插值所有顶点属性信息，那么又如何计算重心坐标系呢？

观察上图，你可以从中看出，重心坐标系值与V0V1P/V1V2P/V2V0P这三个三角形的面积有关。

由向量叉乘的集合意义可知: 叉乘得到的行列式的值，即两个向量所围成的平行四边形面积

结合判断点P是否在三角形内部的函数。我们可以在判断点P是否在三角形内部时，同时计算出P点的重心坐标系值。这是不是非常巧妙！

//改进 edge Funciton，判断点p是否在三角形abc上时，
//同时返回其重心坐标系值
bool edgeFunction(vec2 p, vec2 a, vec2 b, vec2 c, vector<< float>& barycentricCoord ){
//三角形三个点abc构成三个向量ab bc ca
vec2 ab = b - a;
vec2 bc = c - b;
vec2 ca = a - c;
//三角形面积的两倍
float triangleArea = abs(ab * bc);
//待判断点p分别于abc构成三个向量 ap, bp, cp
vec2 ap = p - a;
vec2 bp = p - b;
vec2 cp = p - c;
//得到三次叉乘结果
result1 = ab * ap;
result2 = bc * bp;
result3 = ca * cp;
//计算重心坐标系
barycentricCoord.push_back(abs(result2) / triangleArea);
barycentricCoord.push_back(abs(result3) / triangleArea);
barycentricCoord.push_back(abs(result1) / triangleArea);
//判断是否同号
float threshold = 1e-5;
if(result1 > -threshold && result2 > -threshold && result3 > -threshold)//都为非负数
return true;
if(result1 < threshold && result2 < threshold && result3 ><threshold)//都为非正数
return true;
return false;
}

插值深度

直接插值Z值的问题

根据上文结论可得，对于三角形内部点p，其Z值可以由此插值：
P_z = λ0 * V0_z + λ1 * V1_z + λ2 * V2_z
但实际上直接使用此式子是错误的，因为经过投影变换后，Z值不再满足线性变化。下图可以清楚的展现这种错误：

在投影前，P的Z值为:
P.z = V0.z * (1 - 0.666) + V1.z * 0.666 = 4.001;
投影后，P‘的Z值为
P’.z = V0.z * (1 - 0.8333) + V1.z * 0.8333 = 4.499;

正确的深度插值公式推导

因此正确的深度插值公式为：
1 / P_z = λ0 * 1 / V0_z + λ1 * 1 / V1_z + λ2 * 1 / V2_z

透视校正

同理在插值其他顶点属性（如颜色、法相、纹理坐标）时，如果直接用重心坐标系，也不正确。
因此正确的方法是借助已经正确插值的Z值，属性与Z值满足线性变化关系。

因此正确的顶点属性插值公式为：
Attr = Z * [Attr0 / V0_z * λ0 +Attr1 / V1_z * λ1 + Attr2 / V2_z * λ2 ]
下图展示了使用正确的属性插值与不正确的属性插值的颜色误差。

代码实现

vector2.h

//vector2.h
#pragma once
#include<iostream>using namespace std;template <typename T>
class Vector2 {public:T x, y;float z;public://Vector(){}~Vector2() {}Vector2(T xx = 0, T yy = 0, float zz = 0) :x(xx), y(yy), z(zz) {}Vector2(Vector2& t) { x = t.x; y = t.y; z = t.z; }//标量乘除Vector2 multiplayByScalar(float a, Vector2<T>& result) {result.x = x * a;result.y = y * a;return result;}Vector2 divideByScalar(float a, Vector2<T>& result) {result.x = x / a;result.y = y / a;return result;}//矢量运算Vector2<T> add(Vector2& t, Vector2<T>& result) {result.x = x + t.x;result.y = y + t.y;return result;}Vector2<T> divide(Vector2& t, Vector2<T>& result) {result.x = x - t.x;result.y = y - t.y;return result;}float cross(Vector2& t) {return y * t.x - x * t.y;}//模场float length() {return sqrt(x * x + y * y);}//归一化void normalize() {float l = this->length();if (l < 1e-5) {cout << "向量模长为0.0，不能归一化" << endl;return;}x /= l;y /= l;}//流运算friend ostream& operator<<(ostream& os, const Vector2<T>& t) {os << "(x, y, z):" << "(" << t.x << ", " << t.y << ", " << t.z << ")" << endl;return os;}};

Vector3.h

//Vector3
#pragma once
#include<iostream>using namespace std;template <typename T>
class Vector3 {public:T x, y, z;float w;public://Vector(){}~Vector3() {}Vector3(T xx = 0, T yy = 0, T zz = 0, float ww = 1.0) :x(xx), y(yy), z(zz), w(ww) {}Vector3(Vector3& t) { x = t.x; y = t.y; z = t.z; w = t.w; }//标量乘除Vector3 multiplayByScalar(float a, Vector3<T>& result) {result.x = x * a;result.y = y * a;result.z = z * a;return result;}Vector3 divideByScalar(float a, Vector3<T>& result) {result.x = x / a;result.y = y / a;result.z = z / a;return result;}//矢量运算Vector3<T> add(Vector3& t, Vector3<T>& result) {result.x = x + t.x;result.y = y + t.y;result.z = z + t.z;return result;}Vector3<T> divide(Vector3& t, Vector3<T>& result) {result.x = x - t.x;result.y = y - t.y;result.z = z - t.z;return result;}Vector3<T> cross(Vector3& t, Vector3<T>& result) {result.x = y * t.z - z * t.y;result.y = z * t.x - x * t.z;result.z = x * t.y - y * t.x;return result;}//比较运算bool equal(Vector3<T>& t) {float threshold = 1e-5;if (abs(x - t.x) < threshold && abs(y - t.y) < threshold && abs(z - t.z) < threshold)return true;return false;}//透视除法void perspectiveDivision() {x = x / T(w);y = y / T(w);z = z / T(w);w = 1.0f;}//模场float length() {return sqrt(x * x + y * y + z * z);}//归一化void normalize() {float l = this->length();if (l < 1e-5) {cout << "向量模长为0.0，不能归一化" << endl;return ;}x /= l;y /= l;z /= l;}//流运算friend ostream& operator<<(ostream& os, const Vector3<T>& t) {os << "(x, y, z, w):" << "(" << t.x << ", " << t.y << ", " << t.z << ", " << t.w << ")" << endl;return os;}};

Matrix4.h

//Matrix4.h
#pragma once
#include "Vector3.h"class Matrix4
{public:Matrix4();~Matrix4();Matrix4(float a0, float a1, float a2, float a3,float a4, float a5, float a6, float a7,float a8, float a9, float a10, float a11,float a12, float a13, float a14, float a15);Matrix4(Matrix4& t);//重置为单位矩阵void setIdentityMatrix();//求矩阵的行列式的值float getDeterminant();//求矩阵的逆矩阵Matrix4 invert();//矩阵与向量/点相乘Vector3<float> multiplyByVector(Vector3<float> v);//矩阵与矩阵相乘Matrix4 multiplyByMatrix4(Matrix4& m);//交换矩阵的两行void swapRow(int row1, int row2);//输出矩阵friend ostream& operator<<(ostream& os, Matrix4& t);public:float m[16];};

Matrix4.cpp

//Matrix4.cpp
#include "Matrix4.h"Matrix4::Matrix4()
{m[0] = m[5] = m[10] = m[15] = 1.0f;m[1] = m[2] = m[3] = m[4] = m[6] = m[7] = m[8] = m[9] = m[11] = m[12] = m[13] = m[14] = 0.0f;
}Matrix4::~Matrix4(){}Matrix4::Matrix4(float a0, float a1, float a2, float a3,float a4, float a5, float a6, float a7,float a8, float a9, float a10, float a11,float a12, float a13, float a14, float a15) {m[0] = a0; m[1] = a1; m[2] = a2; m[3] = a3;m[4] = a4; m[5] = a5; m[6] = a6; m[7] = a7;m[8] = a8; m[9] = a9; m[10] = a10; m[11] = a11;m[12] = a12; m[13] = a13; m[14] = a14; m[15] = a15;
}Matrix4::Matrix4(Matrix4& t) {for (int i = 0; i < 16; i++)m[i] = t.m[i];
}//重置为单位矩阵
void Matrix4::setIdentityMatrix() {m[0] = m[5] = m[10] = m[15] = 1.0f;m[1] = m[2] = m[3] = m[4] = m[6] = m[7] = m[8] = m[9] = m[11] = m[12] = m[13] = m[14] = 0.0f;
}//求矩阵的行列式的值
float Matrix4::getDeterminant() {return m[0] * m[5] * m[10] * m[15] + m[1] * m[6] * m[11] * m[12] + m[2] * m[7] * m[8] * m[13] + m[3] * m[4] * m[9] * m[14]- m[3] * m[6] * m[9] * m[12] - m[7] * m[10] * m[13] * m[0] - m[11] * m[14] * m[1] * m[4] - m[15] * m[2] * m[5] * m[8];
}//求矩阵的逆矩阵
/* C++实现矩阵求逆 采用  Gauss-Jordan elimination method
** 原理：利用行变换和新建一个单位阵，当前矩阵通过行变换变为单位阵时，单位阵通过同样的行变换会变成矩阵的逆
** 行变换：1.任意交换两行  2.一行加减除了这行的其他行  3.一行乘除一个数** 步骤：
** 1.通过交换行将对角线上所有值变为非0值
** 2.通过行变换按列从左到右将非对角线上的值变换为0
** 3.通过行缩放，将对角线上的值变为1
*/
Matrix4 Matrix4::invert() {Matrix4 invertM = Matrix4();int row = 0, column = 0;float tempM[16];for (int i = 0; i < 16; i++)tempM[i] = m[i];int matrixSize = 4;float threshold = 1e-5f;//1.通过交换行将对角线上所有值变为非0值for (column = 0; column < matrixSize; column++) {if (abs(tempM[column * matrixSize + column]) <= threshold) {//对角线上值为0 需要交换int swapRow = column;for (int i = 0; i < matrixSize; i++) {if (abs(tempM[i * matrixSize + column]) > threshold && abs(tempM[column * matrixSize + i]) > threshold)swapRow = i;}if (swapRow == column) {cout << "该矩阵没有逆矩阵" << endl;return invertM;}else {//交换for (int i = 0; i < matrixSize; i++) {float temp = tempM[column * matrixSize + i];tempM[column * matrixSize + i] = tempM[swapRow * matrixSize + i];tempM[swapRow * matrixSize + i] = temp;}invertM.swapRow(column, swapRow);}}//2.通过行变换按列从左到右将非对角线上的值变换为0//记录对角线的值float pivotsValue = tempM[column * matrixSize + column];for (row = 0; row < matrixSize; row++) {//如果处于对角线上 则跳过if (column == row)continue;//如果不处于对角线上 则将这一列上其他所有值变为0.0float coeff = tempM[row * matrixSize + column] / pivotsValue;//for (int i = 0; i < matrixSize; i++) {tempM[row * matrixSize + i] -= coeff * tempM[column * matrixSize + i];invertM.m[row * matrixSize + i] -= coeff * invertM.m[column * matrixSize + i];}tempM[row * matrixSize + column] = 0.0f;}}//3.通过行缩放，将对角线上的值变为1for (row = 0; row < matrixSize; row++) {float coeff = 1.0f / tempM[row * matrixSize + row];for (int i = 0; i < matrixSize; i++) {tempM[row * matrixSize + i] *= coeff;invertM.m[row * matrixSize + i] *= coeff;}}return invertM;
}//矩阵与向量/点相乘
/* 矩阵是行主序，右乘 this * v
**
*/Vector3<float> Matrix4::multiplyByVector(Vector3<float> v) {Vector3<float> reusltV;reusltV.x = float(m[0] * v.x + m[1] * v.y + m[2] * v.z + m[3] * v.w);reusltV.y = float(m[4] * v.x + m[5] * v.y + m[6] * v.z + m[7] * v.w);reusltV.z = float(m[8] * v.x + m[9] * v.y + m[10] * v.z + m[11] * v.w);reusltV.w = float(m[12] * v.x + m[13] * v.y + m[14] * v.z + m[15] * v.w);return reusltV;
}//矩阵与矩阵相乘 this * M
Matrix4 Matrix4::multiplyByMatrix4(Matrix4& M) {Matrix4 resultM;resultM.m[0] = m[0] * M.m[0] + m[1] * M.m[4] + m[2] * M.m[8] + m[3] * M.m[12];resultM.m[1] = m[0] * M.m[1] + m[1] * M.m[5] + m[2] * M.m[9] + m[3] * M.m[13];resultM.m[2] = m[0] * M.m[2] + m[1] * M.m[6] + m[2] * M.m[10] + m[3] * M.m[14];resultM.m[3] = m[0] * M.m[3] + m[1] * M.m[7] + m[2] * M.m[11] + m[3] * M.m[15];resultM.m[4] = m[4] * M.m[0] + m[5] * M.m[4] + m[6] * M.m[8] + m[7] * M.m[12];resultM.m[5] = m[4] * M.m[1] + m[5] * M.m[5] + m[6] * M.m[9] + m[7] * M.m[13];resultM.m[6] = m[4] * M.m[2] + m[5] * M.m[6] + m[6] * M.m[10] + m[7] * M.m[14];resultM.m[7] = m[4] * M.m[3] + m[5] * M.m[7] + m[6] * M.m[11] + m[7] * M.m[15];resultM.m[8] = m[8] * M.m[0] + m[9] * M.m[4] + m[10] * M.m[8] + m[11] * M.m[12];resultM.m[9] = m[8] * M.m[1] + m[9] * M.m[5] + m[10] * M.m[9] + m[11] * M.m[13];resultM.m[10] = m[8] * M.m[2] + m[9] * M.m[6] + m[10] * M.m[10] + m[11] * M.m[14];resultM.m[11] = m[8] * M.m[3] + m[9] * M.m[7] + m[10] * M.m[11] + m[11] * M.m[15];resultM.m[12] = m[12] * M.m[0] + m[13] * M.m[4] + m[14] * M.m[8] + m[15] * M.m[12];resultM.m[13] = m[12] * M.m[1] + m[13] * M.m[5] + m[14] * M.m[9] + m[15] * M.m[13];resultM.m[14] = m[12] * M.m[2] + m[13] * M.m[6] + m[14] * M.m[10] + m[15] * M.m[14];resultM.m[15] = m[12] * M.m[3] + m[13] * M.m[7] + m[14] * M.m[11] + m[15] * M.m[15];return resultM;
}//交换矩阵的两行
void Matrix4::swapRow(int row1, int row2) {int matrixSize = 4;for (int i = 0; i < matrixSize; i++) {float temp = m[row1 * matrixSize + i];m[row1 * matrixSize + i] = m[row2 * matrixSize + i];m[row2 * matrixSize + i] = temp;}
}//输出矩阵
ostream& operator<<(ostream& os, Matrix4& t) {os << "[" << t.m[0] << ",\t" << t.m[1] << ",\t" << t.m[2] << ",\t" << t.m[3] << "]" << endl;os << "[" << t.m[4] << ",\t" << t.m[5] << ",\t" << t.m[6] << ",\t" << t.m[7] << "]" << endl;os << "[" << t.m[8] << ",\t" << t.m[9] << ",\t" << t.m[10] << ",\t" << t.m[11] << "]" << endl;os << "[" << t.m[12] << ",\t" << t.m[13] << ",\t" << t.m[14] << ",\t" << t.m[15] << "]" << endl;return os;
}

Camera .h

//Camera .h
#pragma once#include "Vector3.h"
#include "Matrix4.h"
const float PI = 3.14159f;class Camera {public :Vector3<float> position;Vector3<float> direction;Vector3<float> up;float nearClippingPlane;float farClippingPlane;float horizonFov;float verticalFov;Matrix4 cameraCoord;public:Camera();~Camera();Camera(Vector3<float>& p);//设置direction 和 upvoid lookAt(Vector3<float> d, Vector3<float> u);//设置前后裁剪面void setClippingPlane(float n, float f);//设置hf 和 vfvoid setHFandVF(float hf, float vf);//计算视角矩阵函数群Matrix4 getCameraMatrix();Matrix4 getViewMatrix();//计算投影矩阵函数群//Matrix4 getProjectionMatrix_rasterSpace();Matrix4 getProjectionMatrix_NDCSpace();};

Camera.cpp

//Camera.cpp
#include "Camera.h"
Camera::Camera() {position = Vector3<float>(0.0f, 0.0f, 0.0f);direction = Vector3<float>(0.0f, 0.0f, -1.0f, 0.0f);up = Vector3<float>(0.0f, 1.0f, 0.0f, 0.0f);nearClippingPlane = 1.0f;farClippingPlane = 100.0f;horizonFov = verticalFov = 90.0f;
}Camera::~Camera() {}Camera::Camera(Vector3<float>& p) {position = Vector3<float>(p);direction = Vector3<float>(0.0f, 0.0f, -1.0f, 0.0f);up = Vector3<float>(0.0f, 1.0f, 0.0f, 0.0f);nearClippingPlane = 1.0f;farClippingPlane = 100.0f;horizonFov = verticalFov = 90.0f;
}//设置direction 和 up
void Camera::lookAt(Vector3<float> d, Vector3<float> u) {direction = Vector3<float>(d);up = Vector3<float>(u);
}//设置前后裁剪面
void Camera::setClippingPlane(float n, float f) {nearClippingPlane = n;farClippingPlane = f;
}//设置hf 和 vf
void  Camera::setHFandVF(float hf, float vf) {horizonFov = hf;verticalFov = vf;
}//计算相机坐标系
/* 得到相机坐标系基底的矩阵表示
** 视角坐标系Z方向与direction方向相反  W = -normalize(direction)
** 视角坐标系X方向 U = Up * W
** 视角坐标系Y方向 V = W * U
*/
Matrix4 Camera::getCameraMatrix() {Vector3<float> W(1.0, 1.0, 1.0, 0.0);direction.multiplayByScalar(-1.0, W);W.normalize();up.normalize();Vector3<float> U(1.0, 1.0, 1.0, 0.0);up.cross(W, U);U.normalize();Vector3<float> V(1.0, 1.0, 1.0, 0.0);W.cross(U, V);V.normalize();Matrix4 cameraCoord = Matrix4();/*[U.x. V.x, W.x, camera.x ][U.y, V.y, W.y, camera.y ][U.z, V.z, W.z, camera.z ][0.0, 0.0, 0.0, 1.0      ]*/cameraCoord.m[0] = U.x;cameraCoord.m[4] = U.y;cameraCoord.m[8] = U.z;cameraCoord.m[12] = 0.0;cameraCoord.m[1] = V.x;cameraCoord.m[5] = V.y;cameraCoord.m[9] = V.z;cameraCoord.m[13] = 0.0;cameraCoord.m[2] = W.x;cameraCoord.m[6] = W.y;cameraCoord.m[10] = W.z;cameraCoord.m[14] = 0.0;cameraCoord.m[3] = position.x;cameraCoord.m[7] = position.y;cameraCoord.m[11] = position.z;cameraCoord.m[15] = 1.0;this->cameraCoord = cameraCoord;return cameraCoord;
}//计算图形流水线中模型变化（世界坐标系--->视角坐标系）的 视角矩阵
Matrix4 Camera::getViewMatrix() {Matrix4 cameraCoord = this->getCameraMatrix();return cameraCoord.invert();
}计算图形流水线中模型变化（视角坐标系--->rasterSpace）的投影矩阵
//Matrix4 Camera::getProjectionMatrix_rasterSpace() {//
//}//计算图形流水线中模型变化（视角坐标系--->NDCSpace）的投影矩阵
Matrix4 Camera::getProjectionMatrix_NDCSpace() {/*** 投影矩阵的具体推导请看 https://blog.csdn.net/qq_27161673** [ 2N / (r - l), 0           ,  (r + l)/(r - l),  0 ]** [ 0           , 2N / (t - b),  (t + b)/(t - b),  0 ]** [ 0           , 0           ,  (N + F)/(N - F),  2NF / (N - F)   ]** [ 0           , 0           ,  -1             ,  0               ]*///角度转弧度float horizonFov_radian = horizonFov * PI / 180.0f;float verticalFov_radian = verticalFov * PI / 180.0f;/*** r = N * tan(horizonFov_radian / 2);** l = -r;** t = N * tan(verticalFov_radian / 2);** t = -b;*/float r = nearClippingPlane * tan(horizonFov_radian / 2.0f);float l = -r;float t = nearClippingPlane * tan(verticalFov_radian / 2.0f);float b = -t;Matrix4 projectionMatrix = Matrix4();projectionMatrix.m[0] = nearClippingPlane / r;projectionMatrix.m[1] = 0.0f;projectionMatrix.m[2] = 0.0f;projectionMatrix.m[3] = 0.0f;projectionMatrix.m[4] = 0.0f;projectionMatrix.m[5] = nearClippingPlane / t;projectionMatrix.m[6] = 0.0f;projectionMatrix.m[7] = 0.0f;projectionMatrix.m[8] = 0.0f;projectionMatrix.m[9] = 0.0f;projectionMatrix.m[10] = (nearClippingPlane + farClippingPlane) / (nearClippingPlane - farClippingPlane);projectionMatrix.m[11] = 2.0f * nearClippingPlane * farClippingPlane / (nearClippingPlane - farClippingPlane);projectionMatrix.m[12] = 0.0f;projectionMatrix.m[13] = 0.0f;projectionMatrix.m[14] = -1.0f;projectionMatrix.m[15] = 0.0f;return projectionMatrix;
}

main.cpp

/************************************************************************/
/*2019-12-29 自己的光栅化程序                                           */
/************************************************************************//*
** 分成两个阶段
** 一、实现流水线中坐标变化过程：包括模型矩阵、视角矩阵以及投影矩阵。
** 二、实现光栅化：包括隐藏面消除和顶点属性插值
**
** 本程序分为以下几个类
**
*/#include "Camera.h"
#include "Vector2.h"
#include <vector>
#include <fstream>bool edgeFunction(Vector2<float> p, Vector2<float>a, Vector2<float>b, Vector2<float>c, vector<float>& barycentricCoord);int main() {/************************************************************************//* 设定窗口大小                                                                     *//************************************************************************/int screen_width = 400;int screen_height = 400;/************************************************************************//* 设定相机                                                                      *//************************************************************************/// 设定相机位置Vector3<float> cameraPosition = Vector3<float>(0.0f, 0.0f, 3.0f, 1.0f);Camera camera = Camera(cameraPosition);//设定相聚 看向方向direction 和相机上方向upcamera.lookAt(Vector3<float>(0.0f, 0.0f, -1.0f, 0.0f), Vector3<float>(0.0f, 1.0f, 0.0f, 0.0f));//设定视域前/后裁剪面camera.setClippingPlane(1.0f, 10.0f);//设定水平/垂直视域角度camera.setHFandVF(90.0f, 90.0f);/************************************************************************//* 给定一个三角形数据（默认世界坐标系）                                                                     *//************************************************************************/Vector3<float> A_world = Vector3<float>(-2.0f, -1.0f, 0.0f, 1.0f);Vector3<float> B_world = Vector3<float>(2.0f, -1.0f, 1.0f, 1.0f);Vector3<float> C_world = Vector3<float>(-1.0f, 1.0f, -3.0f, 1.0f);cout << "A_world" << A_world << endl;cout << "B_world" << B_world << endl;cout << "C_world" << C_world << endl;Vector3<float> D_world = Vector3<float>(-2.0f, -1.0f, 1.0f, 1.0f);Vector3<float> E_world = Vector3<float>(2.0f, -1.0f, 0.0f, 1.0f);Vector3<float> F_world = Vector3<float>(1.0f, 1.0f, -3.0f, 1.0f);/************************************************************************//* 转换到视角坐标系： p_view = viewMatrix * p_world                 *//************************************************************************/Matrix4 viewMatrix = camera.getViewMatrix();Vector3<float> A_view = viewMatrix.multiplyByVector(A_world);Vector3<float> B_view = viewMatrix.multiplyByVector(B_world);Vector3<float> C_view = viewMatrix.multiplyByVector(C_world);cout << "A_view" << A_view << endl;cout << "B_view" << B_view << endl;cout << "C_view" << C_view << endl;Vector3<float> D_view = viewMatrix.multiplyByVector(D_world);Vector3<float> E_view = viewMatrix.multiplyByVector(E_world);Vector3<float> F_view = viewMatrix.multiplyByVector(F_world);/************************************************************************//* 转换到NDC坐标系 ：p_pro = projectionMatrix * p_view                                              *//************************************************************************/Matrix4 projectionMatrix = camera.getProjectionMatrix_NDCSpace();Vector3<float> A_pro = projectionMatrix.multiplyByVector(A_view);Vector3<float> B_pro = projectionMatrix.multiplyByVector(B_view);Vector3<float> C_pro = projectionMatrix.multiplyByVector(C_view);cout << "透视除法前 A_pro" << A_pro << endl;cout << "透视除法前 B_pro" << B_pro << endl;cout << "透视除法前 C_pro" << C_pro << endl;//透视除法A_pro.perspectiveDivision();B_pro.perspectiveDivision();C_pro.perspectiveDivision();cout << "A_pro" << A_pro << endl;cout << "B_pro" << B_pro << endl;cout << "C_pro" << C_pro << endl;Vector3<float> D_pro = projectionMatrix.multiplyByVector(D_view);Vector3<float> E_pro = projectionMatrix.multiplyByVector(E_view);Vector3<float> F_pro = projectionMatrix.multiplyByVector(F_view);D_pro.perspectiveDivision();E_pro.perspectiveDivision();F_pro.perspectiveDivision();/************************************************************************//* 转换到视图坐标系：p_screen = width/height * p_pro                                                                     *//************************************************************************/float half_screen_width = screen_width / 2;float half_screen_height = screen_height / 2;Vector2<float> A_screen = Vector2<float>(float((A_pro.x + 1.0f) * half_screen_width), float((1.0f - A_pro.y) * half_screen_height), A_view.z);Vector2<float> B_screen = Vector2<float>(float((B_pro.x + 1.0f) * half_screen_width), float((1.0f - B_pro.y) * half_screen_height), B_view.z);Vector2<float> C_screen = Vector2<float>(float((C_pro.x + 1.0f) * half_screen_width), float((1.0f - C_pro.y) * half_screen_height), C_view.z);cout << "A_screen" << A_screen << endl;cout << "B_screen" << B_screen << endl;cout << "C_screen" << C_screen << endl;Vector2<float> D_screen = Vector2<float>(float((D_pro.x + 1.0f) * half_screen_width), float((1.0f - D_pro.y) * half_screen_height), D_view.z);Vector2<float> E_screen = Vector2<float>(float((E_pro.x + 1.0f) * half_screen_width), float((1.0f - E_pro.y) * half_screen_height), E_view.z);Vector2<float> F_screen = Vector2<float>(float((F_pro.x + 1.0f) * half_screen_width), float((1.0f - F_pro.y) * half_screen_height), F_view.z);//点数组vector<Vector2<float>> vertice = vector<Vector2<float>>(6);vertice[0] = (A_screen);vertice[1] = (B_screen);vertice[2] = (C_screen);vertice[3] = (D_screen);vertice[4] = (E_screen);vertice[5] = (F_screen);//索引数组vector<unsigned int> indices = vector<unsigned int>();indices.push_back(0);indices.push_back(1);indices.push_back(2);indices.push_back(3);indices.push_back(4);indices.push_back(5);//点颜色数组Vector3<unsigned char> red = Vector3<unsigned char>(255, 0,0);Vector3<unsigned char> green = Vector3<unsigned char>(255, 255, 0);Vector3<unsigned char> blue = Vector3<unsigned char>(255, 0, 255);Vector3<unsigned char> red1 = Vector3<unsigned char>(0, 0, 255);Vector3<unsigned char> green1 = Vector3<unsigned char>(0, 255, 255);Vector3<unsigned char> blue1 = Vector3<unsigned char>(255, 0, 255);vector<Vector3<unsigned char>> colors = vector<Vector3<unsigned char>>(6);colors[0] = (red);colors[1] = (green);colors[2] = (blue);colors[3] = (red1);colors[4] = (green1);colors[5] = (blue1);/************************************************************************//* 光栅化** 1.判断每个像素点是否在screen space空间的三角形中** 2.如果在 计算其三角形重心坐标，并进行顶点属性插值** 3.将颜色信息写入frame buffer中** *//************************************************************************///创建Z-Buffer 和 FrameBuffer//创建Z-Buffer 记录FrameBuffer对应像素点的深度值int bufferSize = screen_width * screen_height;float* zBuffer = new float[bufferSize];//创建FrameBuffervector<vector<float>> framebuffer(bufferSize); //Z-Buffer初始化为最大值  FrameBuffer初始化为黑色for (int i = 0; i < bufferSize; i++) {zBuffer[i] = FLT_MAX;framebuffer[i] = vector<float>(4);for (int j = 0; j < 4; j++)framebuffer[i][j] = 0.0f;}//遍历所有像素for (int row = 0; row < screen_height; row++) {for (int column = 0; column < screen_width; column++) {//当前像素点坐标Vector2<float> pixelCoord = Vector2<float>(float(column) + 0.5f, float(row) + 0.5f);//遍历所有三角形int triangleNum = indices.size() / 3;for (int triangleIndex = 0; triangleIndex < triangleNum ; triangleIndex++) {//三角形三个点坐标（screen space坐标系下） 默认逆时针方向Vector2<float> a = vertice[indices[triangleIndex * 3]];Vector2<float> b = vertice[indices[triangleIndex * 3 + 1]];Vector2<float> c = vertice[indices[triangleIndex * 3 + 2]];//三角形三个点的颜色Vector3<unsigned char> a_color = colors[indices[triangleIndex * 3]];Vector3<unsigned char> b_color = colors[indices[triangleIndex * 3 + 1]];Vector3<unsigned char> c_color = colors[indices[triangleIndex * 3 + 2]];//记录像素点的重心坐标vector<float> barycentricCoord = vector<float>();//edge Functionif (edgeFunction(pixelCoord, a, b, c, barycentricCoord) == false)continue;else {//插值像素点的Z坐标// 1 / z = λ * 1/ z0 + (1-λ) * 1 /z1pixelCoord.z = 1.0f / (barycentricCoord[0] / a.z + barycentricCoord[1] / b.z + barycentricCoord[2] / c.z);//深度测试 只有Z坐标小于0的物体才能被看到if (pixelCoord.z < 1e-5 && abs(pixelCoord.z) < zBuffer[row * screen_width + column] && (-pixelCoord.z) >= camera.nearClippingPlane && (-pixelCoord.z) <= camera.farClippingPlane) {zBuffer[row * screen_width + column] = abs(pixelCoord.z);//未使用透视矫正/*float R = (barycentricCoord[0] * a_color.x + barycentricCoord[1] * b_color.x  + barycentricCoord[2] * c_color.x );float G = (barycentricCoord[0] * a_color.y + barycentricCoord[1] * b_color.y  + barycentricCoord[2] * c_color.y );float B = (barycentricCoord[0] * a_color.z + barycentricCoord[1] * b_color.z  + barycentricCoord[2] * c_color.z );*///插值像素点处的顶点属性（目前只有颜色)  使用透视矫正// attribute = z * [attribute0 / z0 * λ + attribute1 / z1 * (1-λ)]float R = pixelCoord.z * (barycentricCoord[0] * float(a_color.x) / a.z + barycentricCoord[1] * float(b_color.x) / b.z + barycentricCoord[2] * float(c_color.x) / c.z);float G = pixelCoord.z * (barycentricCoord[0] * float(a_color.y) / a.z + barycentricCoord[1] * float(b_color.y) / b.z + barycentricCoord[2] * float(c_color.y )/ c.z);float B = pixelCoord.z * (barycentricCoord[0] * float(a_color.z) / a.z + barycentricCoord[1] * float(b_color.z) / b.z + barycentricCoord[2] * float(c_color.z) / c.z);//更新frameBuffer中颜色值vector<float> piexlColor = vector<float>(4);piexlColor[0] = R / 255.0f; piexlColor[1] = G / 255.0f; piexlColor[2] = B / 255.0f; piexlColor[3] = 1.0f;framebuffer[row * screen_width + column] = piexlColor;}}}}}//将framebuffer中的数据写入图片中std::ofstream ofs;ofs.open("./output.ppm");ofs << "P3\n" << screen_width << ' ' << screen_height << "\n255\n";for (int i = 0; i < screen_height; i++){for (int j = 0; j < screen_width; j++){float r = framebuffer[i * screen_width + j][0];float g = framebuffer[i * screen_width + j][1];float b = framebuffer[i * screen_width + j][2];int ir = int(255.99 * r);int ig = int(255.99 * g);int ib = int(255.99 * b);ofs << ir << ' ' << ig << ' ' << ib << '\n';}}ofs.close();delete[] zBuffer;return 0;
}/************************************************************************/
/* 判断一个点是否在三角形内部                                           */
/************************************************************************/
bool edgeFunction(Vector2<float> p, Vector2<float>a, Vector2<float>b, Vector2<float>c, vector<float>& barycentricCoord) {//根据旋向和三个点坐标组成三个向量Vector2<float> ab = Vector2<float>(b.x - a.x, b.y - a.y);Vector2<float> bc = Vector2<float>(c.x - b.x, c.y - b.y);Vector2<float> ca = Vector2<float>(a.x - c.x, a.y - c.y);Vector2<float> ap = Vector2<float>(p.x - a.x, p.y - a.y);Vector2<float> bp = Vector2<float>(p.x - b.x, p.y - b.y);Vector2<float> cp = Vector2<float>(p.x - c.x, p.y - c.y);float result1 = ab.cross(ap);float result2 = bc.cross(bp);float result3 = ca.cross(cp);float yuzhi = .1f;if ((result1 > -yuzhi && result2 > -yuzhi && result3 > -yuzhi) || (result1 < yuzhi && result2 < yuzhi && result3 < yuzhi) ) {//计算重心坐标系float tirangle_area = abs(ab.cross(bc));/*cout << tirangle_area << endl;*/barycentricCoord.push_back(abs(result2 / tirangle_area));barycentricCoord.push_back(abs(result3 / tirangle_area));barycentricCoord.push_back(abs(result1 / tirangle_area));return true;}return false;}