记一次从 AviSynth 到 VapourSynth 的迁移（3）：找不到的 Bug 与代码优化

回忆重连

离写这个系列的上一篇 Blog 已经一个月有余，离发布 VapourSynth 版 AreaResize 也三个星期了。一直没有继续写 Blog，有些东西都要忘了。这个时候发现 commit 表情还是挺有用的。

从功能上讲，上一篇 Blog 写到（基于原始 avs 版）修正了颜色问题，也在一定程度上支持高位深输入；但目标尺寸仍必须是“常规”的尺寸，比如 960x540 这样的，在一些稍“不常规”的尺寸上（如 1280x900），8bit 会分屏，32bit 会奔溃。

总之问题集中在尺寸上，尝试好久，包括调换垂直方向与水平方向的处理顺序，以及其他一些细节的尝试，都没有帮助。当时为了找这个 bug，盯着屏幕盯了好像有两三天，也没有结果。

最后没办法， 因为 AreaResize 除了原始的 avs 版，还有一个优化过的 avs 版。优化版作者的本意的加快处理速度，但我把代码迁移到 vs 后，发现上面的 Bug 居然解决了——虽然我真的看不出，两份代码除了速度优化外，其他方面的区别。

神奇的优化

如前面所说，我并没有发现原始 avs 版代码和优化 avs 版代码有除了速度优化之外的区别，但真的迁移了优化版代码后，上面的 Bug 就解决了。

具体的分析，等我看明白了再写吧。代码改到这里，排查了手滑写错了一个循环，除了 700x700 这样奇怪格式外，对于 YUV 格式，在高位深、一般目标尺寸输入上都没有问题了。下面就是支持 RGB 的问题。

在 1 月 2 号晚上和 1 月 3 号一天，都在解决支持 RGB 的问题，其实也就是在学相关的概念和数据结构。如同我在某个 commit 下面的评论中所写：“让我开始阅读 VapourSynth 的源码和官方 plugins 示例（就是 vscore.cpp 中的那些函数）。让我对 vs api 与框架有了近乎全面的认识，至少知道从头到尾是怎么计算的。”其中，数据结构的相关内容写在了这篇 Blog 中。

支持 RGB

写一写支持 RGB 的代码实现。

对于 YUV 格式，写一个循环，依次处理三个平面/通道即可。但我这么套在 RGB 上，不行，崩溃（但其实奔溃还有数据类型的问题，一开始的代码有整型处理，遇到浮点数就崩了）。让后又试着按 stack 格式写了代码，不循环，一个像素点依次处理三个通道，再移动到下一个像素点。像下面这样

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
ps = 3;  // plane_size，用于下面乘3，三个通道一起处理

template <typename T>
static void ResizeVerticalRGB() noexcept
{
    for (int curPixel = 0; curPixel < target_width * ps; curPixel++)
    {
        const T* curSrcpR = srcR + curPixel * ps;
        const T* curSrcpG = srcG + curPixel * ps + 1;
        const T* curSrcpB = srcB + curPixel * ps + 2;
        T* curDstpR = srcR + curPixel * ps;
        T* curDstpG = srcG + curPixel * ps + 1;
        T* curDstpB = srcB + curPixel * ps + 2;
        
        ...
    }
}

