手搓一个小型渲染器（五）

相机移动

在之前的渲染中，我们已经使用了Gouraud Shading，其想法非常简单，模型的制作者给出了模型中每个顶点的法向量，之后我们计算每个顶点的光照强度，然后在三角形中进行插值。

在欧几里得空间中，坐标可以有一个原点和一个基$(O,\vec{i},\vec{j},\vec{k})$组成，假设点$P$的坐标为$(x,y,z)$，则$\overrightarrow{OP}$可以用下面的式子表示：
$$\overrightarrow{OP}=\vec{i}x+\vec{j}y+\vec{k}z=\left[\begin{array}{ccc}\vec{i}& \vec{j}& \vec{k}\end{array}\right]\left[\begin{array}{c}x\\ y\\ z\end{array}\right]$$
现在假设我们有另一个基于点的基$(O^{\prime },{\vec{i}}^{\prime },{\vec{j}}^{\prime },{\vec{k}}^{\prime })$，那么我们如何从之前的基进行变换得到现有的基呢？那么存在一个矩阵$M$,可以得到最终的变换：
$$\left[\begin{array}{ccc}\vec{i}& \vec{j}& \vec{k}\end{array}\right]=\left[\begin{array}{ccc}\overrightarrow{i^{\prime }}& \overrightarrow{j^{\prime }}& \overrightarrow{k^{\prime }}\end{array}\right]\times M$$
其具体过程如下图所示：

那么对于$\overrightarrow{OP}$，我们可以得到：
$$\overrightarrow{OP}=\overrightarrow{O{O}^{\prime }}+\overrightarrow{O’P}=\left[\begin{array}{ccc}\vec{i}& \vec{j}& \vec{k}\end{array}\right]\left[\begin{array}{c}{O}_x^{\prime }\\ O_y^{\prime }\\ O_z^{\prime }\end{array}\right]+\left[\begin{array}{ccc}\overrightarrow{i^{\prime }}& \overrightarrow{j^{\prime }}& \overrightarrow{k^{\prime }}\end{array}\right]\left[\begin{array}{c}{x}^{\prime }\\ y^{\prime }\\ z^{\prime }\end{array}\right]$$
我们可以将上面的式子进行结合
$$\overrightarrow{OP}=\overrightarrow{O{O}^{\prime }}+\overrightarrow{O’P}=\left[\begin{array}{ccc}\vec{i}& \vec{j}& \vec{k}\end{array}\right]\left(\left[\begin{array}{c}{O}_x^{\prime }\\ O_y^{\prime }\\ O_z^{\prime }\end{array}\right]+M\left[\begin{array}{c}{x}^{\prime }\\ y^{\prime }\\ z^{\prime }\end{array}\right]\right)$$
得到下面的式子
$$\left[\begin{array}{c}x\\ y\\ z\end{array}\right]=\left[\begin{array}{c}{O}_x^{\prime }\\ O_y^{\prime }\\ O_z^{\prime }\end{array}\right]+M\left[\begin{array}{c}{x}^{\prime }\\ y^{\prime }\\ z^{\prime }\end{array}\right]⟶\left[\begin{array}{c}{x}^{\prime }\\ y^{\prime }\\ z^{\prime }\end{array}\right]=M^{-1}\left(\left[\begin{array}{c}x\\ y\\ z\end{array}\right]-\left[\begin{array}{c}{O}_x^{\prime }\\ O_y^{\prime }\\ O_z^{\prime }\end{array}\right]\right)$$

创建我们自己的gluLookAt

OpenGL和我们的小型渲染器只能在摄像机位于z轴上的情况下绘制场景。如果我们想移动相机，我们可以移动所有场景，让相机静止不动。我们可以这样考虑我们的问题，我们想将相机放置于点$e$（眼睛）绘制场景，相机应该指向点$c$（中心），这样给定的向量u（向上）在最终渲染中是垂直的。

