【笔记】ejoy2d —— matrix

【笔记】ejoy2d —— matrix

2014 年 1 月 25 日

ejoy2d 是 ejoy 公司开源的 2d 图像引擎,主要作者是业界牛人云风.

矩阵理论基础

在学习这一段代码之前, 先引入一些矩阵的基础理论.

  • 2d图形的变换, 可以用一个6元组 {a, b, c, d, tx, ty} 来表示, 其中, a 和 d 表示缩放, b 和 c 表示旋转, tx 和 ty 表示位移, 这样就能得到一个基础的2d矩阵公式:
![](../data/2014-01-25-ejoy2d-matrix/2.gif)
  • 为了把二维图形的变化统一在一个坐标系里, 引入齐次坐标的概念, 即把一个图形用一个三维矩阵表示, 其中第三列总是 (0, 0, 1), 用来作为坐标系的标准, 这样, 变换6元组就可以用一个 3*3 的矩阵来表示:
这样, 我们就可以得到一个矩阵变换公式:

![](../data/2014-01-25-ejoy2d-matrix/4.gif)
  • 平移: 当

时, 则退化为:

![](../data/2014-01-25-ejoy2d-matrix/6.gif)
  • 缩放: 当

时, 则退化为:

![](../data/2014-01-25-ejoy2d-matrix/8.gif)
  • 旋转: 当

时, 则退化为:

![](../data/2014-01-25-ejoy2d-matrix/10.gif)
  • 归一化: 当

时, 则 M 退化为归一化矩阵:

![](../data/2014-01-25-ejoy2d-matrix/12.gif)
  • 逆矩阵: 所谓逆矩阵, 就是对于矩阵 A, 存在矩阵 B, A B = B A = E, 其中, E 是归一化矩阵, 那么 A 和 B 就互为各自的逆矩阵. 对于 M 而言, 它的逆矩阵可以计算得到:
![](../data/2014-01-25-ejoy2d-matrix/14.gif)

源码分析

