2020-07-06

PyTorchで画像縮小アルゴリズム

pytorch で area average な画像縮小を試してみる．

概要

画像のサイズを変更する画像補間手法は多々存在する．しかしながら，一般によく知られている bilinear や bicubic，nearest-neighbor を使って画像縮小をすると，巨大な画像を小さくする場合ではジャギーが目立つ．したがって，綺麗に画像縮小を行いたい場合は area average を用いるのがよい．本稿では pytorch で画像縮小を行う方法について簡単なコードで見ていく．

area averageとは

area average は，サンプリング点同士の間を等分して得られる領域の平均の画素値を計算するアルゴリズムである． bilinearやbicubic は入力画像中でサンプリング点座標近傍の4点や16点しか見ない．一方で，area average であればサンプリング点がたとえば100画素ぶん離れていれば 100x100 の領域の平均値を出力画素の画素値とする．

f:id:yuuhopro:20200706183208j:plain — アルゴリズムの違い

要するに average pooling である．サンプリング点が離れれば離れるほど， bicubic よりも area を使うほうが良い感じになる．

pytorchの画像リサイズ関数

pytorch では画像をリサイズするための torch.nn.functional.interpolate() が用意されている．使用できるアルゴリズムは，現時点(ver1.5.1)では

nearest
linear
bilinear
bicubic
trilinear
area

であり，おそらく area が area average に該当する． pytorch 最新版documentation

適当に動かして確認する

まずは画像を適当に用意する．3024x3024のデカイ画像を用意した．デカくて載らないので圧縮済みのものを載せておく．

f:id:yuuhopro:20200706180523j:plain — 元画像 3024x3024 掲載のものはjpg品質50に圧縮済み

コード

import numpy as np
import cv2
import torch
import torch.nn.functional as F


if __name__ == '__main__':
    img = cv2.imread('curry.png',cv2.IMREAD_COLOR)

    # pytorchに持ってくる (1,C,H,W)
    tensor = torch.from_numpy(img.astype(np.float32)).permute(2,0,1)[None,:,:,:]
    
    # area average での縮小
    out = F.interpolate(tensor,size=(256,256),mode='area')
    ndarr = out[0].permute(1,2,0).numpy().astype(np.uint8)
    cv2.imwrite('out_area.png',ndarr)

    # bicubic での縮小
    out = F.interpolate(tensor,size=(256,256),mode='bicubic',align_corners=False).clamp(0.0,255.0)
    ndarr = out[0].permute(1,2,0).numpy().astype(np.uint8)
    cv2.imwrite('out_cubic.png',ndarr)

結果

f:id:yuuhopro:20200706180922p:plain — area average 縮小

f:id:yuuhopro:20200706180945p:plain — bicubic 縮小

この写真ではよくわからないかもしれないが，皿のフチの部分がbicubicはジャギっているのがわかる． 128x128 などもっと低解像度に縮小すればより顕著に現れる．

著作権の関係で載せないが，アニメや漫画の画像であればよりわかりやすくこの現象が見てとれる．

結論

深層学習などで画像を縮小して読み込まなければいけない場合， bicubicなどを使うのではなく area average を使うべきである．

2020-06-16

numpyで積分画像を作る

numpyで積分画像を作る際は img.cumsum(axis=0).cumsum(axis=1) が良いよという話。

積分画像 (integral image) とは

数式的表記

画像 ${I(x,y)}$ があるとき、その積分画像 ${ii(x,y)}$ は次のように表される。

${ \displaystyle ii(x,y) = \sum_{w=1}^{x} \sum_{h=1}^{y} I(w,h) }$

となる。ここで画像 ${I(x,y)}$ は左上の画素を原点 ${(1,1)}$ として右方向に ${x}$ 軸を、下方向に ${y}$ 軸をとったものであることに注意。

numpy的表記

わかりやすいように画像を白黒画像(1チャンネル)として、 ${(H,W)-}$ 型ndarrayとする。すなわち img[0,0] は画像左上の画素、 img[H-1,W-1] は右下の画素、img[0,W-1] は右上の画素を表す。

積分画像は次の図のように定義され、積分画像における各画素の値は、元の画像で左上の画素からその画素までの区間の総和となる。つまり ii[y,x] == img[:y+1,:x+1].sum() の関係が成り立つ。

f:id:yuuhopro:20200616205555j:plain — 積分画像の説明図

応用