这意味着我们要在$(c,x^{\prime },y^{\prime },z^{\prime })$中进行渲染，但是我们的模型的基是在$(O,x,y,z)$中给出的，我们所需要的只是计算坐标的变换。以下是一个C++代码，用于计算必要的4x4矩阵ModelView：

mat<4,4> lookat(vec3 eye, vec3 center, vec3 up) {
    vec3 z = (eye - center).normalized();
    vec3 x = cross(up, z).normalized();
    vec3 y = cross(x, z).normalized();
    mat<4, 4> res;
    res = res.identity();
    for (int i = 0; i < 3; i++) {
        res[0][i] = x[i];
        res[1][i] = y[i];
        res[2][i] = z[i];
        res[i][3] = -center[i];
    }
    return res;
}

z’是由向量ce给出的（不要忘记对其进行归一化，稍后会有所帮助）。我们如何计算x’？简单地通过u和z之间的叉积。然后我们计算y’，使其与已经计算的x’和z’正交（让我提醒您，在我们的问题设置中，ce和u不一定正交）。最后一步是将原点转换到观察者e的点，我们的变换矩阵就准备好了。现在，只要在模型帧中获得坐标为（x，y，z，1）的任何点，将其乘以矩阵ModelView，就可以获得相机帧中的坐标！顺便说一句，ModelView这个名称来自OpenGL术语。

ViewPort

在之前我们写过下面这段代码：

screen_coords[j] = Vec2i((v.x+1.)*width/2., (v.y+1.)*height/2.);

它的主要目的是将一个属于[-1.1] * [-1,1]范围内的方形平面中的点映射到(width, height)这个平面，值（v.x+1）在0和2之间变化，（v.x+1）/2在0和1之间变化，并且（v.x+1）*width/2扫描所有图像。因此，我们有效地将双单位正方形映射到图像上。现在我们可以改变一下上面的代码：

mat<4,4> viewport(int x, int y, int w, int h) {
    mat<4, 4> m; 
    m = m.identity();
    m[0][3] = x + w / 2.f;
    m[1][3] = y + h / 2.f;
    m[2][3] = depth / 2.f;

    m[0][0] = w / 2.f;
    m[1][1] = h / 2.f;
    m[2][2] = depth / 2.f
    return m;
}

这意味着双单位立方体[-1,1][-1,1][-1,1]被映射到屏幕立方体[x，x+w][y，y+h][0，d]上。立方体而不是矩形，这是因为使用z缓冲区进行深度计算。这里d是z缓冲区的分辨率。我喜欢它等于255，因为转储z缓冲区的黑白图像以进行调试很简单。

链式坐标变换

所以，让我们总结一下。我们的模型（例如角色）是在它们自己的局部框架（对象坐标）中创建的。它们被插入到以世界坐标表示的场景中。从一个到另一个的转换是用矩阵模型进行的。然后，我们想在相机帧（眼睛坐标）中表达它，这种变换被称为View。然后，我们使用投影矩阵对场景进行变形以创建透视变形（第4课），该矩阵将场景转换为所谓的剪辑坐标。最后，我们绘制场景，将剪辑坐标转换为屏幕坐标的矩阵称为Viewport。

同样，如果我们从.obj文件中读取一个点v，然后在屏幕上绘制它，它将经历以下转换链：

Viewport * Projection * View * Model * v.

我们在具体实现时，使用了下面的代码：

Vec3f v = model->vert(face[j]);
screen_coords[j] =  Vec3f(ViewPort*Projection*ModelView*Matrix(v));

当我们只渲染一个对象时，矩阵Model等于单位矩阵，因此将其与矩阵View合并。

法向量变换

• 如果我们有一个模型，它的法向量是由艺术家给出的，并且这个模型用仿射映射变换，那么法向量要用映射变换，等于原始映射矩阵的逆矩阵的换位