话说这改后，好像也是崩溃（或者是在关闭预览窗口后奔溃），所以放弃了，老老实实地（在process()函数中）先转换为 Interleaved 结构，在送到具体的Ver()中处理。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
if (d->vi->format->colorFamily == cmYUV)
{
    ...
}
else
{
    const T* srcpR = reinterpret_cast<const T*>(vsapi->getReadPtr(src, 0));
    const T* srcpG = reinterpret_cast<const T*>(vsapi->getReadPtr(src, 1));
    const T* srcpB = reinterpret_cast<const T*>(vsapi->getReadPtr(src, 2));
    
    int src_stride = vsapi->getStride(src, 0) / sizeof(T);
    int dst_stride = vsapi->getStirde(dst, 0) / sizeof(T);
    
    T* srcInterleaved = new (std::nothrow) T[d->vi->width * d->vi->height * 3];
    T* dstInterleaved = new (std::nothrow) T[d->target    * d->target     * 3];
    
    // 将待处理的数据转换为Interleaved结构
    for (int y = 0; y < src_height; y++)
    {
        for (int x = 0; x < src_width; x++)
        {
            const unsigned pos = (x + y * src_width) * 3; // 这里要乘3，因为是转成Interleaved
            // 这里的RGB顺序其实无所谓，只要自洽即可，但为了显得专业，按BGR的顺序
            srcInterleaved[pos] = srcpB[x];
            srcInterleaved[pos + 1] = srcpG[x];
            srcInterleaved[pos + 2] = srcpR[x];
        }
        srcpB += src_stride; // src_stride ≈ src_width，完成整个帧的遍历
        srcpG += src_stride;
        srcpR += src_stride;
    }
    
    // 核心处理
    Func(dst, src, (const T*)srcInterleaved, dstInterleaved, src_stride, dst_stride, d, vsapi);
    
    T* VS_RESTRICT dstpR = reinterpret_cast<T*>(vsapi->getWritePtr(dst, 0));
    T* VS_RESTRICT dstpG = reinterpret_cast<T*>(vsapi->getWritePtr(dst, 1));
    T* VS_RESTRICT dstpB = reinterpret_cast<T*>(vsapi->getWritePtr(dst, 2));
    
    // 将处理后的数据转换为Planar结构
    for (int y = 0; y < target_height; y++)
    {
        for (int x = 0; x < target_width; x++)
        {
            const unsigned pos = (x + y * target_width) * 3;
            dstpB[x] = dstInterleaved[pos];
            dstpG[x] = dstInterleaved[pos + 1];
            dstpR[x] = dstInterleaved[pos + 2];
        }
        dstpB += dst_stride;
        dstpG += dst_stride;
        dstpR += dst_stride;
    }  // 此时，dst已经是Planar结构
    
    delete[] srcInterleaved;
    delete[] dstInterleaved;
}

至此，数据结构的问题已经解决，但对于 RGB 输入仍不能正常处理。

然后修正了一个蛋疼的 Bug：我在写peak的时候，把1 << bitsPerSample写成了1 << sampleType，所以一直有问题、崩溃。改正之后，8bit RGB 输入终于能看了。而 16bit 和 32bit 由于查表的范围问题，还会崩溃。

然后将gamma_LUT声明为 double。到这里，输出结果是图像中的一部分带有蓝点红点绿点（提高gamma_LUT的精度，在单个帧中，会让这些点变少，而正常图像的部分增多）。

顺带，说一个 C++ 的基础知识，对于auto声明，应该是确定了后就不能改了，而不能在不同变量类型间反复横跳。

修正查表（LUT）的相关问题

查表，即 Look up table，在代码中一般缩写为 LUT。

论数组与指针的区别

1 月 12 号解决了上述（8bit）下图像带有红点绿点的问题。解决方法是把gamma_LUT的声明由数组变成指针。（开发版 commit：a5dc85d）

1
2
double gamma_LUT[RGB_PIXEL_RANGE_EXTENDED];  // old
double gamma_LUT = new (std::nothrow) duoble[RGB_PIXEL_RANGE_EXTENDED];  // new

至于为什么从数组变到指针就没有问题了，这个…我猜是数组一些地方溢出了？但我尝试过加[0, 255]限制，也不能彻底解决啊。

增加 16bit 支持

这个更偏向于体力活，没有太多要写的。只是需要明白一个事情，在 8bit 处理中用到的宏常量RGB_PIXEL_RANGE_EXTENDED，是为了让 8bit 的计算不那么僵硬，增加中间数据的过渡，就和我们要把8bit的图像素材扩展到 10bit 或 16bit 处理是一个道理。这里定义的宏常量值为25501，相当于把 8bit 的计算范围扩大了 100 倍。

