2次元アフィン変換の高速転送(2) - melpon日記 - HaskellもC++もまともに扱えないへたれのページ

2次元アフィン変換の高速転送(1)

前回でアフィン変換を行うための行列クラスを作ったので、今回は実際に転送を行ってみようと思います。

オレオレ用語

転送元と転送先と書くと、パッと見分かりづらいので、それぞれ src, dst と書くことにします。

とりあえず転送してみる。

まず何も考えずに、1ピクセル単位ずつ行列を掛けて転送してみましょう。
元の画像が

こんなのだとして、これを 30 度回転させて 1.5 倍に拡大してみます。

matrix mat;
mat.scale(fixed(3, 2), fixed(3, 2)); // 縦横を 1.5 倍にする
mat.rotate(30); // 30 度回転
for (int y = 0; y < src->width(); y++)
{
    for (int x = 0; x < src->height(); x++)
    {
        fixed_point p = mat.transform(fixed_point(x, y)); // 座標変換
        dst->set_pixel(p.x, p.y, src->get_pixel(x, y));
    }
}

こんな感じになりましたが、このコードには 2 つの大きな問題があります。
1つ目はすぐに分かります。これを実行してみると、表示される画像は以下のようになります。

穴あきになってますね。src のループで dst にピクセルを転送しているので、拡大する場合はどうしても足りない部分が出てきてしまうのです。
これを何とかしないと使い物になりません。

2つ目は、毎回座標変換を行っているので遅いということです。これも何とかしないと使い物になりません。
もちろん set_pixel や get_pixel も遅いのですが、これはポインタを使って書くのが面倒で簡略化しただけなので考えないことにします。

画像が穴あきになってしまう問題

これが起こってしまう原因は、の src ピクセル数より dst のピクセル数の方が多いからです。
src を主体としたループで dst に転送するのは問題があるのです。

なのでこれは src 主体のループではなく、dst を主体としたループに変更してやります。
dst のあるピクセルは src のどのピクセルに対応しているのかというのを計算するのです。
dst 主体のループにするためには、アフィン変換した後の四角形から最小の矩形を出してやる必要があります。
つまり、

この赤いラインの矩形を出す必要があります。これはアフィン変換後の四角形から最小、最大の位置を探すだけです。

matrix mat;
matrix.scale(1.5f, 1.5f);
matrix.rotate(30);

fixed_point ps[] = {
    fixed_point(0, 0),
    fixed_point(0, src->height()),
    fixed_point(src->width(), src->height()),
    fixed_point(src->width(), 0),
};

ps[0] = mat.transform(ps[0]);
ps[1] = mat.transform(ps[1]);
ps[2] = mat.transform(ps[2]);
ps[3] = mat.transform(ps[3]);

// アフィン変換した四角形から、その四角形を包含する最小の矩形を計算
fixed_point lt(
    std::min(std::min(std::min(ps[0].x, ps[1].x), ps[2].x), ps[3].x),
    std::min(std::min(std::min(ps[0].y, ps[1].y), ps[2].y), ps[3].y));
fixed_point rb(
    std::max(std::max(std::max(ps[0].x, ps[1].x), ps[2].x), ps[3].x),
    std::max(std::max(std::max(ps[0].y, ps[1].y), ps[2].y), ps[3].y));

これで dst 主体のループが書けます。

for (int y = lt.y; y <= rb.y; y++)
{
    for (int x = lt.x; x <= rb.x; x++)
    {
        fixed_point p = ... // x, y から対応する src の位置を求める
        dst->set_pixel(x, y, src->get_pixel(p.x, p.y));
    }
}

あとは ... となっている部分（コメントの部分）をどうやって埋めるかです。
src の座標から dst の座標を求めるには、変換行列を $T$ として、
$\LARGE T \Large\begin{pmatrix}x \\ y \\ 1\end{pmatrix}=\begin{pmatrix}x^\prime \\ y^\prime \\ 1\end{pmatrix}$
という計算を行いました。今は dst の座標から src の座標を求めたいので、この逆の計算、つまり $x^\prime$ , $y^\prime$ から $x$ , $y$ を求めたいのです。
ここで使うのが逆行列です。逆行列は $T^{-1}$ と書き、以下のような性質があります。
$\Large T\cdot T^{-1}=T^{-1}\cdot T=E$ ( $E$ は恒等行列)
これを使って、以下のようにして $x^\prime$ , $y^\prime$ から $x$ , $y$ を求めることができます。
$\Large T^{-1}\begin{pmatrix}x^\prime \\ y^\prime \\ 1\end{pmatrix}=T^{-1}T \begin{pmatrix}x \\ y \\ 1\end{pmatrix}=E\begin{pmatrix}x \\ y \\ 1\end{pmatrix}=\begin{pmatrix}x \\ y \\ 1\end{pmatrix}$
これで先ほどのプログラムを完成させることができます。

// 逆行列を求める
matrix inv(mat);
inv.invert();
for (int y = lt.y; y <= rb.y; y++)
{
    for (int x = lt.x; x <= rb.x; x++)
    {
        fixed_point p = inv.transform(fixed_point(x, y)); // x, y から対応する src の位置を求める
        dst->set_pixel(x, y, src->get_pixel(p.x, p.y));
    }
}

図で表すと、以下のようになります。

少し分かりづらいですが、赤い半透明の塗りつぶしが最初の矩形（プログラム中で ps を初期化した段階）です。
これをアフィン変換した四角形が青い半透明の塗りつぶしで、それを内包する最小の矩形が赤い矩形です。
そして、その赤い矩形に逆行列を掛けたものが青い四角形です。赤い塗りつぶしの矩形をちゃんと覆っているのが分かります。
赤い矩形のあるピクセルが、青い四角形のどのピクセルに対応しているかを考えると分かりやすいかもしれません。

あと図を見れば明らかですが、赤い矩形に行列を掛けた後のピクセルは src の領域（青い半透明の塗りつぶし）の外側に行くことがあるので、ちゃんとクリッピングしてやる必要があります。このプログラムでは、get_pixel で src の領域の外側だった場合は透明のピクセルを返していると思ってください。

これで無事、以下のような画像が表示されます。

線形補間等はしていないので綺麗ではないですが、穴あきの問題は解決されていることが分かります。

あとは速度の問題だけです。
→書きました。

2次元アフィン変換の高速転送(3)