对于上面这句话，我们可以具体举一个例子说明，假设2D平面上存在一个三角形，三角形三个点坐标为（0，0），（1，0），（0，1），我们可以得到一个向量n（1，1），垂直于斜边，现在我们对这三个坐标进行仿射变换，将y坐标拉伸两倍，如果我们也将向量n的y坐标拉伸，那么最终的结果向量n将不再垂直于斜边。

对于法向量的处理，我们不能对原先的法向量进行简单的变换，因为他们可能不再是正常的，因此我们需要计算变换后模型的新的法向量。回到3D空间中，我们现在有一个法向量n=(A,B,C)，这个法向量穿过原点，因此我们有下式：
$$Ax+By+Cz=0$$

现在我们将n的坐标扩展为齐次坐标，可以得到下面这个式子：
$$\left[\begin{array}{cccc}A& B& C& 0\end{array}\right]\times \left[\begin{array}{c}x\\ y\\ z\\ 1\end{array}\right]=0$$
其中由于n为向量，因此用0填充，而xyz为点则用1来填充。

让我们在两者之间插入一个单位矩阵（M乘以M的倒数等于单位）：
$$\left(\left[\begin{array}{cccc}A& B& C& 0\end{array}\right]\times M^{-1}\right)\times \left(M\times \left[\begin{array}{c}x\\ y\\ z\\ 1\end{array}\right]\right)=0$$
右括号中的表达式-表示对象的变换点。在左边-表示变换对象的法线向量。在标准约定中，我们通常将坐标写为列，因此我们可以将前面的表达式重写如下
$${\left((M^T{)}^{-1}\times \left[\begin{array}{c}A\\ B\\ C\\ 0\end{array}\right]\right)}^T\times \left(M\times \left[\begin{array}{c}x\\ y\\ z\\ 1\end{array}\right]\right)=0$$
左括号告诉我们，通过应用仿射映射的逆转置矩阵，可以从旧法线计算变换对象的法线。如果我们的变换矩阵M是均匀缩放、旋转和平移的组合，那么M等于它的逆转置，因为在这种情况下，逆转置和转置是相互抵消的。但由于我们的矩阵包括透视变形，通常这种技巧没有帮助。

main.cpp

#include <vector>
#include <iostream>
#include <cmath>
#include <limits>
#include "tgaimage.h"
#include "model.h"
#include "geometry.h"

const int width = 800;
const int height = 800;
const int depth = 1024;

double* zbuffer = NULL;
vec3 light_dir = {0, 1, 0};
//vec3 light_dir = light_dir.normalized();
vec3 eye = { 3, 3, 3 };
vec3 center{ 0, 0, 0 };

mat<4,4> viewport(int x, int y, int w, int h) {
    mat<4, 4> m; 
    m = m.identity();
    m[0][3] = x + w / 2.f;
    m[1][3] = y + h / 2.f;
    m[2][3] = depth / 2.f;

    m[0][0] = w / 2.f;
    m[1][1] = h / 2.f;
    m[2][2] = depth / 2.f;
    return m;
}

mat<4,4> lookat(vec3 eye, vec3 center, vec3 up) {
    vec3 z = (eye - center).normalized();
    vec3 x = cross(up, z).normalized();
    vec3 y = cross(x, z).normalized();
    mat<4, 4> res;
    res = res.identity();
    for (int i = 0; i < 3; i++) {
        res[0][i] = x[i];
        res[1][i] = y[i];
        res[2][i] = z[i];
        res[i][3] = -center[i];
    }
    return res;
}