これがいったい何の役に立つかというと、主に画像中の矩形領域の画素の総和を計算するのに使える。 (矩形はここでは全ての辺が ${x}$ 軸または ${y}$ 軸に平行な長方形(正方形)のことをいう。 AABB: Axis Aligned Bounding Box )

具体的に考えてみる。

f:id:yuuhopro:20200616203207j:plain — 画像中の矩形

この画像中の矩形箇所(図中斜線領域)の画素値の総和を求める際に、通常であれば二重のfor文で全ての画素値を足し合わせることになる。

つまり

def area_sum(img,x1,x2,y1,y2):
    result = 0
    for h in range(y1,y2):
        for w in range(x1,x2):
            result += img[h,w]
    return result

または、numpyの機能を使うと

def area_sum(img,x1,x2,y1,y2):
    return img[y1:y2,x1:x2].sum()

となる。いずれにしてもアルゴリズム上では同じ計算量である。 ${x_2-x_1=N}$ 、 ${y_2-y_1=M}$ とすると計算量は ${O(NM)}$ である。

一方で積分画像を予め用意しておいた場合は、

def area_sum(ii, x1,x2,y1,y2):
    return ii[y2,x2] + ii[y1,x1] - ii[y2,x1] - ii[y1,x2]

という ${O(1)}$ のアルゴリズムで求まる。これは、図形的に表すと以下のようになる。

f:id:yuuhopro:20200616212435j:plain — 矩形領域の計算方法

計算で使用している全ての項が、左上から数えた画素値の和つまり積分画像に記憶されている値になっていることがミソである。これにより簡単に画像中の任意の矩形範囲の画素値の和を高速に求めることができるようになった。これは求めたい範囲が複数個ある場合に非常に有効になる。一箇所しか求めない場合は積分画像を作成する際の計算量でメリットが打ち消されてしまうが、調べたい矩形が何千箇所もの数になる場合は圧倒的な差である。

積分画像の作成

積分画像作成の愚直な実装では以下のようになる。

ii = np.empty((H,W))
for h in range(H):
    for w in range(W):
        ii[h,w] = img[:h+1,:w+1].sum()

当然このような素人のアルゴリズムは採用すべきではない。通常であれば累積和により求めることが妥当である。つまり以下のようにすると .sum() の部分の二乗オーダーの計算が省略でき、高速である。

ii = np.empty((H,W))
ii[0,0] = img[0,0]
for w in range(1,W):
    ii[0,w] = ii[0,w-1] + img[0,w]
for h in range(1,H):
    ii[h,0] = ii[h-1,0] + img[h,0]
    for w in range(1,W):
        ii[h,w] = ii[h,w-1] + img[h,w]

しかしながら、pythonを使う以上はfor文に依る低速化は無視できない。そこで numpy の cumsum() 関数により、各列での累積和を計算し、その値を用いて各行の累積和を計算する。

まず内部的な動作を示す。以下のように各列ごとの累積和を求め、その後に各行ごとの累積和を求めるようにすると自動的に積分画像ができる。

ii = np.empty((H,W))

# cumsum(axis=0)
for w in range(W):
    ii[0,w] = img[0,w]
    for h in range(1,H):
        ii[h,w] = ii[h-1,w] + img[h,w]

# cumsum(axis=1)
for h in range(H):
    ii[h,0] = ii[h,0] + img[h,0]
    for w in range(1,W):
        ii[h,w] = ii[h,w-1] + img[h,w]

これは単純に

ii = img.cumsum(axis=0).cumsum(axis=1)

によって計算可能であり、完全にfor文を排して高速に計算できる。

おまけ cv2.integral, cv2.integral2 について

画像処理、画像認識に使えるライブラリの一つである cv2(OpenCV) は積分画像計算用の関数を備えている。

ii = cv2.integral(img)

で積分画像が求まるが、numpyの cumsum() を用いたものと違い、最初に 0 が挿入される。

つまり [ 1,1,1,1,1 ] という配列に対してそれぞれの関数を適用すると、

cumsum()適用: 出力 [1,2,3,4,5]
integral()適用: 出力[0,1,2,3,4,5]

となる。したがって integral() では入力サイズと出力サイズは異なる。

また、 cv2.integral2() を利用すれば、

ii, ii2 = cv2.integral2(img)

のように戻り値は二つであり、それぞれは

ii : cv2.integral( img )
ii2 : cv2.integral( img**2 )

となっている。つまり通常の積分画像に加えて、同時に画素値を2乗した場合の積分画像も計算してくれる。

これは、おそらく矩形範囲の画素値の分散などの計算に便利である。

カオスノカケラ

めもちょう