钾肥喵的窝

我在 CODING 部署的 Hexo 博客

0%

[GAMES101 作业1]旋转与投影详解

发表于 更新于 分类于

任务清单

  1. 实现三维空间中绕 \(z\) 轴旋转, 构造变换矩阵
  2. 构造透视投影变换矩阵
  3. (提高作业)绕任意过原点的轴旋转, 构造变换矩阵

任务详解

任务1

对应到代码中是实现函数:

1
Eigen::Matrix4f get_model_matrix(float rotation_angle);

需要注意的是函数参数中的的旋转角度是角度制而非弧度制, 需要手动转换后再使用.

绕 \(z\) 轴旋转还是很简单的, 套公式就行了:

\[ \mathbf{R}_{z}(\alpha) = \begin{bmatrix} \cos{\alpha} & -\sin{\alpha} & 0 & 0 \\ \sin{\alpha} & \cos{\alpha} & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix} \]

代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
Eigen::Matrix4f get_model_matrix(float rotation_angle)
{
    Eigen::Matrix4f model = Eigen::Matrix4f::Identity();

    float r = rotation_angle / 180.0 * acos(-1);
    Eigen::Matrix4f translate;
    translate <<    cos(r),    -sin(r), 0.0,    0.0,
                    sin(r),    cos(r),  0.0,    0.0,
                    0.0,       0.0,     1.0,    0.0,
                    0.0,       0.0,     0.0,    1.0;

    model = translate * model;

    return model;
}

任务2

对应到代码中是实现函数:

1
Eigen::Matrix4f get_projection_matrix(float eye_fov, float aspect_ratio, float zNear, float zFar);

参数分别是: 垂直方向的视场角, 长宽比, 近处距离, 远处距离.

因为是透视投影, 我们默认中心点是 \((0, 0, z)\), 因此省略掉平移操作, 变换矩阵也能简化.

用 \(n, f, h, w\) 分别表示近处距离, 远处距离, 视场高度, 视场宽度, 有变换矩阵:

\[ \begin{bmatrix} \frac{2n}{w} & 0 & 0 & 0 \\ 0 & \frac{2n}{h} & 0 & 0 \\ 0 & 0 & \frac{n+f}{f-n} & \frac{-2nf}{f-n} \\ 0 & 0 & 1 & 0 \end{bmatrix}, \qquad 其中\begin{cases} n = \mathrm{zNear}, \\ f = \mathrm{zFar}, \\ h = 2n\cdot\tan{\frac{\mathrm{eye\_fov}}{2}}, \\ w = h\cdot\mathrm{aspect\_ratio}, \end{cases} \]

需要注意直接使用上述矩阵进行变换得到的是一个倒三角, 这是因为在右手系中我们通常是向 \(-z\) 方向看, 考虑到这点, 我们需要对符号进行一些处理, 最终的变换矩阵如下(符号定义同上):

\[ \begin{bmatrix} \frac{2n}{w} & 0 & 0 & 0 \\ 0 & \frac{2n}{h} & 0 & 0 \\ 0 & 0 & \frac{n+f}{n-f} & \frac{-2nf}{n-f} \\ 0 & 0 & -1 & 0 \end{bmatrix} \]

更正: 此处推导存在问题, 已于 作业2 中修正并给出具体过程.

代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
Eigen::Matrix4f get_projection_matrix(float eye_fov, float aspect_ratio,
    float zNear, float zFar)
{
    // Students will implement this function

    Eigen::Matrix4f projection = Eigen::Matrix4f::Identity();

    // TODO: Implement this function
    // Create the projection matrix for the given parameters.
    // Then return it.
    float half_fov = eye_fov / 360 * acos(-1);

    float n = zNear;
    float f = zFar;
    float h = 2 * n * tan(half_fov);
    float w = h * aspect_ratio;

    projection <<   -2*n/w, 0,      0,              0,
                    0,      -2*n/h, 0,              0,
                    0,      0,      -(n+f)/(n-f),   -2*n*f/(n-f),
                    0,      0,      1,              0;

    return projection;
}

任务3

实现函数:

1
Eigen::Matrix4f get_rotation(Vector3f axis, float angle);

直接套用Rodrigues' rotation formula就可以了, 简单推导可见附录.

代码(注意要把矩阵扩展成齐次坐标形式)如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
Eigen::Matrix4f get_rotation(Vector3f axis, float angle) 
{
    float r = angle / 180.0 * acos(-1);

    Eigen::Matrix3f i;

    i << 0, -axis.z(), axis.y(),
        axis.z(), 0, -axis.x(),
        -axis.y(), axis.x(), 0;

    Eigen::Matrix3f rotation = cos(r) * Eigen::Matrix3f::Identity() 
        + (1 - cos(r)) * axis * axis.transpose()
        + sin(r) * i;

    Eigen::Matrix4f translate = Eigen::Matrix4f::Identity();

    translate.block<3, 3>(0, 0) = rotation;

    return translate;
}