void triangle(vec3 t0, vec3 t1, vec3 t2, float ity0, float ity1, float ity2, TGAImage& image, double* zbuffer) {
    if (t0.y == t1.y && t0.y == t2.y) return; // i dont care about degenerate triangles
    if (t0.y > t1.y) { std::swap(t0, t1); std::swap(ity0, ity1); }
    if (t0.y > t2.y) { std::swap(t0, t2); std::swap(ity0, ity2); }
    if (t1.y > t2.y) { std::swap(t1, t2); std::swap(ity1, ity2); }

    int total_height = t2.y - t0.y;
    for (int i = 0; i < total_height; i++) {
        bool second_half = i > t1.y - t0.y || t1.y == t0.y;
        int segment_height = second_half ? t2.y - t1.y : t1.y - t0.y;
        float alpha = (float)i / total_height;
        float beta = (float)(i - (second_half ? t1.y - t0.y : 0)) / segment_height; // be careful: with above conditions no division by zero here
        vec3 A = t0 + vec3(t2 - t0) * alpha;
        vec3 B = second_half ? t1 + vec3(t2 - t1) * beta : t0 + vec3(t1 - t0) * beta;
        float ityA = ity0 + (ity2 - ity0) * alpha;
        float ityB = second_half ? ity1 + (ity2 - ity1) * beta : ity0 + (ity1 - ity0) * beta;
        if (A.x > B.x) { std::swap(A, B); std::swap(ityA, ityB); }
        for (int j = A.x; j <= B.x; j++) {
            float phi = B.x == A.x ? 1. : (float)(j - A.x) / (B.x - A.x);
            vec3    P = vec3(A) + vec3(B - A) * phi;
            float ityP = ityA + (ityB - ityA) * phi;
            int idx = P.x + P.y * width;
            if (P.x >= width || P.y >= height || P.x < 0 || P.y < 0) continue;
            if (zbuffer[idx] < P.z) {
                zbuffer[idx] = P.z;
                TGAColor color{ (std::uint8_t)255 * ityP, (std::uint8_t)255 * ityP, (std::uint8_t)255 * ityP, (std::uint8_t)255 * ityP };
                image.set(P.x, P.y, color );
            }
        }
    }
}


vec3 barycentric(vec3 A, vec3 B, vec3 C, vec3 P) {
    vec3 s[2];
    for (int i = 2; i--; ) {
        s[i][0] = C[i] - A[i];
        s[i][1] = B[i] - A[i];
        s[i][2] = A[i] - P[i];
    }
    vec3 u = cross(s[0], s[1]);
    if (std::abs(u[2]) > 1e-2) // dont forget that u[2] is integer. If it is zero then triangle ABC is degenerate
        return vec3{ 1.f - (u.x + u.y) / u.z, u.y / u.z, u.x / u.z };
    return vec3{ -1, 1, 1 }; // in this case generate negative coordinates, it will be thrown away by the rasterizator
}

/*
*   带z-buffer的三角形渲染算法
*/
void triangle(vec3* pts, vec2* puv, vec3* origin, Model* model, double* zbuffer, TGAImage& image)
{
    //确定三角形的包围框
    vec2 bboxmin{ std::numeric_limits<double>::max(), std::numeric_limits<double>::max() };
    vec2 bboxmax{ std::numeric_limits<double>::min(), std::numeric_limits<double>::min() };
    vec2 bboxclamp{ ((double)image.width() - 1.), ((double)image.height() - 1.) };

    //遍历三角形顶点，确定四边形边框
    for (int i = 0; i < 3; i++)
    {
        bboxmin.x = std::max(0., std::min(bboxmin.x, pts[i].x));
        bboxmin.y = std::max(0., std::min(bboxmin.y, pts[i].y));
        bboxmax.x = std::min(bboxclamp.x, std::max(bboxmax.x, pts[i].x));
        bboxmax.y = std::min(bboxclamp.y, std::max(bboxmax.y, pts[i].y));
    }

    //遍历框内像素，确定像素点是否在三角形内部（使用重心坐标判断）
    vec3 p{};
    for (p.x = bboxmin.x; p.x <= bboxmax.x; p.x++) {
        for (p.y = bboxmin.y; p.y <= bboxmax.y; p.y++) {
            
            vec3 bc_screen = barycentric(pts[0], pts[1], pts[2], p);
            if (bc_screen.x < 0 || bc_screen.y < 0 || bc_screen.z < 0) continue;
            //利用重心坐标计算三角形内点p的z坐标值
            p.z = 0;

            //for (int i = 0; i < 3; i++) p.z += pts[i][2] * bc_screen[i];
            p.z = 1. / (bc_screen[0] / pts[0].z + bc_screen[1] / pts[1].z + bc_screen[2] / pts[2].z);

            vec2 color_coords{};
            for (int i = 0; i < 3; i++) {
                color_coords.x += puv[i].x * (double)model->diffuse().width() *  bc_screen[i] / pts[i].z;
                color_coords.y += puv[i].y * (double)model->diffuse().height() * bc_screen[i] / pts[i].z;
            }
            //for (int i = 0; i < 3; i++) {
            //    color_coords.x += puv[i].x * (double)model->diffuse().width() * bc_screen[i];
            //    color_coords.y += puv[i].y * (double)model->diffuse().height() * bc_screen[i];
            //}
            TGAColor color = model->diffuse().get((int)(color_coords.x * p.z), (int)(color_coords.y * p.z));

            //更新z-buffer,并且只有更新的时候，需要渲染更近的像素颜色值
            //double到int会有精度损失，因此不能涵盖到所有zbuffer
            if (zbuffer[(int)p.x + (int)p.y * image.width()] <= p.z) {
                zbuffer[(int)p.x + (int)p.y * image.width()] = p.z;
                image.set((int)p.x, (int)p.y, color);
            }
        }
    }
}