简单回顾一下理论, 下面就可以开始 ejoy2d 的矩阵源码分析, lib/matrix.c 和 lib/matrix.h 封装了矩阵的操作, lib/lmatrix.c 和 lib/lmatrix.h是导出给lua的接口.

  • 六元组表示的矩阵(matrix) /* * | m[0] m[1] 0 | * | m[2] m[3] 0 | * | m[4] m[5] 1 | */ struct matrix { int m[6]; };
  • 矩阵归一化(matrix identity) matrix_identity(struct matrix *mm) { int *mat = mm->m; mat[0] = 1024; mat[1] = 0; mat[2] = 0; mat[3] = 1024; mat[4] = 0; mat[5] = 0; } 为了避免浮点运算, 这里以 1024 为基处理, 来保持整数运算, 实际上, 现在的变换矩阵 M 演变为:
  • 矩阵的乘法(matrix multiply) static inline void matrix_mul(struct matrix *mm, const struct matrix *mm1, const struct matrix *mm2) { int *m = mm->m; const int *m1 = mm1->m; const int *m2 = mm2->m; m[0] = (m1[0] * m2[0] + m1[1] * m2[2]) /1024; m[1] = (m1[0] * m2[1] + m1[1] * m2[3]) /1024; m[2] = (m1[2] * m2[0] + m1[3] * m2[2]) /1024; m[3] = (m1[2] * m2[1] + m1[3] * m2[3]) /1024; m[4] = (m1[4] * m2[0] + m1[5] * m2[2]) /1024 + m2[4]; m[5] = (m1[4] * m2[1] + m1[5] * m2[3]) /1024 + m2[5]; } 这里就是 3*3 的矩阵乘法, 唯一做了变化的就是在乘完之后除 1024 来保持scale.
  • 矩阵的逆运算(matrix inverse) int matrix_inverse(const struct matrix *mm, struct matrix *mo) { const int *m = mm->m; int *o = mo->m; if (m[1] == 0 && m[2] == 0) { return _inverse_scale(m,o); } if (m[0] == 0 && m[3] == 0) { return _inverse_rot(m,o); } int t = m[0] * m[3] - m[1] * m[2] ; if (t == 0) return 1; o[0] = (int32_t)((int64_t)m[3] * (1024 * 1024) / t); o[1] = (int32_t)(- (int64_t)m[1] * (1024 * 1024) / t); o[2] = (int32_t)(- (int64_t)m[2] * (1024 * 1024) / t); o[3] = (int32_t)((int64_t)m[0] * (1024 * 1024) / t); o[4] = - (m[4] * o[0] + m[5] * o[2]) / 1024; o[5] = - (m[4] * o[1] + m[5] * o[3]) / 1024; return 0; } 这里有几点需要注意的:
    1. 因为归一化时做了 1024 的线性变换, 所以这里有 1024 的参数.
    2. 当 b = c = 0 时, 退化为线性变换_inverse_scale, 这里做了简化处理提高效率.
    3. 当 a = d = 0 时, 退化为旋转变换_inverse_rot, 同样简化处理.
  • 矩阵的缩放(matrix scale) static inline void scale_mat(int *m, int sx, int sy) { if (sx != 1024) { m[0] = m[0] * sx / 1024; m[2] = m[2] * sx / 1024; m[4] = m[4] * sx / 1024; } if (sy != 1024) { m[1] = m[1] * sy / 1024; m[3] = m[3] * sy / 1024; m[5] = m[5] * sy / 1024; } } 之前说到, 当 M 退化为只有 a 和 d 时, 可以看做一个缩放(scale). 这里的 a 对应就是代码中的 sx, b 对应的就是代码中的 sy. 此时演进成如下的形式:
  • 矩阵的旋转(matrix rotate) static inline void rot_mat(int *m, int d) { if (d==0) return; int cosd = icosd(d); int sind = isind(d); int m0_cosd = m[0] * cosd; int m0_sind = m[0] * sind; int m1_cosd = m[1] * cosd; int m1_sind = m[1] * sind; int m2_cosd = m[2] * cosd; int m2_sind = m[2] * sind; int m3_cosd = m[3] * cosd; int m3_sind = m[3] * sind; int m4_cosd = m[4] * cosd; int m4_sind = m[4] * sind; int m5_cosd = m[5] * cosd; int m5_sind = m[5] * sind; m[0] = (m0_cosd - m1_sind) /1024; m[1] = (m0_sind + m1_cosd) /1024; m[2] = (m2_cosd - m3_sind) /1024; m[3] = (m2_sind + m3_cosd) /1024; m[4] = (m4_cosd - m5_sind) /1024; m[5] = (m4_sind + m5_cosd) /1024; } 在这种情况下, M 退化为旋转矩阵, 此处代码的原理也很容易理解:
代码中的 d 是经过变换的角度, 这个在lua的接口中有体现, 源码截取如下:

``` C
double r = luaL_checknumber(L,2);
matrix_rot(m, r * (1024.0 / 360.0));
```

代码中对角度的计算做了一个 cos 表, cos 和 sin 的计算都是通过查表来得到的, 代码中做 sin 计算时的 64, 是 ![](../data/2014-01-25-ejoy2d-matrix/18.gif) 经过 1024 换算后得到的.
``` C
static inline int
icost(int dd) {
    static int t[256] = {
    ...
    };
    if (dd < 0) {
        dd = 256 - (-dd % 256);
    } else {
        dd %= 256;
    }

    return t[dd];
}

static inline int
icosd(int d) {
    int dd = d/4;
    return icost(dd);
}

static inline int
isind(int d) {
    int dd = 64 - d/4;
    return icost(dd);
}
```
  • 矩阵的SRT变换 ejoy2d.matrix中用一个 srt(scale, rotate, translate) 结构体来封装了变换矩阵, 对矩阵的 srt 操作依次就是: scale_mat, rot_mat 以及 translate (直接增加 M 中的 tx 和 ty). struct srt { int offx; int offy; int scalex; int scaley; int rot; }; void matrix_srt(struct matrix *mm, const struct srt *srt) { scale_mat(mm->m, srt->scalex, srt->scaley); rot_mat(mm->m, srt->rot); mm->m[4] += srt->offx; mm->m[5] += srt->offy; } void matrix_sr(struct matrix *mat, int sx, int sy, int d) { int *m = mat->m; int cosd = icosd(d); int sind = isind(d); int m0_cosd = sx * cosd; int m0_sind = sx * sind; int m3_cosd = sy * cosd; int m3_sind = sy * sind; m[0] = m0_cosd /1024; m[1] = m0_sind /1024; m[2] = -m3_sind /1024; m[3] = m3_cosd /1024; } void matrix_rs(struct matrix *mat, int sx, int sy, int d) { int *m = mat->m; int cosd = icosd(d); int sind = isind(d); int m0_cosd = sx * cosd; int m0_sind = sx * sind; int m3_cosd = sy * cosd; int m3_sind = sy * sind; m[0] = m0_cosd /1024; m[1] = m3_sind /1024; m[2] = -m0_sind /1024; m[3] = m3_cosd /1024; } 其中函数 matrix_sr() 是先 scale 再 rotate, 没有 transform 的简化演进, b = c = 0:
与之相对的, 函数 matrix_rs() 是先 rotate 再 scale, 没有 transform 的演进:

![](../data/2014-01-25-ejoy2d-matrix/20.gif)
  • 矩阵的 Lua 接口 int ejoy2d_matrix(lua_State *L) { luaL_Reg l[] = { { "new", lnew }, { "scale", lscale }, { "trans", ltrans }, { "rot", lrot }, { "sr", lsr }, { "rs", lrs }, { "inverse", linverse }, { "mul", lmul }, { "tostring", ltostring }, { "identity", lidentity}, { "export", lexport }, { "import", limport }, { NULL, NULL }, }; luaL_newlib(L,l); return 1; } 这其中要注意的一点是, 如果 lnew() 没有输入参数的话, 则返回归一矩阵(scale 1024).

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

发表于

我来说两句

0 条评论
登录 后参与评论

相关文章

来自专栏机器之心

教程 | 利用TensorFlow和神经网络来处理文本分类问题

3397
来自专栏编程软文

前端慌不慌?用深度学习自动生成HTML代码

5066
来自专栏小鹏的专栏

trick—Batch Normalization

深度学习中 Batch Normalization为什么效果好? 这里分五部分简单解释一下Batch Normalization (BN)。 1. What ...

3128
来自专栏大数据挖掘DT机器学习

【Keras】完整实现‘交通标志’分类、‘票据’分类两个项目,让你掌握深度学习图像分类

我们一般用深度学习做图片分类的入门教材都是MNIST或者CIFAR-10,因为数据都是别人准备好的,有的甚至是一个函数就把所有数据都load进来了,所以跑起来...

5205
来自专栏Soul Joy Hub

自动微分(Automatic Differentiation)简介

http://blog.csdn.net/aws3217150/article/details/70214422 现代深度学习系统中(比如MXNet, Tens...

5243
来自专栏机器学习从入门到成神

XGboost数据比赛实战之调参篇(完整流程)

这一篇博客的内容是在上一篇博客Scikit中的特征选择,XGboost进行回归预测,模型优化的实战的基础上进行调参优化的,所以在阅读本篇博客之前,请先移步看一下...

2.9K8
来自专栏mathor

matlab—影像分析进阶

在这一章里面我们要做的事情全部都围绕两个问题,一个图像当中有多少个xxx,他们的大小是多少,举个例子

2002
来自专栏人工智能LeadAI

基于Spark /Tensorflow使用CNN处理NLP的尝试

01 前言 关于CNN如何和NLP结合,其实是被这篇文章(http://www.wildml.com/2015/11/understanding-convolu...

3826
来自专栏Gaussic

使用TensorFlow训练循环神经网络语言模型

读了将近一个下午的TensorFlow Recurrent Neural Network教程,翻看其在PTB上的实现,感觉晦涩难懂,因此参考了部分代码,自己写了...

2123
来自专栏人工智能LeadAI

用CNN做句子分类:CNN Sentence Classification (with Theano code)

01 Intro 本篇文章来细说CNN在NLP中的一大应用————句子分类。通过Yoon Kim的论文介绍一个应用,分析代码,并重构代码。 重构后的代码放在gi...

6665

扫码关注云+社区