至于怎么显示在图形中, 对代码略作修改就可以了, 在绕 \(z\) 轴旋转的变换矩阵之前乘上这个矩阵就可以了. 

1
r.set_model(get_rotation(axis, angle2) * get_model_matrix(angle));

附录

Rodrigues' rotation formula的简单推导

符号约定: 用 \(\boldsymbol{k}, \boldsymbol{v}, \boldsymbol{v'}, \alpha\) 表示旋转轴(单位向量), 初始向量, 旋转后的向量, 旋转角(逆时针).

特殊情况

首先考虑 \(\boldsymbol{k}\perp\boldsymbol{v}\) 的特殊情况. 参考绕坐标轴旋转的推导, 我们同样可以将 \(\boldsymbol{v'}\) 分解到两个相互垂直的方向上.

首先想到的就是其中一个方向可以用 \(\boldsymbol{v}\) 表示, 那么另一个方向如何找呢, 很简单, \(\boldsymbol{k} \times \boldsymbol{v}\) 就是要找的另一个方向(注意叉乘的顺序, 这样叉乘得到的方向符合右手系, 便于后续操作).

由此就得到了特殊情况下的旋转公式:

\[ \boldsymbol{v'} = \cos{\alpha}\cdot\boldsymbol{v} + \sin{\alpha}\cdot(\boldsymbol{k} \times \boldsymbol{v}) \]

一般情况

一般情况可以由特殊情况推广得到, 我们可以拆分 \(\boldsymbol{v}\) 得到 \(\boldsymbol{v_{1}}\perp\boldsymbol{k}\) 以及 \(\boldsymbol{v_{2}}\parallel\boldsymbol{k}\). 显然, 在绕 \(\boldsymbol{k}\) 旋转的过程中, \(\boldsymbol{v_{2}}\) 分量不受影响, \(\boldsymbol{v_{1}}\) 变为 \(\boldsymbol{v_{1}'}\). 因此公式为:

\[ \boldsymbol{v'} = \cos{\alpha}\cdot\boldsymbol{v_{1}} + \sin{\alpha}\cdot(\boldsymbol{k} \times \boldsymbol{v_{1}}) + \boldsymbol{v_{2}} \]

从而将问题转换为如何求 \(\boldsymbol{v}\) 在 \(\boldsymbol{k}\) 及其法平面上的投影:

\[ \begin{cases} \boldsymbol{v_{1}} = \boldsymbol{v} - \boldsymbol{v_{2}} \\ \boldsymbol{v_{2}} = (\boldsymbol{k}\cdot\boldsymbol{v})\boldsymbol{k} \end{cases} \]

代入易得:

\[ \boldsymbol{v'} = \cos{\alpha}\cdot\boldsymbol{v} + \sin{\alpha}\cdot(\boldsymbol{k} \times \boldsymbol{v}) + (1-\cos{\alpha})(\boldsymbol{k}\cdot\boldsymbol{v})\boldsymbol{k} \]

矩阵形式

\[ \boldsymbol{v'} = \mathbf{R}\boldsymbol{v} \]

点乘的处理

即处理 \((1-\cos{\alpha})(\boldsymbol{k}\cdot\boldsymbol{v})\boldsymbol{k}\), 标量系数可以先不管, 同时做一下变形: \(\boldsymbol{k}(\boldsymbol{k}^{T}\cdot\boldsymbol{v})\), 把它们当作三个矩阵相乘, 用结合律很容易得到 \((\boldsymbol{k}\cdot\boldsymbol{k}^{T})\boldsymbol{v}\)

叉乘的处理

写出坐标表示就很容易得到变换矩阵了.

\[ \begin{align*} \boldsymbol{c} &= \boldsymbol{a} \times \boldsymbol{b} \\ &= \begin{bmatrix} a_{y}b_{z}-a_{z}b_{y} \\ a_{z}b_{x}-a_{x}b_{z} \\ a_{x}b_{y}-a_{y}b_{x} \end{bmatrix} \\ &= \begin{bmatrix} 0 &-a_{z} &a_{y} \\ a_{z} &0 &-a_{x} \\ -a_{y} &a_{x} &0 \end{bmatrix} \begin{bmatrix} b_{x} \\ b_{y} \\ b_{z} \end{bmatrix} \\ &= \mathbf{R_{a}}b \end{align*} \]

变换矩阵

综上, 我们很容易得到变换矩阵:

\[ \mathbf{R} = \cos{\alpha}\cdot\mathbf{I} + (1-\cos{\alpha})(\boldsymbol{k}\cdot\boldsymbol{k^{T}}) + \sin{\alpha}\cdot\begin{bmatrix} 0 &-k_{z} &k_{y} \\ k_{z} &0 &-k_{x} \\ -k_{y} &k_{x} &0 \end{bmatrix} \]