事实上，log2(25501)介于 14 与 15 之间。所以对于喂给 AreaResize 的 RGB 输入，只有 16bit（和一般不会有人用的 15bit）能实现比 8bit 更“精细”的处理。（但话说回来，就算 8bit 的中间数据再精细，终归还要化为 8bit，所以 10bit 的优势还是在的？）

增加 32bit 支持与gamma参数

1 月 13 号做的事情包括增加 16bit 与 32bit 支持，增加gamma参数。

先思考一下为要查表。就 AreaResize 而言，由于对 RGB 格式增加了 Gamma 校正（但我不太明白原作者为什么不给 YUV 格式增加 Gamma 校正？是一般的滤镜都不做的惯例吗？），计算量增大，不如把 Gamma 校正后的颜色值用查表的形式给出。

但是，对 32bit 浮点数数据，由于数据过于精细，计算表格本身的计算量过大，已经不适合再查表了。而且查表的初衷是为了做 Gamma 校正，而无论是我查的一些资料，还是我的实际测试，在 32bit 下，由于精度足够高，做不做 Gamma 校正都没什么区别。

所以在 32bit 的处理上，没有加 Gamma 校正，也没有做成查表的形式，直接是用到的时候再计算。

另外，给 8~16bit RGB 增加了 Gamma 校正参数，让使用者能够自行调整 Gamma 校正的数值。因为 Gamma 校正具体用什么值，本身似乎就是需要商榷的。

至此，面上的问题都解决了。然后整理了一下代码，就发布了。

未来要做的

bool 型的函数

其实我在迁移的过程中，一直都回避了一个问题。在 avs 版代码中，核心处理函数Ver()和Hor()都是bool型，在调用过程中会借此给一些特殊情况单独处理，以增强鲁棒性（在最开始，我以为这个判断是防止内存溢出的）。但我为了快点迁移代码，把这个问题暂时忽略了。

代码优化了什么

虽然说，优化前后的两版 avs 代码，都不难读懂，但这种优化的思路，值得我去学习。这是把实际项目和算法题联系起来的桥梁。况且，为什么换成优化后的代码，就没有“挑剔输入尺寸”的 Bug，这个问题我还不知道。

一些优化

没什么用的指令集优化

在 AreaResize 中，接受 8bit 和 16bit RGB 输入时，要进行查表操作。联想到 JincResize 的经验，能否在计算 table 的时候进行指令集优化，以加快速度。

于是风风火火写了相关代码。但代码写完，冷静下来才注意到，AreaResize 所需的表格很简单。相比 JincResize 中由两个向量（数组）得到 table，AreaResize 中计算 table 仅需一些算数乘除法，好像不需要我再专门优化，因为由于很简单…或许编译器本身就给优化了？

特意进行了稍长一点的测试，发现 8bit 下优化之后的代码稍快，而 16bit 则是优化前的稍快。所以…好像没什么必要吧。

或许画蛇添足的 Parallel

因为偶然看到《OpenCL 2.0 异构计算》这本书，读到了 C++ AMP（已经成为 C++11 标准的一个扩展），明白了parallel_for的含义，可以简单理解为将 for 循环并行化。

由于有不同的实现方式，所以多少让人感到凌乱。这里用了微软的ppl.h（已经整合到了 msvc 中），至于 gcc…先放着吧…Orz

1
2
3
4
5
6
7
8
9
10
#if defined(_MSC_VER)
#include <ppl.h>

Concurrency::parallel_for(0, (int)range, [&](int i)
// for (int i = 0; i < range; i++)
{
    ...
});

#endif

测试发现加了上述优化后，速度反而降了…Orz…然后我才注意到，VapourSynth 本身已经有fmParallel参数，应该是已经包含了 Parallel 相关优化？所以我这样做是画蛇添足？

目前 Github 上 repo 的情况，Parallel 部分提交到了主分支，而 AVX 部分则是提交到了新建的 dev 分支。

总之，这应该是我走出前人已有算法，做出属于我的微小贡献的第一步吧。