mat<4, 1> v2m(vec3 coords) {
    mat<4, 1> matrix{};
    matrix[0][0] = coords.x;
    matrix[1][0] = coords.y;
    matrix[2][0] = coords.z;
    matrix[3][0] = 1;
    return matrix;
}

vec3 m2v(mat<4, 1> matrix) {
    return vec3{ matrix[0][0] / matrix[3][0], matrix[1][0] / matrix[3][0], matrix[2][0] / matrix[3][0] };
}

int main() {

    Model* model = new Model("obj/african_head/african_head.obj");
    

    zbuffer = new double[width * height];
    for (int i = 0; i < width * height; i++) {
        zbuffer[i] = std::numeric_limits<double>::min();
    }

    { // draw the model
        mat<4, 4> ModelView = lookat(eye, center, vec3{ 0, 1, 0 });
        mat<4, 4> Projection;
        Projection = Projection.identity();
        mat<4, 4> ViewPort = viewport(width / 8, height / 8, width * 3 / 4, height * 3 / 4);
        Projection[3][2] = -1.f / (eye - center).norm();

        std::cerr << ModelView << std::endl;
        std::cerr << Projection << std::endl;
        std::cerr << ViewPort << std::endl;
        mat<4, 4> z = ViewPort * Projection * ModelView;
        std::cerr << z << std::endl;

        TGAImage image(width, height, TGAImage::RGB);
        //第几个三角形
        for (int i = 0; i < model->nfaces(); i++) {
            vec3 world_coords[3]{};
            vec3 origin_coords[3]{};
            vec3 screen_coords[3]{};
            vec2 uv_coords[3]{};
            float intensity[3]{};
            for (int j = 0; j < 3; j++) {
                world_coords[j] = model->vert(i, j);
                uv_coords[j] = model->uv(i, j);
                screen_coords[j] = m2v(ViewPort * Projection * ModelView * v2m(world_coords[j]));
                origin_coords[j] = m2v(Projection * ModelView * v2m(world_coords[j]));
                //intensity[j] = model->normal(i, j) * light_dir;
            }
            //triangle(screen_coords[0], screen_coords[1], screen_coords[2], intensity[0], intensity[1], intensity[2], image, zbuffer);
            triangle(screen_coords, uv_coords, origin_coords, model, zbuffer, image);
        }
        image.flip_vertically(); // i want to have the origin at the left bottom corner of the image
        image.write_tga_file("output.tga");
    }

    delete model;
    delete[] zbuffer;
    return 0;
}